Actual source code: ex13.c
2: static char help[] = "Scatters from a sequential vector to a parallel vector. \n\
3: uses block index sets\n\n";
5: #include <petscvec.h>
7: int main(int argc, char **argv)
8: {
9: PetscInt bs = 1, n = 5, i, low;
10: PetscInt ix0[3] = {5, 7, 9}, iy0[3] = {1, 2, 4}, ix1[3] = {2, 3, 4}, iy1[3] = {0, 1, 3};
11: PetscMPIInt size, rank;
12: PetscScalar *array;
13: Vec x, y;
14: IS isx, isy;
15: VecScatter ctx;
16: PetscViewer sviewer;
19: PetscInitialize(&argc, &argv, (char *)0, help);
20: MPI_Comm_size(PETSC_COMM_WORLD, &size);
21: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
25: PetscOptionsGetInt(NULL, NULL, "-bs", &bs, NULL);
26: n = bs * n;
28: /* Create vector x over shared memory */
29: VecCreate(PETSC_COMM_WORLD, &x);
30: VecSetSizes(x, n, PETSC_DECIDE);
31: VecSetFromOptions(x);
33: VecGetOwnershipRange(x, &low, NULL);
34: VecGetArray(x, &array);
35: for (i = 0; i < n; i++) array[i] = (PetscScalar)(i + low);
36: VecRestoreArray(x, &array);
38: /* Create a sequential vector y */
39: VecCreateSeq(PETSC_COMM_SELF, n, &y);
40: VecGetArray(y, &array);
41: for (i = 0; i < n; i++) array[i] = (PetscScalar)(i + 100 * rank);
42: VecRestoreArray(y, &array);
44: /* Create two index sets */
45: if (rank == 0) {
46: ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix0, PETSC_COPY_VALUES, &isx);
47: ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy0, PETSC_COPY_VALUES, &isy);
48: } else {
49: ISCreateBlock(PETSC_COMM_SELF, bs, 3, ix1, PETSC_COPY_VALUES, &isx);
50: ISCreateBlock(PETSC_COMM_SELF, bs, 3, iy1, PETSC_COPY_VALUES, &isy);
51: }
53: if (rank == 10) {
54: PetscPrintf(PETSC_COMM_SELF, "\n[%d] isx:\n", rank);
55: ISView(isx, PETSC_VIEWER_STDOUT_SELF);
56: }
58: VecScatterCreate(y, isy, x, isx, &ctx);
59: VecScatterSetFromOptions(ctx);
61: /* Test forward vecscatter */
62: VecSet(x, 0.0);
63: VecScatterBegin(ctx, y, x, ADD_VALUES, SCATTER_FORWARD);
64: VecScatterEnd(ctx, y, x, ADD_VALUES, SCATTER_FORWARD);
65: VecView(x, PETSC_VIEWER_STDOUT_WORLD);
67: /* Test reverse vecscatter */
68: VecScale(x, -1.0);
69: VecScatterBegin(ctx, x, y, ADD_VALUES, SCATTER_REVERSE);
70: VecScatterEnd(ctx, x, y, ADD_VALUES, SCATTER_REVERSE);
71: PetscViewerGetSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
72: if (rank == 1) VecView(y, sviewer);
73: PetscViewerRestoreSubViewer(PETSC_VIEWER_STDOUT_WORLD, PETSC_COMM_SELF, &sviewer);
75: /* Free spaces */
76: VecScatterDestroy(&ctx);
77: ISDestroy(&isx);
78: ISDestroy(&isy);
79: VecDestroy(&x);
80: VecDestroy(&y);
81: PetscFinalize();
82: return 0;
83: }
85: /*TEST
87: test:
88: nsize: 3
90: TEST*/