1 2 static char help[] = "Tests VecView() contour plotting for 2d DMDAs.\n\n"; 3 4 /* 5 MATLAB must be installed to configure PETSc to have MATLAB engine. 6 Unless you have specific important reasons for using the MATLAB engine, we do not 7 recommend it. If you want to use MATLAB for visualization and maybe a little post processing 8 then you can use the socket viewer and send the data to MATLAB via that. 9 10 VecView() on DMDA vectors first puts the Vec elements into global natural ordering before printing (or plotting) 11 them. In 2d 5 by 2 DMDA this means the numbering is 12 13 5 6 7 8 9 14 0 1 2 3 4 15 16 Now the default split across 2 processors with the DM is (by rank) 17 18 0 0 0 1 1 19 0 0 0 1 1 20 21 So the global PETSc ordering is 22 23 3 4 5 8 9 24 0 1 2 6 7 25 26 Use the options 27 -da_grid_x <nx> - number of grid points in x direction, if M < 0 28 -da_grid_y <ny> - number of grid points in y direction, if N < 0 29 -da_processors_x <MX> number of processors in x directio 30 -da_processors_y <MY> number of processors in x direction 31 */ 32 33 #include <petscdm.h> 34 #include <petscdmda.h> 35 36 int main(int argc, char **argv) 37 { 38 PetscMPIInt rank; 39 PetscInt M = 10, N = 8; 40 PetscBool flg = PETSC_FALSE; 41 DM da; 42 PetscViewer viewer; 43 Vec local, global; 44 PetscScalar value; 45 DMBoundaryType bx = DM_BOUNDARY_NONE, by = DM_BOUNDARY_NONE; 46 DMDAStencilType stype = DMDA_STENCIL_BOX; 47 #if defined(PETSC_HAVE_MATLAB) 48 PetscViewer mviewer; 49 PetscMPIInt size; 50 #endif 51 52 PetscFunctionBeginUser; 53 PetscCall(PetscInitialize(&argc, &argv, (char *)0, help)); 54 PetscCall(PetscViewerDrawOpen(PETSC_COMM_WORLD, 0, "", 300, 0, 300, 300, &viewer)); 55 #if defined(PETSC_HAVE_MATLAB) 56 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size)); 57 if (size == 1) PetscCall(PetscViewerMatlabOpen(PETSC_COMM_WORLD, "tmp.mat", FILE_MODE_WRITE, &mviewer)); 58 #endif 59 60 PetscCall(PetscOptionsGetBool(NULL, NULL, "-star_stencil", &flg, NULL)); 61 if (flg) stype = DMDA_STENCIL_STAR; 62 63 /* Create distributed array and get vectors */ 64 PetscCall(DMDACreate2d(PETSC_COMM_WORLD, bx, by, stype, M, N, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL, &da)); 65 PetscCall(DMSetFromOptions(da)); 66 PetscCall(DMSetUp(da)); 67 PetscCall(DMCreateGlobalVector(da, &global)); 68 PetscCall(DMCreateLocalVector(da, &local)); 69 70 value = -3.0; 71 PetscCall(VecSet(global, value)); 72 PetscCall(DMGlobalToLocalBegin(da, global, INSERT_VALUES, local)); 73 PetscCall(DMGlobalToLocalEnd(da, global, INSERT_VALUES, local)); 74 75 PetscCallMPI(MPI_Comm_rank(PETSC_COMM_WORLD, &rank)); 76 value = rank + 1; 77 PetscCall(VecScale(local, value)); 78 PetscCall(DMLocalToGlobalBegin(da, local, ADD_VALUES, global)); 79 PetscCall(DMLocalToGlobalEnd(da, local, ADD_VALUES, global)); 80 81 flg = PETSC_FALSE; 82 PetscCall(PetscOptionsGetBool(NULL, NULL, "-view_global", &flg, NULL)); 83 if (flg) { /* view global vector in natural ordering */ 84 PetscCall(VecView(global, PETSC_VIEWER_STDOUT_WORLD)); 85 } 86 PetscCall(DMView(da, viewer)); 87 PetscCall(VecView(global, viewer)); 88 #if defined(PETSC_HAVE_MATLAB) 89 if (size == 1) { 90 PetscCall(DMView(da, mviewer)); 91 PetscCall(VecView(global, mviewer)); 92 } 93 #endif 94 95 /* Free memory */ 96 #if defined(PETSC_HAVE_MATLAB) 97 if (size == 1) PetscCall(PetscViewerDestroy(&mviewer)); 98 #endif 99 PetscCall(PetscViewerDestroy(&viewer)); 100 PetscCall(VecDestroy(&local)); 101 PetscCall(VecDestroy(&global)); 102 PetscCall(DMDestroy(&da)); 103 PetscCall(PetscFinalize()); 104 return 0; 105 } 106 107 /*TEST 108 109 test: 110 requires: x 111 nsize: 2 112 args: -nox 113 114 TEST*/ 115