Actual source code: ex21.c

  1: static const char help[] = "Test DMCreateInjection() for mapping coordinates in 3D";

  3: #include <petscvec.h>
  4: #include <petscmat.h>
  5: #include <petscdm.h>
  6: #include <petscdmda.h>

  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;

 16:   PetscFunctionBeginUser;
 17:   bx = DM_BOUNDARY_NONE;
 18:   by = DM_BOUNDARY_NONE;
 19:   bz = DM_BOUNDARY_NONE;

 21:   periodicity = 0;

 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:   }

 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));

 37:   PetscCall(DMCoarsen(daf, MPI_COMM_NULL, &dac));

 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));

 42:   {
 43:     DM         cdaf, cdac;
 44:     Vec        coordsc, coordsf, coordsf2;
 45:     Mat        inject;
 46:     VecScatter vscat;
 47:     Mat        interp;
 48:     PetscReal  norm;

 50:     PetscCall(DMGetCoordinateDM(dac, &cdac));
 51:     PetscCall(DMGetCoordinateDM(daf, &cdaf));

 53:     PetscCall(DMGetCoordinates(dac, &coordsc));
 54:     PetscCall(DMGetCoordinates(daf, &coordsf));

 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));

 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:   }

 73:   if (0) {
 74:     PetscCall(DMCreateGlobalVector(dac, &ac));

 76:     PetscCall(DMCreateGlobalVector(daf, &af));

 78:     PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "dac_7.vtu", &vv));
 79:     PetscCall(VecView(ac, vv));
 80:     PetscCall(PetscViewerDestroy(&vv));

 82:     PetscCall(PetscViewerASCIIOpen(PETSC_COMM_WORLD, "daf_7.vtu", &vv));
 83:     PetscCall(VecView(af, vv));
 84:     PetscCall(PetscViewerDestroy(&vv));
 85:     PetscCall(VecDestroy(&ac));
 86:     PetscCall(VecDestroy(&af));
 87:   }
 88:   PetscCall(DMDestroy(&dac));
 89:   PetscCall(DMDestroy(&daf));
 90:   PetscFunctionReturn(PETSC_SUCCESS);
 91: }

 93: int main(int argc, char **argv)
 94: {
 95:   PetscInt mx, my, mz;

 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: }

110: /*TEST

112:       test:
113:          nsize: 5
114:          args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_x 5
115:          output_file: output/empty.out

117:       test:
118:          suffix: 2
119:          nsize: 5
120:          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_x 5
121:          output_file: output/empty.out

123:       test:
124:          suffix: 3
125:          nsize: 5
126:          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_x 5
127:          output_file: output/empty.out

129:       test:
130:          suffix: 4
131:          nsize: 5
132:          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_x 5
133:          output_file: output/empty.out

135:       test:
136:          suffix: 5
137:          nsize: 5
138:          args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_y 5
139:          output_file: output/empty.out

141:       test:
142:          suffix: 6
143:          nsize: 5
144:          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_y 5
145:          output_file: output/empty.out

147:       test:
148:          suffix: 7
149:          nsize: 5
150:          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_y 5
151:          output_file: output/empty.out

153:       test:
154:          suffix: 8
155:          nsize: 5
156:          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_y 5
157:          output_file: output/empty.out

159:       test:
160:          suffix: 9
161:          nsize: 5
162:          args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_z 5
163:          output_file: output/empty.out

165:       test:
166:          suffix: 10
167:          nsize: 5
168:          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_z 5
169:          output_file: output/empty.out

171:       test:
172:          suffix: 11
173:          nsize: 5
174:          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_z 5
175:          output_file: output/empty.out

177:       test:
178:          suffix: 12
179:          nsize: 5
180:          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_z 5
181:          output_file: output/empty.out

183:       test:
184:          suffix: 13
185:          nsize: 5
186:          args: -mx 30 -my 30 -mz 30 -periodic 0
187:          output_file: output/empty.out

189:       test:
190:          suffix: 14
191:          nsize: 5
192:          args: -mx 29 -my 30 -mz 30 -periodic 1
193:          output_file: output/empty.out

195:       test:
196:          suffix: 15
197:          nsize: 5
198:          args: -mx 30 -my 29 -mz 30 -periodic 2
199:          output_file: output/empty.out

201:       test:
202:          suffix: 16
203:          nsize: 5
204:          args: -mx 30 -my 30 -mz 29 -periodic 3
205:          output_file: output/empty.out

207: TEST*/