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 #include <petscviewerhdf5.h> 7 8 int main(int argc,char **argv) 9 { 10 PetscMPIInt rank,size; 11 PetscInt N = 6,M=8,P=5,dof=1; 12 PetscInt stencil_width=1,pt=0,st=0; 13 PetscBool flg2,flg3,isbinary,mpiio; 14 DMBoundaryType bx = DM_BOUNDARY_NONE,by = DM_BOUNDARY_NONE,bz = DM_BOUNDARY_NONE; 15 DMDAStencilType stencil_type = DMDA_STENCIL_STAR; 16 DM da,da2; 17 Vec global1,global2; 18 PetscScalar mone = -1.0; 19 PetscReal norm; 20 PetscViewer viewer; 21 PetscRandom rdm; 22 #if defined(PETSC_HAVE_HDF5) 23 PetscBool ishdf5; 24 #endif 25 26 PetscCall(PetscInitialize(&argc,&argv,(char*)0,help)); 27 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD,&rank)); 28 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD,&size)); 29 30 PetscCall(PetscOptionsGetInt(NULL,NULL,"-M",&M,NULL)); 31 PetscCall(PetscOptionsGetInt(NULL,NULL,"-N",&N,NULL)); 32 PetscCall(PetscOptionsGetInt(NULL,NULL,"-P",&P,NULL)); 33 PetscCall(PetscOptionsGetInt(NULL,NULL,"-dof",&dof,NULL)); 34 PetscCall(PetscOptionsGetInt(NULL,NULL,"-stencil_width",&stencil_width,NULL)); 35 PetscCall(PetscOptionsGetInt(NULL,NULL,"-periodic",&pt,NULL)); 36 if (pt == 1) bx = DM_BOUNDARY_PERIODIC; 37 if (pt == 2) by = DM_BOUNDARY_PERIODIC; 38 if (pt == 4) {bx = DM_BOUNDARY_PERIODIC; by = DM_BOUNDARY_PERIODIC;} 39 40 PetscCall(PetscOptionsGetInt(NULL,NULL,"-stencil_type",&st,NULL)); 41 stencil_type = (DMDAStencilType) st; 42 43 PetscCall(PetscOptionsHasName(NULL,NULL,"-oned",&flg2)); 44 PetscCall(PetscOptionsHasName(NULL,NULL,"-twod",&flg2)); 45 PetscCall(PetscOptionsHasName(NULL,NULL,"-threed",&flg3)); 46 47 PetscCall(PetscOptionsHasName(NULL,NULL,"-binary",&isbinary)); 48 #if defined(PETSC_HAVE_HDF5) 49 PetscCall(PetscOptionsHasName(NULL,NULL,"-hdf5",&ishdf5)); 50 #endif 51 PetscCall(PetscOptionsHasName(NULL,NULL,"-mpiio",&mpiio)); 52 if (flg2) { 53 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,&da)); 54 } else if (flg3) { 55 PetscCall(DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,stencil_type,M,N,P,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,0,&da)); 56 } else { 57 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da)); 58 } 59 PetscCall(DMSetFromOptions(da)); 60 PetscCall(DMSetUp(da)); 61 62 PetscCall(DMCreateGlobalVector(da,&global1)); 63 PetscCall(PetscObjectSetName((PetscObject)global1,"Test_Vec")); 64 PetscCall(PetscRandomCreate(PETSC_COMM_WORLD,&rdm)); 65 PetscCall(PetscRandomSetFromOptions(rdm)); 66 PetscCall(VecSetRandom(global1,rdm)); 67 if (isbinary) { 68 if (mpiio) { 69 PetscCall(PetscOptionsSetValue(NULL,"-viewer_binary_mpiio","")); 70 } 71 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer)); 72 #if defined(PETSC_HAVE_HDF5) 73 } else if (ishdf5) { 74 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD,"temp",FILE_MODE_WRITE,&viewer)); 75 #endif 76 } else SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid Viewer : Run with -binary or -hdf5 option"); 77 PetscCall(VecView(global1,viewer)); 78 PetscCall(PetscViewerDestroy(&viewer)); 79 80 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Global vector written to temp file is \n")); 81 PetscCall(VecView(global1,PETSC_VIEWER_STDOUT_WORLD)); 82 83 if (flg2) { 84 PetscCall(DMDACreate2d(PETSC_COMM_WORLD,bx,by,stencil_type,M,N,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,&da2)); 85 } else if (flg3) { 86 PetscCall(DMDACreate3d(PETSC_COMM_WORLD,bx,by,bz,stencil_type,M,N,P,PETSC_DECIDE,PETSC_DECIDE,PETSC_DECIDE,dof,stencil_width,0,0,0,&da2)); 87 } else { 88 PetscCall(DMDACreate1d(PETSC_COMM_WORLD,bx,M,dof,stencil_width,NULL,&da2)); 89 } 90 PetscCall(DMSetFromOptions(da2)); 91 PetscCall(DMSetUp(da2)); 92 93 if (isbinary) { 94 PetscCall(PetscViewerBinaryOpen(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer)); 95 #if defined(PETSC_HAVE_HDF5) 96 } else if (ishdf5) { 97 PetscCall(PetscViewerHDF5Open(PETSC_COMM_WORLD,"temp",FILE_MODE_READ,&viewer)); 98 #endif 99 } else SETERRQ(PETSC_COMM_WORLD,PETSC_ERR_SUP,"Invalid Viewer : Run with -binary or -hdf5 option"); 100 101 PetscCall(DMCreateGlobalVector(da2,&global2)); 102 PetscCall(PetscObjectSetName((PetscObject)global2,"Test_Vec")); 103 PetscCall(VecLoad(global2,viewer)); 104 PetscCall(PetscViewerDestroy(&viewer)); 105 106 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"Global vector read from temp file is \n")); 107 PetscCall(VecView(global2,PETSC_VIEWER_STDOUT_WORLD)); 108 PetscCall(VecAXPY(global2,mone,global1)); 109 PetscCall(VecNorm(global2,NORM_MAX,&norm)); 110 if (norm != 0.0) { 111 PetscCall(PetscPrintf(PETSC_COMM_WORLD,"ex23: Norm of difference %g should be zero\n",(double)norm)); 112 PetscCall(PetscPrintf(PETSC_COMM_WORLD," Number of processors %d\n",size)); 113 PetscCall(PetscPrintf(PETSC_COMM_WORLD," M,N,P,dof %D %D %D %D\n",M,N,P,dof)); 114 PetscCall(PetscPrintf(PETSC_COMM_WORLD," stencil_width %D stencil_type %d periodic %d\n",stencil_width,(int)stencil_type,(int)pt)); 115 PetscCall(PetscPrintf(PETSC_COMM_WORLD," dimension %d\n",1 + (int) flg2 + (int) flg3)); 116 } 117 118 PetscCall(PetscRandomDestroy(&rdm)); 119 PetscCall(DMDestroy(&da)); 120 PetscCall(DMDestroy(&da2)); 121 PetscCall(VecDestroy(&global1)); 122 PetscCall(VecDestroy(&global2)); 123 PetscCall(PetscFinalize()); 124 return 0; 125 } 126