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