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