1 2 static char help[] = "Tests VecView()/VecLoad() for DMDA vectors (this tests DMDAGlobalToNatural()).\n\n"; 3 4 #include <petscdm.h> 5 #include <petscdmda.h> 6 7 int main(int argc,char **argv) 8 { 9 PetscMPIInt size; 10 PetscInt N = 6,m=PETSC_DECIDE,n=PETSC_DECIDE,p=PETSC_DECIDE,M=8,dof=1,stencil_width=1,P=5,pt = 0,st = 0; 11 PetscBool flg2,flg3,native = PETSC_FALSE; 12 DMBoundaryType bx = DM_BOUNDARY_NONE,by = DM_BOUNDARY_NONE,bz = DM_BOUNDARY_NONE; 13 DMDAStencilType stencil_type = DMDA_STENCIL_STAR; 14 DM da; 15 Vec global1,global2,global3,global4; 16 PetscScalar mone = -1.0; 17 PetscReal norm; 18 PetscViewer viewer; 19 PetscRandom rdm; 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,"-P",&P,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 PetscCall(PetscOptionsGetBool(NULL,NULL,"-native",&native,NULL)); 29 if (pt == 1) bx = DM_BOUNDARY_PERIODIC; 30 if (pt == 2) by = DM_BOUNDARY_PERIODIC; 31 if (pt == 3) {bx = DM_BOUNDARY_PERIODIC; by = DM_BOUNDARY_PERIODIC;} 32 if (pt == 4) bz = DM_BOUNDARY_PERIODIC; 33 34 PetscCall(PetscOptionsGetInt(NULL,NULL,"-stencil_type",&st,NULL)); 35 stencil_type = (DMDAStencilType) st; 36 37 PetscCall(PetscOptionsHasName(NULL,NULL,"-one",&flg2)); 38 PetscCall(PetscOptionsHasName(NULL,NULL,"-two",&flg2)); 39 PetscCall(PetscOptionsHasName(NULL,NULL,"-three",&flg3)); 40 if (flg2) { 41 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,m,n,dof,stencil_width,0,0,&da)); 42 } else if (flg3) { 43 PetscCall(DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,stencil_type,M,N,P,m,n,p,dof,stencil_width,0,0,0,&da)); 44 } else { 45 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da)); 46 } 47 PetscCall(DMSetFromOptions(da)); 48 PetscCall(DMSetUp(da)); 49 50 PetscCall(DMCreateGlobalVector(da,&global1)); 51 PetscCall(PetscRandomCreate(PETSC_COMM_WORLD,&rdm)); 52 PetscCall(PetscRandomSetFromOptions(rdm)); 53 PetscCall(DMCreateGlobalVector(da,&global2)); 54 PetscCall(DMCreateGlobalVector(da,&global3)); 55 PetscCall(DMCreateGlobalVector(da,&global4)); 56 57 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer)); 58 if (native) PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE)); 59 PetscCall(VecSetRandom(global1,rdm)); 60 PetscCall(VecView(global1,viewer)); 61 PetscCall(VecSetRandom(global3,rdm)); 62 PetscCall(VecView(global3,viewer)); 63 if (native) PetscCall(PetscViewerPopFormat(viewer)); 64 PetscCall(PetscViewerDestroy(&viewer)); 65 66 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer)); 67 if (native) PetscCall(PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE)); 68 PetscCall(VecLoad(global2,viewer)); 69 PetscCall(VecLoad(global4,viewer)); 70 if (native) PetscCall(PetscViewerPopFormat(viewer)); 71 PetscCall(PetscViewerDestroy(&viewer)); 72 73 if (native) { 74 Vec filenative; 75 PetscBool same; 76 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer)); 77 PetscCall(DMDACreateNaturalVector(da,&filenative)); 78 /* DMDA "natural" Vec does not commandeer VecLoad. The following load will only work when run on the same process 79 * layout, where as the standard VecView/VecLoad (using DMDA and not PETSC_VIEWER_NATIVE) can be read on a different 80 * number of processors. */ 81 PetscCall(VecLoad(filenative,viewer)); 82 PetscCall(VecEqual(global2,filenative,&same)); 83 if (!same) { 84 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ex23: global vector does not match contents of file\n")); 85 PetscCall(VecView(global2,0)); 86 PetscCall(VecView(filenative,0)); 87 } 88 PetscCall(PetscViewerDestroy(&viewer)); 89 PetscCall(VecDestroy(&filenative)); 90 } 91 92 PetscCall(VecAXPY(global2,mone,global1)); 93 PetscCall(VecNorm(global2,NORM_MAX,&norm)); 94 if (norm != 0.0) { 95 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 96 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ex23: Norm of difference %g should be zero\n",(double)norm)); 97 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of processors %d\n",size)); 98 PetscCall(PetscPrintf(PETSC_COMM_WORLD," M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",M,N,P,dof)); 99 PetscCall(PetscPrintf(PETSC_COMM_WORLD," stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n",stencil_width,(int)stencil_type,(int)bx)); 100 PetscCall(PetscPrintf(PETSC_COMM_WORLD," dimension %d\n",1 + (int) flg2 + (int) flg3)); 101 } 102 PetscCall(VecAXPY(global4,mone,global3)); 103 PetscCall(VecNorm(global4,NORM_MAX,&norm)); 104 if (norm != 0.0) { 105 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 106 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ex23: Norm of difference %g should be zero\n",(double)norm)); 107 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of processors %d\n",size)); 108 PetscCall(PetscPrintf(PETSC_COMM_WORLD," M,N,P,dof %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT " %" PetscInt_FMT "\n",M,N,P,dof)); 109 PetscCall(PetscPrintf(PETSC_COMM_WORLD," stencil_width %" PetscInt_FMT " stencil_type %d periodic %d\n",stencil_width,(int)stencil_type,(int)bx)); 110 PetscCall(PetscPrintf(PETSC_COMM_WORLD," dimension %d\n",1 + (int) flg2 + (int) flg3)); 111 } 112 113 PetscCall(PetscRandomDestroy(&rdm)); 114 PetscCall(DMDestroy(&da)); 115 PetscCall(VecDestroy(&global1)); 116 PetscCall(VecDestroy(&global2)); 117 PetscCall(VecDestroy(&global3)); 118 PetscCall(VecDestroy(&global4)); 119 PetscCall(PetscFinalize()); 120 return 0; 121 } 122 123 /*TEST 124 125 test: 126 nsize: {{1 3}} 127 args: -one -dof {{1 2 3}} -stencil_type {{0 1}} 128 129 test: 130 suffix: 3 131 nsize: {{2 4}} 132 args: -two -dof {{1 3}} -stencil_type {{0 1}} 133 134 test: 135 suffix: 4 136 nsize: {{1 4}} 137 args: -three -dof {{2 3}} -stencil_type {{0 1}} 138 139 test: 140 suffix: 2 141 nsize: 2 142 args: -two -native 143 144 TEST*/ 145