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 PetscCall(PetscInitialize(&argc,&argv,0,help)); 100 mx = 2; 101 my = 2; 102 mz = 2; 103 PetscCall(PetscOptionsGetInt(NULL,NULL,"-mx", &mx, 0)); 104 PetscCall(PetscOptionsGetInt(NULL,NULL,"-my", &my, 0)); 105 PetscCall(PetscOptionsGetInt(NULL,NULL,"-mz", &mz, 0)); 106 PetscCall(test1_DAInjection3d(mx,my,mz)); 107 PetscCall(PetscFinalize()); 108 return 0; 109 } 110 111 /*TEST 112 113 test: 114 nsize: 5 115 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_x 5 116 117 test: 118 suffix: 2 119 nsize: 5 120 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_x 5 121 122 test: 123 suffix: 3 124 nsize: 5 125 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_x 5 126 127 test: 128 suffix: 4 129 nsize: 5 130 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_x 5 131 132 test: 133 suffix: 5 134 nsize: 5 135 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_y 5 136 137 test: 138 suffix: 6 139 nsize: 5 140 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_y 5 141 142 test: 143 suffix: 7 144 nsize: 5 145 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_y 5 146 147 test: 148 suffix: 8 149 nsize: 5 150 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_y 5 151 152 test: 153 suffix: 9 154 nsize: 5 155 args: -mx 30 -my 30 -mz 30 -periodic 0 -da_processors_z 5 156 157 test: 158 suffix: 10 159 nsize: 5 160 args: -mx 29 -my 30 -mz 30 -periodic 1 -da_processors_z 5 161 162 test: 163 suffix: 11 164 nsize: 5 165 args: -mx 30 -my 29 -mz 30 -periodic 2 -da_processors_z 5 166 167 test: 168 suffix: 12 169 nsize: 5 170 args: -mx 30 -my 30 -mz 29 -periodic 3 -da_processors_z 5 171 172 test: 173 suffix: 13 174 nsize: 5 175 args: -mx 30 -my 30 -mz 30 -periodic 0 176 177 test: 178 suffix: 14 179 nsize: 5 180 args: -mx 29 -my 30 -mz 30 -periodic 1 181 182 test: 183 suffix: 15 184 nsize: 5 185 args: -mx 30 -my 29 -mz 30 -periodic 2 186 187 test: 188 suffix: 16 189 nsize: 5 190 args: -mx 30 -my 30 -mz 29 -periodic 3 191 192 TEST*/ 193