Actual source code: ex51.c
2: static char help[] = "Tests DMDAGlobalToNaturalAllCreate() using contour plotting for 2d DMDAs.\n\n";
4: #include <petscdm.h>
5: #include <petscdmda.h>
6: #include <petscdraw.h>
8: int main(int argc, char **argv)
9: {
10: PetscInt i, j, M = 10, N = 8, m = PETSC_DECIDE, n = PETSC_DECIDE;
11: PetscMPIInt rank;
12: PetscBool flg = PETSC_FALSE;
13: DM da;
14: PetscViewer viewer;
15: Vec localall, global;
16: PetscScalar value, *vlocal;
17: DMBoundaryType bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE;
18: DMDAStencilType stype = DMDA_STENCIL_BOX;
19: VecScatter tolocalall, fromlocalall;
20: PetscInt start, end;
23: PetscInitialize(&argc, &argv, (char *)0, help);
24: PetscViewerDrawOpen(PETSC_COMM_WORLD, 0, "", 300, 0, 300, 300, &viewer);
26: /* Read options */
27: PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL);
28: PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL);
29: PetscOptionsGetInt(NULL, NULL, "-m", &m, NULL);
30: PetscOptionsGetInt(NULL, NULL, "-n", &n, NULL);
31: PetscOptionsGetBool(NULL, NULL, "-star_stencil", &flg, NULL);
32: if (flg) stype = DMDA_STENCIL_STAR;
34: /* Create distributed array and get vectors */
35: DMDACreate2d(PETSC_COMM_WORLD, bx, by, stype, M, N, m, n, 1, 1, NULL, NULL, &da);
36: DMSetFromOptions(da);
37: DMSetUp(da);
39: DMCreateGlobalVector(da, &global);
40: VecCreateSeq(PETSC_COMM_SELF, M * N, &localall);
42: MPI_Comm_rank(PETSC_COMM_WORLD, &rank);
43: VecGetOwnershipRange(global, &start, &end);
44: for (i = start; i < end; i++) {
45: value = 5.0 * rank;
46: VecSetValues(global, 1, &i, &value, INSERT_VALUES);
47: }
48: VecView(global, viewer);
50: /*
51: Create Scatter from global DMDA parallel vector to local vector that
52: contains all entries
53: */
54: DMDAGlobalToNaturalAllCreate(da, &tolocalall);
55: DMDANaturalAllToGlobalCreate(da, &fromlocalall);
57: VecScatterBegin(tolocalall, global, localall, INSERT_VALUES, SCATTER_FORWARD);
58: VecScatterEnd(tolocalall, global, localall, INSERT_VALUES, SCATTER_FORWARD);
60: VecGetArray(localall, &vlocal);
61: for (j = 0; j < N; j++) {
62: for (i = 0; i < M; i++) *vlocal++ += i + j * M;
63: }
64: VecRestoreArray(localall, &vlocal);
66: /* scatter back to global vector */
67: VecScatterBegin(fromlocalall, localall, global, INSERT_VALUES, SCATTER_FORWARD);
68: VecScatterEnd(fromlocalall, localall, global, INSERT_VALUES, SCATTER_FORWARD);
70: VecView(global, viewer);
72: /* Free memory */
73: VecScatterDestroy(&tolocalall);
74: VecScatterDestroy(&fromlocalall);
75: PetscViewerDestroy(&viewer);
76: VecDestroy(&localall);
77: VecDestroy(&global);
78: DMDestroy(&da);
79: PetscFinalize();
80: return 0;
81: }
83: /*TEST
85: build:
86: requires: !complex
88: test:
89: nsize: 3
91: TEST*/