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