Actual source code: ex7.c
2: static char help[] = "Tests DMLocalToLocalxxx() for DMDA.\n\n";
4: #include <petscdmda.h>
6: int main(int argc, char **argv)
7: {
8: PetscMPIInt rank;
9: PetscInt M = 8, dof = 1, stencil_width = 1, i, start, end, P = 5, N = 6, m = PETSC_DECIDE, n = PETSC_DECIDE, p = PETSC_DECIDE, pt = 0, st = 0;
10: PetscBool flg = PETSC_FALSE, flg2, flg3;
11: DMBoundaryType periodic;
12: DMDAStencilType stencil_type;
13: DM da;
14: Vec local, global, local_copy;
15: PetscScalar value;
16: PetscReal norm, work;
17: PetscViewer viewer;
18: char filename[64];
19: FILE *file;
21: PetscFunctionBeginUser;
22: PetscCall(PetscInitialize(&argc, &argv, (char *)0, help));
23: PetscCall(PetscOptionsGetInt(NULL, NULL, "-M", &M, NULL));
24: PetscCall(PetscOptionsGetInt(NULL, NULL, "-N", &N, NULL));
25: PetscCall(PetscOptionsGetInt(NULL, NULL, "-dof", &dof, NULL));
26: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_width", &stencil_width, NULL));
27: PetscCall(PetscOptionsGetInt(NULL, NULL, "-periodic", &pt, NULL));
29: periodic = (DMBoundaryType)pt;
31: PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL));
33: stencil_type = (DMDAStencilType)st;
35: PetscCall(PetscOptionsHasName(NULL, NULL, "-grid2d", &flg2));
36: PetscCall(PetscOptionsHasName(NULL, NULL, "-grid3d", &flg3));
37: if (flg2) {
38: PetscCall(DMDACreate2d(PETSC_COMM_WORLD, periodic, periodic, stencil_type, M, N, m, n, dof, stencil_width, NULL, NULL, &da));
39: } else if (flg3) {
40: PetscCall(DMDACreate3d(PETSC_COMM_WORLD, periodic, periodic, periodic, stencil_type, M, N, P, m, n, p, dof, stencil_width, NULL, NULL, NULL, &da));
41: } else {
42: PetscCall(DMDACreate1d(PETSC_COMM_WORLD, periodic, M, dof, stencil_width, NULL, &da));
43: }
44: PetscCall(DMSetFromOptions(da));
45: PetscCall(DMSetUp(da));
47: PetscCall(DMCreateGlobalVector(da, &global));
48: PetscCall(DMCreateLocalVector(da, &local));
49: PetscCall(VecDuplicate(local, &local_copy));
51: /* zero out vectors so that ghostpoints are zero */
52: value = 0;
53: PetscCall(VecSet(local, value));
54: PetscCall(VecSet(local_copy, value));
56: PetscCall(VecGetOwnershipRange(global, &start, &end));
57: for (i = start; i < end; i++) {
58: value = i + 1;
59: PetscCall(VecSetValues(global, 1, &i, &value, INSERT_VALUES));
60: }
61: PetscCall(VecAssemblyBegin(global));
62: PetscCall(VecAssemblyEnd(global));
64: PetscCall(DMGlobalToLocalBegin(da, global, INSERT_VALUES, local));
65: PetscCall(DMGlobalToLocalEnd(da, global, INSERT_VALUES, local));
67: PetscCall(DMLocalToLocalBegin(da, local, INSERT_VALUES, local_copy));
68: PetscCall(DMLocalToLocalEnd(da, local, INSERT_VALUES, local_copy));
70: PetscCall(PetscOptionsGetBool(NULL, NULL, "-save", &flg, NULL));
71: if (flg) {
72: PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank));
73: PetscCall(PetscSNPrintf(filename, PETSC_STATIC_ARRAY_LENGTH(filename), "local.%d", rank));
74: PetscCall(PetscViewerASCIIOpen(PETSC_COMM_SELF, filename, &viewer));
75: PetscCall(PetscViewerASCIIGetPointer(viewer, &file));
76: PetscCall(VecView(local, viewer));
77: fprintf(file, "Vector with correct ghost points\n");
78: PetscCall(VecView(local_copy, viewer));
79: PetscCall(PetscViewerDestroy(&viewer));
80: }
82: PetscCall(VecAXPY(local_copy, -1.0, local));
83: PetscCall(VecNorm(local_copy, NORM_MAX, &work));
84: PetscCall(MPIU_Allreduce(&work, &norm, 1, MPIU_REAL, MPIU_MAX, PETSC_COMM_WORLD));
85: PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Norm of difference %g should be zero\n", (double)norm));
87: PetscCall(VecDestroy(&local_copy));
88: PetscCall(VecDestroy(&local));
89: PetscCall(VecDestroy(&global));
90: PetscCall(DMDestroy(&da));
91: PetscCall(PetscFinalize());
92: return 0;
93: }
95: /*TEST
97: test:
98: nsize: 8
99: args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic
101: test:
102: suffix: 2
103: nsize: 8
104: args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic -grid2d
105: output_file: output/ex7_1.out
107: test:
108: suffix: 3
109: nsize: 8
110: args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic -grid3d
111: output_file: output/ex7_1.out
113: TEST*/