Actual source code: ex23.c
2: static char help[] = "Tests VecView()/VecLoad() for DMDA vectors (this tests DMDAGlobalToNatural()).\n\n";
4: #include <petscdm.h>
5: #include <petscdmda.h>
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;
21: PetscInitialize(&argc,&argv,(char*)0,help);
22: PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL);
23: PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL);
24: PetscOptionsGetInt(NULL,NULL,"-P",&P,NULL);
25: PetscOptionsGetInt(NULL,NULL,"-dof",&dof,NULL);
26: PetscOptionsGetInt(NULL,NULL,"-stencil_width",&stencil_width,NULL);
27: PetscOptionsGetInt(NULL,NULL,"-periodic",&pt,NULL);
28: 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;
34: PetscOptionsGetInt(NULL,NULL,"-stencil_type",&st,NULL);
35: stencil_type = (DMDAStencilType) st;
37: PetscOptionsHasName(NULL,NULL,"-one",&flg2);
38: PetscOptionsHasName(NULL,NULL,"-two",&flg2);
39: PetscOptionsHasName(NULL,NULL,"-three",&flg3);
40: if (flg2) {
41: DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,m,n,dof,stencil_width,0,0,&da);
42: } else if (flg3) {
43: 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: DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da);
46: }
47: DMSetFromOptions(da);
48: DMSetUp(da);
50: DMCreateGlobalVector(da,&global1);
51: PetscRandomCreate(PETSC_COMM_WORLD,&rdm);
52: PetscRandomSetFromOptions(rdm);
53: DMCreateGlobalVector(da,&global2);
54: DMCreateGlobalVector(da,&global3);
55: DMCreateGlobalVector(da,&global4);
57: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer);
58: if (native) PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE);
59: VecSetRandom(global1,rdm);
60: VecView(global1,viewer);
61: VecSetRandom(global3,rdm);
62: VecView(global3,viewer);
63: if (native) PetscViewerPopFormat(viewer);
64: PetscViewerDestroy(&viewer);
66: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer);
67: if (native) PetscViewerPushFormat(viewer,PETSC_VIEWER_NATIVE);
68: VecLoad(global2,viewer);
69: VecLoad(global4,viewer);
70: if (native) PetscViewerPopFormat(viewer);
71: PetscViewerDestroy(&viewer);
73: if (native) {
74: Vec filenative;
75: PetscBool same;
76: PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer);
77: 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: VecLoad(filenative,viewer);
82: VecEqual(global2,filenative,&same);
83: if (!same) {
84: PetscPrintf(PETSC_COMM_WORLD,"ex23: global vector does not match contents of file\n");
85: VecView(global2,0);
86: VecView(filenative,0);
87: }
88: PetscViewerDestroy(&viewer);
89: VecDestroy(&filenative);
90: }
92: VecAXPY(global2,mone,global1);
93: VecNorm(global2,NORM_MAX,&norm);
94: if (norm != 0.0) {
95: MPI_Comm_size(PETSC_COMM_WORLD,&size);
96: PetscPrintf(PETSC_COMM_WORLD,"ex23: Norm of difference %g should be zero\n",(double)norm);
97: PetscPrintf(PETSC_COMM_WORLD," Number of processors %d\n",size);
98: PetscPrintf(PETSC_COMM_WORLD," M,N,P,dof %D %D %D %D\n",M,N,P,dof);
99: PetscPrintf(PETSC_COMM_WORLD," stencil_width %D stencil_type %d periodic %d\n",stencil_width,(int)stencil_type,(int)bx);
100: PetscPrintf(PETSC_COMM_WORLD," dimension %d\n",1 + (int) flg2 + (int) flg3);
101: }
102: VecAXPY(global4,mone,global3);
103: VecNorm(global4,NORM_MAX,&norm);
104: if (norm != 0.0) {
105: MPI_Comm_size(PETSC_COMM_WORLD,&size);
106: PetscPrintf(PETSC_COMM_WORLD,"ex23: Norm of difference %g should be zero\n",(double)norm);
107: PetscPrintf(PETSC_COMM_WORLD," Number of processors %d\n",size);
108: PetscPrintf(PETSC_COMM_WORLD," M,N,P,dof %D %D %D %D\n",M,N,P,dof);
109: PetscPrintf(PETSC_COMM_WORLD," stencil_width %D stencil_type %d periodic %d\n",stencil_width,(int)stencil_type,(int)bx);
110: PetscPrintf(PETSC_COMM_WORLD," dimension %d\n",1 + (int) flg2 + (int) flg3);
111: }
113: PetscRandomDestroy(&rdm);
114: DMDestroy(&da);
115: VecDestroy(&global1);
116: VecDestroy(&global2);
117: VecDestroy(&global3);
118: VecDestroy(&global4);
119: PetscFinalize();
120: return 0;
121: }
123: /*TEST
125: test:
126: nsize: {{1 3}}
127: args: -one -dof {{1 2 3}} -stencil_type {{0 1}}
129: test:
130: suffix: 3
131: nsize: {{2 4}}
132: args: -two -dof {{1 3}} -stencil_type {{0 1}}
134: test:
135: suffix: 4
136: nsize: {{1 4}}
137: args: -three -dof {{2 3}} -stencil_type {{0 1}}
139: test:
140: suffix: 2
141: nsize: 2
142: args: -two -native
144: TEST*/