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