xref: /petsc/src/dm/tests/ex21.c (revision ebead697dbf761eb322f829370bbe90b3bd93fa3)
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(0);
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 
118       test:
119          suffix: 2
120          nsize: 5
121          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_x 5
122 
123       test:
124          suffix: 3
125          nsize: 5
126          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_x 5
127 
128       test:
129          suffix: 4
130          nsize: 5
131          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_x 5
132 
133       test:
134          suffix: 5
135          nsize: 5
136          args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_y 5
137 
138       test:
139          suffix: 6
140          nsize: 5
141          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_y 5
142 
143       test:
144          suffix: 7
145          nsize: 5
146          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_y 5
147 
148       test:
149          suffix: 8
150          nsize: 5
151          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_y 5
152 
153       test:
154          suffix: 9
155          nsize: 5
156          args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_z 5
157 
158       test:
159          suffix: 10
160          nsize: 5
161          args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_z 5
162 
163       test:
164          suffix: 11
165          nsize: 5
166          args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_z 5
167 
168       test:
169          suffix: 12
170          nsize: 5
171          args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_z 5
172 
173       test:
174          suffix: 13
175          nsize: 5
176          args: -mx 30 -my 30 -mz 30 -periodic 0
177 
178       test:
179          suffix: 14
180          nsize: 5
181          args: -mx 29 -my 30 -mz 30 -periodic 1
182 
183       test:
184          suffix: 15
185          nsize: 5
186          args: -mx 30 -my 29 -mz 30 -periodic 2
187 
188       test:
189          suffix: 16
190          nsize: 5
191          args: -mx 30 -my 30 -mz 29 -periodic 3
192 
193 TEST*/
194