Actual source code: ex33.c
2: static char help[] = "Tests the routines VecScatterCreateToAll(), VecScatterCreateToZero()\n\n";
4: #include <petscvec.h>
6: int main(int argc, char **argv)
7: {
8: PetscInt n = 3, i, len, start, end;
9: PetscMPIInt size, rank;
10: PetscScalar value, *yy;
11: Vec x, y, z, y_t;
12: VecScatter toall, tozero;
15: PetscInitialize(&argc, &argv, (char *)0, help);
16: MPI_Comm_size(PETSC_COMM_WORLD, &size);
17: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
19: /* create two vectors */
20: VecCreateMPI(PETSC_COMM_WORLD, PETSC_DECIDE, size * n, &x);
22: /* each processor inserts its values */
24: VecGetOwnershipRange(x, &start, &end);
25: for (i = start; i < end; i++) {
26: value = (PetscScalar)i;
27: VecSetValues(x, 1, &i, &value, INSERT_VALUES);
28: }
29: VecAssemblyBegin(x);
30: VecAssemblyEnd(x);
31: VecView(x, PETSC_VIEWER_STDOUT_WORLD);
33: VecScatterCreateToAll(x, &toall, &y);
34: VecScatterBegin(toall, x, y, INSERT_VALUES, SCATTER_FORWARD);
35: VecScatterEnd(toall, x, y, INSERT_VALUES, SCATTER_FORWARD);
36: VecScatterDestroy(&toall);
38: /* Cannot view the above vector with VecView(), so place it in an MPI Vec
39: and do a VecView() */
40: VecGetArray(y, &yy);
41: VecGetLocalSize(y, &len);
42: VecCreateMPIWithArray(PETSC_COMM_WORLD, 1, len, PETSC_DECIDE, yy, &y_t);
43: VecView(y_t, PETSC_VIEWER_STDOUT_WORLD);
44: VecDestroy(&y_t);
45: VecRestoreArray(y, &yy);
47: VecScatterCreateToAll(x, &tozero, &z);
48: VecScatterBegin(tozero, x, z, INSERT_VALUES, SCATTER_FORWARD);
49: VecScatterEnd(tozero, x, z, INSERT_VALUES, SCATTER_FORWARD);
50: VecScatterDestroy(&tozero);
51: if (rank == 0) VecView(z, PETSC_VIEWER_STDOUT_SELF);
52: VecDestroy(&z);
54: VecScatterCreateToZero(x, &tozero, &z);
55: VecScatterBegin(tozero, x, z, INSERT_VALUES, SCATTER_FORWARD);
56: VecScatterEnd(tozero, x, z, INSERT_VALUES, SCATTER_FORWARD);
57: VecScatterDestroy(&tozero);
58: VecDestroy(&z);
60: VecDestroy(&x);
61: VecDestroy(&y);
63: PetscFinalize();
64: return 0;
65: }
67: /*TEST
69: test:
70: nsize: 4
72: TEST*/