1 static const char help[] = "Test DMCreateInjection() for mapping coordinates in 3D"; 2 3 #include <petscvec.h> 4 #include <petscmat.h> 5 #include <petscdm.h> 6 #include <petscdmda.h> 7 8 PetscErrorCode test1_DAInjection3d(PetscInt mx, PetscInt my, PetscInt mz) { 9 DM dac, daf; 10 PetscViewer vv; 11 Vec ac, af; 12 PetscInt periodicity; 13 DMBoundaryType bx, by, bz; 14 15 PetscFunctionBeginUser; 16 bx = DM_BOUNDARY_NONE; 17 by = DM_BOUNDARY_NONE; 18 bz = DM_BOUNDARY_NONE; 19 20 periodicity = 0; 21 22 PetscCall(PetscOptionsGetInt(NULL, NULL, "-periodic", &periodicity, NULL)); 23 if (periodicity == 1) { 24 bx = DM_BOUNDARY_PERIODIC; 25 } else if (periodicity == 2) { 26 by = DM_BOUNDARY_PERIODIC; 27 } else if (periodicity == 3) { 28 bz = DM_BOUNDARY_PERIODIC; 29 } 30 31 PetscCall(DMDACreate3d(PETSC_COMM_WORLD, bx, by, bz, DMDA_STENCIL_BOX, mx + 1, my + 1, mz + 1, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, /* 1 dof */ 32 1, /* stencil = 1 */ NULL, NULL, NULL, &daf)); 33 PetscCall(DMSetFromOptions(daf)); 34 PetscCall(DMSetUp(daf)); 35 36 PetscCall(DMCoarsen(daf, MPI_COMM_NULL, &dac)); 37 38 PetscCall(DMDASetUniformCoordinates(dac, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0)); 39 PetscCall(DMDASetUniformCoordinates(daf, -1.0, 1.0, -1.0, 1.0, -1.0, 1.0)); 40 41 { 42 DM cdaf, cdac; 43 Vec coordsc, coordsf, coordsf2; 44 Mat inject; 45 VecScatter vscat; 46 Mat interp; 47 PetscReal norm; 48 49 PetscCall(DMGetCoordinateDM(dac, &cdac)); 50 PetscCall(DMGetCoordinateDM(daf, &cdaf)); 51 52 PetscCall(DMGetCoordinates(dac, &coordsc)); 53 PetscCall(DMGetCoordinates(daf, &coordsf)); 54 55 PetscCall(DMCreateInjection(cdac, cdaf, &inject)); 56 PetscCall(MatScatterGetVecScatter(inject, &vscat)); 57 PetscCall(VecScatterBegin(vscat, coordsf, coordsc, INSERT_VALUES, SCATTER_FORWARD)); 58 PetscCall(VecScatterEnd(vscat, coordsf, coordsc, INSERT_VALUES, SCATTER_FORWARD)); 59 PetscCall(MatDestroy(&inject)); 60 61 PetscCall(DMCreateInterpolation(cdac, cdaf, &interp, NULL)); 62 PetscCall(VecDuplicate(coordsf, &coordsf2)); 63 PetscCall(MatInterpolate(interp, coordsc, coordsf2)); 64 PetscCall(VecAXPY(coordsf2, -1.0, coordsf)); 65 PetscCall(VecNorm(coordsf2, NORM_MAX, &norm)); 66 /* The fine coordinates are only reproduced in certain cases */ 67 if (!bx && !by && !bz && norm > PETSC_SQRT_MACHINE_EPSILON) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Norm %g\n", (double)norm)); 68 PetscCall(VecDestroy(&coordsf2)); 69 PetscCall(MatDestroy(&interp)); 70 } 71 72 if (0) { 73 PetscCall(DMCreateGlobalVector(dac, &ac)); 74 PetscCall(VecZeroEntries(ac)); 75 76 PetscCall(DMCreateGlobalVector(daf, &af)); 77 PetscCall(VecZeroEntries(af)); 78 79 PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "dac_7.vtu", &vv)); 80 PetscCall(VecView(ac, vv)); 81 PetscCall(PetscViewerDestroy(&vv)); 82 83 PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "daf_7.vtu", &vv)); 84 PetscCall(VecView(af, vv)); 85 PetscCall(PetscViewerDestroy(&vv)); 86 PetscCall(VecDestroy(&ac)); 87 PetscCall(VecDestroy(&af)); 88 } 89 PetscCall(DMDestroy(&dac)); 90 PetscCall(DMDestroy(&daf)); 91 PetscFunctionReturn(0); 92 } 93 94 int main(int argc, char **argv) { 95 PetscInt mx, my, mz; 96 97 PetscFunctionBeginUser; 98 PetscCall(PetscInitialize(&argc, &argv, 0, help)); 99 mx = 2; 100 my = 2; 101 mz = 2; 102 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mx", &mx, 0)); 103 PetscCall(PetscOptionsGetInt(NULL, NULL, "-my", &my, 0)); 104 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mz", &mz, 0)); 105 PetscCall(test1_DAInjection3d(mx, my, mz)); 106 PetscCall(PetscFinalize()); 107 return 0; 108 } 109 110 /*TEST 111 112 test: 113 nsize: 5 114 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_x 5 115 116 test: 117 suffix: 2 118 nsize: 5 119 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_x 5 120 121 test: 122 suffix: 3 123 nsize: 5 124 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_x 5 125 126 test: 127 suffix: 4 128 nsize: 5 129 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_x 5 130 131 test: 132 suffix: 5 133 nsize: 5 134 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_y 5 135 136 test: 137 suffix: 6 138 nsize: 5 139 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_y 5 140 141 test: 142 suffix: 7 143 nsize: 5 144 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_y 5 145 146 test: 147 suffix: 8 148 nsize: 5 149 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_y 5 150 151 test: 152 suffix: 9 153 nsize: 5 154 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_z 5 155 156 test: 157 suffix: 10 158 nsize: 5 159 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_z 5 160 161 test: 162 suffix: 11 163 nsize: 5 164 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_z 5 165 166 test: 167 suffix: 12 168 nsize: 5 169 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_z 5 170 171 test: 172 suffix: 13 173 nsize: 5 174 args: -mx 30 -my 30 -mz 30 -periodic 0 175 176 test: 177 suffix: 14 178 nsize: 5 179 args: -mx 29 -my 30 -mz 30 -periodic 1 180 181 test: 182 suffix: 15 183 nsize: 5 184 args: -mx 30 -my 29 -mz 30 -periodic 2 185 186 test: 187 suffix: 16 188 nsize: 5 189 args: -mx 30 -my 30 -mz 29 -periodic 3 190 191 TEST*/ 192