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