1 2 static char help[] = "Tests DMLocalToLocalxxx() for DMDA.\n\n"; 3 4 #include <petscdmda.h> 5 6 int main(int argc,char **argv) 7 { 8 PetscMPIInt rank; 9 PetscInt M=8,dof=1,stencil_width=1,i,start,end,P=5,N = 6,m=PETSC_DECIDE,n=PETSC_DECIDE,p=PETSC_DECIDE,pt = 0,st = 0; 10 PetscBool flg = PETSC_FALSE,flg2,flg3; 11 DMBoundaryType periodic; 12 DMDAStencilType stencil_type; 13 DM da; 14 Vec local,global,local_copy; 15 PetscScalar value; 16 PetscReal norm,work; 17 PetscViewer viewer; 18 char filename[64]; 19 FILE *file; 20 21 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 22 PetscCall(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 23 PetscCall(PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL)); 24 PetscCall(PetscOptionsGetInt(NULL,NULL,"-dof",&dof,NULL)); 25 PetscCall(PetscOptionsGetInt(NULL,NULL,"-stencil_width",&stencil_width,NULL)); 26 PetscCall(PetscOptionsGetInt(NULL,NULL,"-periodic",&pt,NULL)); 27 28 periodic = (DMBoundaryType) pt; 29 30 PetscCall(PetscOptionsGetInt(NULL,NULL,"-stencil_type",&st,NULL)); 31 32 stencil_type = (DMDAStencilType) st; 33 34 PetscCall(PetscOptionsHasName(NULL,NULL,"-grid2d",&flg2)); 35 PetscCall(PetscOptionsHasName(NULL,NULL,"-grid3d",&flg3)); 36 if (flg2) { 37 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,periodic,periodic,stencil_type,M,N,m,n,dof,stencil_width,NULL,NULL,&da)); 38 } else if (flg3) { 39 PetscCall(DMDACreate3d(PETSC_COMM_WORLD,periodic,periodic,periodic,stencil_type,M,N,P,m,n,p,dof,stencil_width,NULL,NULL,NULL,&da)); 40 } else { 41 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,periodic,M,dof,stencil_width,NULL,&da)); 42 } 43 PetscCall(DMSetFromOptions(da)); 44 PetscCall(DMSetUp(da)); 45 46 PetscCall(DMCreateGlobalVector(da,&global)); 47 PetscCall(DMCreateLocalVector(da,&local)); 48 PetscCall(VecDuplicate(local,&local_copy)); 49 50 /* zero out vectors so that ghostpoints are zero */ 51 value = 0; 52 PetscCall(VecSet(local,value)); 53 PetscCall(VecSet(local_copy,value)); 54 55 PetscCall(VecGetOwnershipRange(global,&start,&end)); 56 for (i=start; i<end; i++) { 57 value = i + 1; 58 PetscCall(VecSetValues(global,1,&i,&value,INSERT_VALUES)); 59 } 60 PetscCall(VecAssemblyBegin(global)); 61 PetscCall(VecAssemblyEnd(global)); 62 63 PetscCall(DMGlobalToLocalBegin(da,global,INSERT_VALUES,local)); 64 PetscCall(DMGlobalToLocalEnd(da,global,INSERT_VALUES,local)); 65 66 PetscCall(DMLocalToLocalBegin(da,local,INSERT_VALUES,local_copy)); 67 PetscCall(DMLocalToLocalEnd(da,local,INSERT_VALUES,local_copy)); 68 69 PetscCall(PetscOptionsGetBool(NULL,NULL,"-save",&flg,NULL)); 70 if (flg) { 71 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 72 sprintf(filename,"local.%d",rank); 73 PetscCall(PetscViewerASCIIOpen(PETSC_COMM_SELF,filename,&viewer)); 74 PetscCall(PetscViewerASCIIGetPointer(viewer,&file)); 75 PetscCall(VecView(local,viewer)); 76 fprintf(file,"Vector with correct ghost points\n"); 77 PetscCall(VecView(local_copy,viewer)); 78 PetscCall(PetscViewerDestroy(&viewer)); 79 } 80 81 PetscCall(VecAXPY(local_copy,-1.0,local)); 82 PetscCall(VecNorm(local_copy,NORM_MAX,&work)); 83 PetscCallMPI(MPI_Allreduce(&work,&norm,1,MPIU_REAL,MPIU_MAX,PETSC_COMM_WORLD)); 84 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Norm of difference %g should be zero\n",(double)norm)); 85 86 PetscCall(VecDestroy(&local_copy)); 87 PetscCall(VecDestroy(&local)); 88 PetscCall(VecDestroy(&global)); 89 PetscCall(DMDestroy(&da)); 90 PetscCall(PetscFinalize()); 91 return 0; 92 } 93 94 /*TEST 95 96 test: 97 nsize: 8 98 args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic 99 100 test: 101 suffix: 2 102 nsize: 8 103 args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic -grid2d 104 output_file: output/ex7_1.out 105 106 test: 107 suffix: 3 108 nsize: 8 109 args: -dof 3 -stencil_width 2 -M 50 -N 50 -periodic -grid3d 110 output_file: output/ex7_1.out 111 112 TEST*/ 113