1 2 static char help[] = "Tests DMLocalToLocalxxx() for DMDA.\n\n"; 3 4 #include <petscdmda.h> 5 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; 20 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)); 28 29 periodic = (DMBoundaryType)pt; 30 31 PetscCall(PetscOptionsGetInt(NULL, NULL, "-stencil_type", &st, NULL)); 32 33 stencil_type = (DMDAStencilType)st; 34 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)); 46 47 PetscCall(DMCreateGlobalVector(da, &global)); 48 PetscCall(DMCreateLocalVector(da, &local)); 49 PetscCall(VecDuplicate(local, &local_copy)); 50 51 /* zero out vectors so that ghostpoints are zero */ 52 value = 0; 53 PetscCall(VecSet(local, value)); 54 PetscCall(VecSet(local_copy, value)); 55 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)); 63 64 PetscCall(DMGlobalToLocalBegin(da, global, INSERT_VALUES, local)); 65 PetscCall(DMGlobalToLocalEnd(da, global, INSERT_VALUES, local)); 66 67 PetscCall(DMLocalToLocalBegin(da, local, INSERT_VALUES, local_copy)); 68 PetscCall(DMLocalToLocalEnd(da, local, INSERT_VALUES, local_copy)); 69 70 PetscCall(PetscOptionsGetBool(NULL, NULL, "-save", &flg, NULL)); 71 if (flg) { 72 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 73 sprintf(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 } 81 82 PetscCall(VecAXPY(local_copy, -1.0, local)); 83 PetscCall(VecNorm(local_copy, NORM_MAX, &work)); 84 PetscCallMPI(MPI_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)); 86 87 PetscCall(VecDestroy(&local_copy)); 88 PetscCall(VecDestroy(&local)); 89 PetscCall(VecDestroy(&global)); 90 PetscCall(DMDestroy(&da)); 91 PetscCall(PetscFinalize()); 92 return 0; 93 } 94 95 /*TEST 96 97 test: 98 nsize: 8 99 args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic 100 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 106 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 112 113 TEST*/ 114