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