xref: /petsc/src/snes/tutorials/ex19.c (revision e08b1d6d0faae6eca507e20c9d3498f81719d047)
1 
2 static char help[] = "Nonlinear driven cavity with multigrid in 2d.\n \
3   \n\
4 The 2D driven cavity problem is solved in a velocity-vorticity formulation.\n\
5 The flow can be driven with the lid or with bouyancy or both:\n\
6   -lidvelocity &ltlid&gt, where &ltlid&gt = dimensionless velocity of lid\n\
7   -grashof &ltgr&gt, where &ltgr&gt = dimensionless temperature gradent\n\
8   -prandtl &ltpr&gt, where &ltpr&gt = dimensionless thermal/momentum diffusity ratio\n\
9  -contours : draw contour plots of solution\n\n";
10 /* in HTML, '&lt' = '<' and '&gt' = '>' */
11 
12 /*
13       See src/ksp/ksp/tutorials/ex45.c
14 */
15 
16 /*F-----------------------------------------------------------------------
17 
18     We thank David E. Keyes for contributing the driven cavity discretization within this example code.
19 
20     This problem is modeled by the partial differential equation system
21 
22 \begin{eqnarray}
23         - \triangle U - \nabla_y \Omega & = & 0  \\
24         - \triangle V + \nabla_x\Omega & = & 0  \\
25         - \triangle \Omega + \nabla \cdot ([U*\Omega,V*\Omega]) - GR* \nabla_x T & = & 0  \\
26         - \triangle T + PR* \nabla \cdot ([U*T,V*T]) & = & 0
27 \end{eqnarray}
28 
29     in the unit square, which is uniformly discretized in each of x and y in this simple encoding.
30 
31     No-slip, rigid-wall Dirichlet conditions are used for $ [U,V]$.
32     Dirichlet conditions are used for Omega, based on the definition of
33     vorticity: $ \Omega = - \nabla_y U + \nabla_x V$, where along each
34     constant coordinate boundary, the tangential derivative is zero.
35     Dirichlet conditions are used for T on the left and right walls,
36     and insulation homogeneous Neumann conditions are used for T on
37     the top and bottom walls.
38 
39     A finite difference approximation with the usual 5-point stencil
40     is used to discretize the boundary value problem to obtain a
41     nonlinear system of equations.  Upwinding is used for the divergence
42     (convective) terms and central for the gradient (source) terms.
43 
44     The Jacobian can be either
45       * formed via finite differencing using coloring (the default), or
46       * applied matrix-free via the option -snes_mf
47         (for larger grid problems this variant may not converge
48         without a preconditioner due to ill-conditioning).
49 
50   ------------------------------------------------------------------------F*/
51 
52 /*
53    Include "petscdmda.h" so that we can use distributed arrays (DMDAs).
54    Include "petscsnes.h" so that we can use SNES solvers.  Note that this
55    file automatically includes:
56      petscsys.h       - base PETSc routines   petscvec.h - vectors
57      petscmat.h - matrices
58      petscis.h     - index sets            petscksp.h - Krylov subspace methods
59      petscviewer.h - viewers               petscpc.h  - preconditioners
60      petscksp.h   - linear solvers
61 */
62 #if defined(PETSC_APPLE_FRAMEWORK)
63 #import <PETSc/petscsnes.h>
64 #import <PETSc/petscdmda.h>
65 #else
66 #include <petscsnes.h>
67 #include <petscdm.h>
68 #include <petscdmda.h>
69 #endif
70 
71 /*
72    User-defined routines and data structures
73 */
74 typedef struct {
75   PetscScalar u,v,omega,temp;
76 } Field;
77 
78 PetscErrorCode FormFunctionLocal(DMDALocalInfo*,Field**,Field**,void*);
79 
80 typedef struct {
81   PetscReal   lidvelocity,prandtl,grashof;  /* physical parameters */
82   PetscBool   draw_contours;                /* flag - 1 indicates drawing contours */
83 } AppCtx;
84 
85 extern PetscErrorCode FormInitialGuess(AppCtx*,DM,Vec);
86 extern PetscErrorCode NonlinearGS(SNES,Vec,Vec,void*);
87 
88 int main(int argc,char **argv)
89 {
90   AppCtx         user;                /* user-defined work context */
91   PetscInt       mx,my,its;
92   MPI_Comm       comm;
93   SNES           snes;
94   DM             da;
95   Vec            x;
96 
97   PetscFunctionBeginUser;
98   PetscCall(PetscInitialize(&argc,&argv,(char*)0,help));
99   comm = PETSC_COMM_WORLD;
100   PetscCall(SNESCreate(comm,&snes));
101 
102   /*
103       Create distributed array object to manage parallel grid and vectors
104       for principal unknowns (x) and governing residuals (f)
105   */
106   PetscCall(DMDACreate2d(PETSC_COMM_WORLD,DM_BOUNDARY_NONE,DM_BOUNDARY_NONE,DMDA_STENCIL_STAR,4,4,PETSC_DECIDE,PETSC_DECIDE,4,1,0,0,&da));
107   PetscCall(DMSetFromOptions(da));
108   PetscCall(DMSetUp(da));
109   PetscCall(SNESSetDM(snes,(DM)da));
110   PetscCall(SNESSetNGS(snes, NonlinearGS, (void*)&user));
111 
112   PetscCall(DMDAGetInfo(da,0,&mx,&my,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE,PETSC_IGNORE));
113   /*
114      Problem parameters (velocity of lid, prandtl, and grashof numbers)
115   */
116   user.lidvelocity = 1.0/(mx*my);
117   user.prandtl     = 1.0;
118   user.grashof     = 1.0;
119 
120   PetscCall(PetscOptionsGetReal(NULL,NULL,"-lidvelocity",&user.lidvelocity,NULL));
121   PetscCall(PetscOptionsGetReal(NULL,NULL,"-prandtl",&user.prandtl,NULL));
122   PetscCall(PetscOptionsGetReal(NULL,NULL,"-grashof",&user.grashof,NULL));
123   PetscCall(PetscOptionsHasName(NULL,NULL,"-contours",&user.draw_contours));
124 
125   PetscCall(DMDASetFieldName(da,0,"x_velocity"));
126   PetscCall(DMDASetFieldName(da,1,"y_velocity"));
127   PetscCall(DMDASetFieldName(da,2,"Omega"));
128   PetscCall(DMDASetFieldName(da,3,"temperature"));
129 
130   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
131      Create user context, set problem data, create vector data structures.
132      Also, compute the initial guess.
133      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
134 
135   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
136      Create nonlinear solver context
137 
138      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
139   PetscCall(DMSetApplicationContext(da,&user));
140   PetscCall(DMDASNESSetFunctionLocal(da,INSERT_VALUES,(PetscErrorCode (*)(DMDALocalInfo*,void*,void*,void*))FormFunctionLocal,&user));
141   PetscCall(SNESSetFromOptions(snes));
142   PetscCall(PetscPrintf(comm,"lid velocity = %g, prandtl # = %g, grashof # = %g\n",(double)user.lidvelocity,(double)user.prandtl,(double)user.grashof));
143 
144   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
145      Solve the nonlinear system
146      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
147   PetscCall(DMCreateGlobalVector(da,&x));
148   PetscCall(FormInitialGuess(&user,da,x));
149 
150   PetscCall(SNESSolve(snes,NULL,x));
151 
152   PetscCall(SNESGetIterationNumber(snes,&its));
153   PetscCall(PetscPrintf(comm,"Number of SNES iterations = %" PetscInt_FMT "\n", its));
154 
155   /*
156      Visualize solution
157   */
158   if (user.draw_contours) PetscCall(VecView(x,PETSC_VIEWER_DRAW_WORLD));
159 
160   /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
161      Free work space.  All PETSc objects should be destroyed when they
162      are no longer needed.
163      - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
164   PetscCall(VecDestroy(&x));
165   PetscCall(DMDestroy(&da));
166   PetscCall(SNESDestroy(&snes));
167   PetscCall(PetscFinalize());
168   return 0;
169 }
170 
171 /* ------------------------------------------------------------------- */
172 
173 /*
174    FormInitialGuess - Forms initial approximation.
175 
176    Input Parameters:
177    user - user-defined application context
178    X - vector
179 
180    Output Parameter:
181    X - vector
182 */
183 PetscErrorCode FormInitialGuess(AppCtx *user,DM da,Vec X)
184 {
185   PetscInt       i,j,mx,xs,ys,xm,ym;
186   PetscReal      grashof,dx;
187   Field          **x;
188 
189   PetscFunctionBeginUser;
190   grashof = user->grashof;
191 
192   PetscCall(DMDAGetInfo(da,0,&mx,0,0,0,0,0,0,0,0,0,0,0));
193   dx   = 1.0/(mx-1);
194 
195   /*
196      Get local grid boundaries (for 2-dimensional DMDA):
197        xs, ys   - starting grid indices (no ghost points)
198        xm, ym   - widths of local grid (no ghost points)
199   */
200   PetscCall(DMDAGetCorners(da,&xs,&ys,NULL,&xm,&ym,NULL));
201 
202   /*
203      Get a pointer to vector data.
204        - For default PETSc vectors, VecGetArray() returns a pointer to
205          the data array.  Otherwise, the routine is implementation dependent.
206        - You MUST call VecRestoreArray() when you no longer need access to
207          the array.
208   */
209   PetscCall(DMDAVecGetArrayWrite(da,X,&x));
210 
211   /*
212      Compute initial guess over the locally owned part of the grid
213      Initial condition is motionless fluid and equilibrium temperature
214   */
215   for (j=ys; j<ys+ym; j++) {
216     for (i=xs; i<xs+xm; i++) {
217       x[j][i].u     = 0.0;
218       x[j][i].v     = 0.0;
219       x[j][i].omega = 0.0;
220       x[j][i].temp  = (grashof>0)*i*dx;
221     }
222   }
223 
224   /*
225      Restore vector
226   */
227   PetscCall(DMDAVecRestoreArrayWrite(da,X,&x));
228   PetscFunctionReturn(0);
229 }
230 
231 PetscErrorCode FormFunctionLocal(DMDALocalInfo *info,Field **x,Field **f,void *ptr)
232 {
233   AppCtx         *user = (AppCtx*)ptr;
234   PetscInt       xints,xinte,yints,yinte,i,j;
235   PetscReal      hx,hy,dhx,dhy,hxdhy,hydhx;
236   PetscReal      grashof,prandtl,lid;
237   PetscScalar    u,uxx,uyy,vx,vy,avx,avy,vxp,vxm,vyp,vym;
238 
239   PetscFunctionBeginUser;
240   grashof = user->grashof;
241   prandtl = user->prandtl;
242   lid     = user->lidvelocity;
243 
244   /*
245      Define mesh intervals ratios for uniform grid.
246 
247      Note: FD formulae below are normalized by multiplying through by
248      local volume element (i.e. hx*hy) to obtain coefficients O(1) in two dimensions.
249 
250   */
251   dhx   = (PetscReal)(info->mx-1);  dhy = (PetscReal)(info->my-1);
252   hx    = 1.0/dhx;                   hy = 1.0/dhy;
253   hxdhy = hx*dhy;                 hydhx = hy*dhx;
254 
255   xints = info->xs; xinte = info->xs+info->xm; yints = info->ys; yinte = info->ys+info->ym;
256 
257   /* Test whether we are on the bottom edge of the global array */
258   if (yints == 0) {
259     j     = 0;
260     yints = yints + 1;
261     /* bottom edge */
262     for (i=info->xs; i<info->xs+info->xm; i++) {
263       f[j][i].u     = x[j][i].u;
264       f[j][i].v     = x[j][i].v;
265       f[j][i].omega = x[j][i].omega + (x[j+1][i].u - x[j][i].u)*dhy;
266       f[j][i].temp  = x[j][i].temp-x[j+1][i].temp;
267     }
268   }
269 
270   /* Test whether we are on the top edge of the global array */
271   if (yinte == info->my) {
272     j     = info->my - 1;
273     yinte = yinte - 1;
274     /* top edge */
275     for (i=info->xs; i<info->xs+info->xm; i++) {
276       f[j][i].u     = x[j][i].u - lid;
277       f[j][i].v     = x[j][i].v;
278       f[j][i].omega = x[j][i].omega + (x[j][i].u - x[j-1][i].u)*dhy;
279       f[j][i].temp  = x[j][i].temp-x[j-1][i].temp;
280     }
281   }
282 
283   /* Test whether we are on the left edge of the global array */
284   if (xints == 0) {
285     i     = 0;
286     xints = xints + 1;
287     /* left edge */
288     for (j=info->ys; j<info->ys+info->ym; j++) {
289       f[j][i].u     = x[j][i].u;
290       f[j][i].v     = x[j][i].v;
291       f[j][i].omega = x[j][i].omega - (x[j][i+1].v - x[j][i].v)*dhx;
292       f[j][i].temp  = x[j][i].temp;
293     }
294   }
295 
296   /* Test whether we are on the right edge of the global array */
297   if (xinte == info->mx) {
298     i     = info->mx - 1;
299     xinte = xinte - 1;
300     /* right edge */
301     for (j=info->ys; j<info->ys+info->ym; j++) {
302       f[j][i].u     = x[j][i].u;
303       f[j][i].v     = x[j][i].v;
304       f[j][i].omega = x[j][i].omega - (x[j][i].v - x[j][i-1].v)*dhx;
305       f[j][i].temp  = x[j][i].temp - (PetscReal)(grashof>0);
306     }
307   }
308 
309   /* Compute over the interior points */
310   for (j=yints; j<yinte; j++) {
311     for (i=xints; i<xinte; i++) {
312 
313       /*
314        convective coefficients for upwinding
315       */
316       vx  = x[j][i].u; avx = PetscAbsScalar(vx);
317       vxp = .5*(vx+avx); vxm = .5*(vx-avx);
318       vy  = x[j][i].v; avy = PetscAbsScalar(vy);
319       vyp = .5*(vy+avy); vym = .5*(vy-avy);
320 
321       /* U velocity */
322       u         = x[j][i].u;
323       uxx       = (2.0*u - x[j][i-1].u - x[j][i+1].u)*hydhx;
324       uyy       = (2.0*u - x[j-1][i].u - x[j+1][i].u)*hxdhy;
325       f[j][i].u = uxx + uyy - .5*(x[j+1][i].omega-x[j-1][i].omega)*hx;
326 
327       /* V velocity */
328       u         = x[j][i].v;
329       uxx       = (2.0*u - x[j][i-1].v - x[j][i+1].v)*hydhx;
330       uyy       = (2.0*u - x[j-1][i].v - x[j+1][i].v)*hxdhy;
331       f[j][i].v = uxx + uyy + .5*(x[j][i+1].omega-x[j][i-1].omega)*hy;
332 
333       /* Omega */
334       u             = x[j][i].omega;
335       uxx           = (2.0*u - x[j][i-1].omega - x[j][i+1].omega)*hydhx;
336       uyy           = (2.0*u - x[j-1][i].omega - x[j+1][i].omega)*hxdhy;
337       f[j][i].omega = uxx + uyy + (vxp*(u - x[j][i-1].omega) + vxm*(x[j][i+1].omega - u))*hy +
338                       (vyp*(u - x[j-1][i].omega) + vym*(x[j+1][i].omega - u))*hx -
339                       .5*grashof*(x[j][i+1].temp - x[j][i-1].temp)*hy;
340 
341       /* Temperature */
342       u            = x[j][i].temp;
343       uxx          = (2.0*u - x[j][i-1].temp - x[j][i+1].temp)*hydhx;
344       uyy          = (2.0*u - x[j-1][i].temp - x[j+1][i].temp)*hxdhy;
345       f[j][i].temp =  uxx + uyy  + prandtl*((vxp*(u - x[j][i-1].temp) + vxm*(x[j][i+1].temp - u))*hy +
346                                             (vyp*(u - x[j-1][i].temp) + vym*(x[j+1][i].temp - u))*hx);
347     }
348   }
349 
350   /*
351      Flop count (multiply-adds are counted as 2 operations)
352   */
353   PetscCall(PetscLogFlops(84.0*info->ym*info->xm));
354   PetscFunctionReturn(0);
355 }
356 
357 /*
358     Performs sweeps of point block nonlinear Gauss-Seidel on all the local grid points
359 */
360 PetscErrorCode NonlinearGS(SNES snes, Vec X, Vec B, void *ctx)
361 {
362   DMDALocalInfo  info;
363   Field          **x,**b;
364   Vec            localX, localB;
365   DM             da;
366   PetscInt       xints,xinte,yints,yinte,i,j,k,l;
367   PetscInt       max_its,tot_its;
368   PetscInt       sweeps;
369   PetscReal      rtol,atol,stol;
370   PetscReal      hx,hy,dhx,dhy,hxdhy,hydhx;
371   PetscReal      grashof,prandtl,lid;
372   PetscScalar    u,uxx,uyy,vx,vy,avx,avy,vxp,vxm,vyp,vym;
373   PetscScalar    fu, fv, fomega, ftemp;
374   PetscScalar    dfudu;
375   PetscScalar    dfvdv;
376   PetscScalar    dfodu, dfodv, dfodo;
377   PetscScalar    dftdu, dftdv, dftdt;
378   PetscScalar    yu=0, yv=0, yo=0, yt=0;
379   PetscScalar    bjiu, bjiv, bjiomega, bjitemp;
380   PetscBool      ptconverged;
381   PetscReal      pfnorm,pfnorm0,pynorm,pxnorm;
382   AppCtx         *user = (AppCtx*)ctx;
383 
384   PetscFunctionBeginUser;
385   grashof = user->grashof;
386   prandtl = user->prandtl;
387   lid     = user->lidvelocity;
388   tot_its = 0;
389   PetscCall(SNESNGSGetTolerances(snes,&rtol,&atol,&stol,&max_its));
390   PetscCall(SNESNGSGetSweeps(snes,&sweeps));
391   PetscCall(SNESGetDM(snes,(DM*)&da));
392   PetscCall(DMGetLocalVector(da,&localX));
393   if (B) {
394     PetscCall(DMGetLocalVector(da,&localB));
395   }
396   /*
397      Scatter ghost points to local vector, using the 2-step process
398         DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
399   */
400   PetscCall(DMGlobalToLocalBegin(da,X,INSERT_VALUES,localX));
401   PetscCall(DMGlobalToLocalEnd(da,X,INSERT_VALUES,localX));
402   if (B) {
403     PetscCall(DMGlobalToLocalBegin(da,B,INSERT_VALUES,localB));
404     PetscCall(DMGlobalToLocalEnd(da,B,INSERT_VALUES,localB));
405   }
406   PetscCall(DMDAGetLocalInfo(da,&info));
407   PetscCall(DMDAVecGetArrayWrite(da,localX,&x));
408   if (B) {
409     PetscCall(DMDAVecGetArrayRead(da,localB,&b));
410   }
411   /* looks like a combination of the formfunction / formjacobian routines */
412   dhx   = (PetscReal)(info.mx-1);dhy   = (PetscReal)(info.my-1);
413   hx    = 1.0/dhx;               hy    = 1.0/dhy;
414   hxdhy = hx*dhy;                hydhx = hy*dhx;
415 
416   xints = info.xs; xinte = info.xs+info.xm; yints = info.ys; yinte = info.ys+info.ym;
417 
418   /* Set the boundary conditions on the momentum equations */
419   /* Test whether we are on the bottom edge of the global array */
420   if (yints == 0) {
421     j     = 0;
422     /* bottom edge */
423     for (i=info.xs; i<info.xs+info.xm; i++) {
424 
425       if (B) {
426         bjiu = b[j][i].u;
427         bjiv = b[j][i].v;
428       } else {
429         bjiu = 0.0;
430         bjiv = 0.0;
431       }
432       x[j][i].u = 0.0 + bjiu;
433       x[j][i].v = 0.0 + bjiv;
434     }
435   }
436 
437   /* Test whether we are on the top edge of the global array */
438   if (yinte == info.my) {
439     j     = info.my - 1;
440     /* top edge */
441     for (i=info.xs; i<info.xs+info.xm; i++) {
442       if (B) {
443         bjiu = b[j][i].u;
444         bjiv = b[j][i].v;
445       } else {
446         bjiu = 0.0;
447         bjiv = 0.0;
448       }
449       x[j][i].u = lid + bjiu;
450       x[j][i].v = bjiv;
451     }
452   }
453 
454   /* Test whether we are on the left edge of the global array */
455   if (xints == 0) {
456     i     = 0;
457     /* left edge */
458     for (j=info.ys; j<info.ys+info.ym; j++) {
459       if (B) {
460         bjiu = b[j][i].u;
461         bjiv = b[j][i].v;
462       } else {
463         bjiu = 0.0;
464         bjiv = 0.0;
465       }
466       x[j][i].u = 0.0 + bjiu;
467       x[j][i].v = 0.0 + bjiv;
468     }
469   }
470 
471   /* Test whether we are on the right edge of the global array */
472   if (xinte == info.mx) {
473     i     = info.mx - 1;
474     /* right edge */
475     for (j=info.ys; j<info.ys+info.ym; j++) {
476       if (B) {
477         bjiu = b[j][i].u;
478         bjiv = b[j][i].v;
479       } else {
480         bjiu = 0.0;
481         bjiv = 0.0;
482       }
483       x[j][i].u = 0.0 + bjiu;
484       x[j][i].v = 0.0 + bjiv;
485     }
486   }
487 
488   for (k=0; k < sweeps; k++) {
489     for (j=info.ys; j<info.ys + info.ym; j++) {
490       for (i=info.xs; i<info.xs + info.xm; i++) {
491         ptconverged = PETSC_FALSE;
492         pfnorm0     = 0.0;
493         fu          = 0.0;
494         fv          = 0.0;
495         fomega      = 0.0;
496         ftemp       = 0.0;
497         /*  Run Newton's method on a single grid point */
498         for (l = 0; l < max_its && !ptconverged; l++) {
499           if (B) {
500             bjiu     = b[j][i].u;
501             bjiv     = b[j][i].v;
502             bjiomega = b[j][i].omega;
503             bjitemp  = b[j][i].temp;
504           } else {
505             bjiu     = 0.0;
506             bjiv     = 0.0;
507             bjiomega = 0.0;
508             bjitemp  = 0.0;
509           }
510 
511           if (i != 0 && i != info.mx - 1 && j != 0 && j != info.my-1) {
512             /* U velocity */
513             u     = x[j][i].u;
514             uxx   = (2.0*u - x[j][i-1].u - x[j][i+1].u)*hydhx;
515             uyy   = (2.0*u - x[j-1][i].u - x[j+1][i].u)*hxdhy;
516             fu    = uxx + uyy - .5*(x[j+1][i].omega-x[j-1][i].omega)*hx - bjiu;
517             dfudu = 2.0*(hydhx + hxdhy);
518             /* V velocity */
519             u     = x[j][i].v;
520             uxx   = (2.0*u - x[j][i-1].v - x[j][i+1].v)*hydhx;
521             uyy   = (2.0*u - x[j-1][i].v - x[j+1][i].v)*hxdhy;
522             fv    = uxx + uyy + .5*(x[j][i+1].omega-x[j][i-1].omega)*hy - bjiv;
523             dfvdv = 2.0*(hydhx + hxdhy);
524             /*
525              convective coefficients for upwinding
526              */
527             vx  = x[j][i].u; avx = PetscAbsScalar(vx);
528             vxp = .5*(vx+avx); vxm = .5*(vx-avx);
529             vy  = x[j][i].v; avy = PetscAbsScalar(vy);
530             vyp = .5*(vy+avy); vym = .5*(vy-avy);
531             /* Omega */
532             u      = x[j][i].omega;
533             uxx    = (2.0*u - x[j][i-1].omega - x[j][i+1].omega)*hydhx;
534             uyy    = (2.0*u - x[j-1][i].omega - x[j+1][i].omega)*hxdhy;
535             fomega = uxx + uyy +  (vxp*(u - x[j][i-1].omega) + vxm*(x[j][i+1].omega - u))*hy +
536                      (vyp*(u - x[j-1][i].omega) + vym*(x[j+1][i].omega - u))*hx -
537                      .5*grashof*(x[j][i+1].temp - x[j][i-1].temp)*hy - bjiomega;
538             /* convective coefficient derivatives */
539             dfodo = 2.0*(hydhx + hxdhy) + ((vxp - vxm)*hy + (vyp - vym)*hx);
540             if (PetscRealPart(vx) > 0.0) dfodu = (u - x[j][i-1].omega)*hy;
541             else dfodu = (x[j][i+1].omega - u)*hy;
542 
543             if (PetscRealPart(vy) > 0.0) dfodv = (u - x[j-1][i].omega)*hx;
544             else dfodv = (x[j+1][i].omega - u)*hx;
545 
546             /* Temperature */
547             u     = x[j][i].temp;
548             uxx   = (2.0*u - x[j][i-1].temp - x[j][i+1].temp)*hydhx;
549             uyy   = (2.0*u - x[j-1][i].temp - x[j+1][i].temp)*hxdhy;
550             ftemp =  uxx + uyy  + prandtl*((vxp*(u - x[j][i-1].temp) + vxm*(x[j][i+1].temp - u))*hy + (vyp*(u - x[j-1][i].temp) + vym*(x[j+1][i].temp - u))*hx) - bjitemp;
551             dftdt = 2.0*(hydhx + hxdhy) + prandtl*((vxp - vxm)*hy + (vyp - vym)*hx);
552             if (PetscRealPart(vx) > 0.0) dftdu = prandtl*(u - x[j][i-1].temp)*hy;
553             else dftdu = prandtl*(x[j][i+1].temp - u)*hy;
554 
555             if (PetscRealPart(vy) > 0.0) dftdv = prandtl*(u - x[j-1][i].temp)*hx;
556             else dftdv = prandtl*(x[j+1][i].temp - u)*hx;
557 
558             /* invert the system:
559              [ dfu / du     0        0        0    ][yu] = [fu]
560              [     0    dfv / dv     0        0    ][yv]   [fv]
561              [ dfo / du dfo / dv dfo / do     0    ][yo]   [fo]
562              [ dft / du dft / dv     0    dft / dt ][yt]   [ft]
563              by simple back-substitution
564            */
565             yu = fu / dfudu;
566             yv = fv / dfvdv;
567             yo = (fomega - (dfodu*yu + dfodv*yv)) / dfodo;
568             yt = (ftemp - (dftdu*yu + dftdv*yv)) / dftdt;
569 
570             x[j][i].u     = x[j][i].u - yu;
571             x[j][i].v     = x[j][i].v - yv;
572             x[j][i].temp  = x[j][i].temp - yt;
573             x[j][i].omega = x[j][i].omega - yo;
574           }
575           if (i == 0) {
576             fomega        = x[j][i].omega - (x[j][i+1].v - x[j][i].v)*dhx - bjiomega;
577             ftemp         = x[j][i].temp - bjitemp;
578             yo            = fomega;
579             yt            = ftemp;
580             x[j][i].omega = x[j][i].omega - fomega;
581             x[j][i].temp  = x[j][i].temp - ftemp;
582           }
583           if (i == info.mx - 1) {
584             fomega        = x[j][i].omega - (x[j][i].v - x[j][i-1].v)*dhx - bjiomega;
585             ftemp         = x[j][i].temp - (PetscReal)(grashof>0) - bjitemp;
586             yo            = fomega;
587             yt            = ftemp;
588             x[j][i].omega = x[j][i].omega - fomega;
589             x[j][i].temp  = x[j][i].temp - ftemp;
590           }
591           if (j == 0) {
592             fomega        = x[j][i].omega + (x[j+1][i].u - x[j][i].u)*dhy - bjiomega;
593             ftemp         = x[j][i].temp-x[j+1][i].temp - bjitemp;
594             yo            = fomega;
595             yt            = ftemp;
596             x[j][i].omega = x[j][i].omega - fomega;
597             x[j][i].temp  = x[j][i].temp - ftemp;
598           }
599           if (j == info.my - 1) {
600             fomega        = x[j][i].omega + (x[j][i].u - x[j-1][i].u)*dhy - bjiomega;
601             ftemp         = x[j][i].temp-x[j-1][i].temp - bjitemp;
602             yo            = fomega;
603             yt            = ftemp;
604             x[j][i].omega = x[j][i].omega - fomega;
605             x[j][i].temp  = x[j][i].temp - ftemp;
606           }
607           tot_its++;
608           pfnorm = PetscRealPart(fu*fu + fv*fv + fomega*fomega + ftemp*ftemp);
609           pfnorm = PetscSqrtReal(pfnorm);
610           pynorm = PetscRealPart(yu*yu + yv*yv + yo*yo + yt*yt);
611           pynorm = PetscSqrtReal(pynorm);
612           pxnorm = PetscRealPart(x[j][i].u*x[j][i].u + x[j][i].v*x[j][i].v + x[j][i].omega*x[j][i].omega + x[j][i].temp*x[j][i].temp);
613           pxnorm = PetscSqrtReal(pxnorm);
614           if (l == 0) pfnorm0 = pfnorm;
615           if (rtol*pfnorm0 >pfnorm || atol > pfnorm || pxnorm*stol > pynorm) ptconverged = PETSC_TRUE;
616         }
617       }
618     }
619   }
620   PetscCall(DMDAVecRestoreArrayWrite(da,localX,&x));
621   if (B) {
622     PetscCall(DMDAVecRestoreArrayRead(da,localB,&b));
623   }
624   PetscCall(DMLocalToGlobalBegin(da,localX,INSERT_VALUES,X));
625   PetscCall(DMLocalToGlobalEnd(da,localX,INSERT_VALUES,X));
626   PetscCall(PetscLogFlops(tot_its*(84.0 + 41.0 + 26.0)));
627   PetscCall(DMRestoreLocalVector(da,&localX));
628   if (B) {
629     PetscCall(DMRestoreLocalVector(da,&localB));
630   }
631   PetscFunctionReturn(0);
632 }
633 
634 /*TEST
635 
636    test:
637       nsize: 2
638       args: -da_refine 3 -snes_monitor_short -pc_type mg -ksp_type fgmres -pc_mg_type full
639       requires: !single
640 
641    test:
642       suffix: 10
643       nsize: 3
644       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -pc_fieldsplit_type symmetric_multiplicative -snes_view -da_refine 1 -ksp_type fgmres
645       requires: !single
646 
647    test:
648       suffix: 11
649       nsize: 4
650       requires: pastix
651       args: -snes_monitor_short -pc_type redundant -dm_mat_type mpiaij -redundant_pc_factor_mat_solver_type pastix -pc_redundant_number 2 -da_refine 4 -ksp_type fgmres
652 
653    test:
654       suffix: 12
655       nsize: 12
656       requires: pastix
657       args: -snes_monitor_short -pc_type redundant -dm_mat_type mpiaij -redundant_pc_factor_mat_solver_type pastix -pc_redundant_number 5 -da_refine 4 -ksp_type fgmres
658 
659    test:
660       suffix: 13
661       nsize: 3
662       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -pc_fieldsplit_type multiplicative -snes_view -da_refine 1 -ksp_type fgmres -snes_mf_operator
663       requires: !single
664 
665    test:
666       suffix: 14
667       nsize: 4
668       args: -snes_monitor_short -pc_type mg -dm_mat_type baij -mg_coarse_pc_type bjacobi -da_refine 3 -ksp_type fgmres
669       requires: !single
670 
671    test:
672       suffix: 14_ds
673       nsize: 4
674       args: -snes_converged_reason -pc_type mg -dm_mat_type baij -mg_coarse_pc_type bjacobi -da_refine 3 -ksp_type fgmres -mat_fd_type ds
675       output_file: output/ex19_2.out
676       requires: !single
677 
678    test:
679       suffix: 17
680       args: -snes_monitor_short -ksp_pc_side right
681       requires: !single
682 
683    test:
684       suffix: 18
685       args: -snes_monitor_ksp draw::draw_lg -ksp_pc_side right
686       requires: x !single
687 
688    test:
689       suffix: 19
690       nsize: 2
691       args: -da_refine 3 -snes_monitor_short -pc_type mg -ksp_type fgmres -pc_mg_type full -snes_type newtontrdc
692       requires: !single
693 
694    test:
695       suffix: 20
696       nsize: 2
697       args: -da_refine 3 -snes_monitor_short -pc_type mg -ksp_type fgmres -pc_mg_type full -snes_type newtontrdc -snes_trdc_use_cauchy false
698       requires: !single
699 
700    test:
701       suffix: 2
702       nsize: 4
703       args: -da_refine 3 -snes_converged_reason -pc_type mg -mat_fd_type ds
704       requires: !single
705 
706    test:
707       suffix: 2_bcols1
708       nsize: 4
709       args: -da_refine 3 -snes_converged_reason -pc_type mg -mat_fd_type ds -mat_fd_coloring_bcols
710       output_file: output/ex19_2.out
711       requires: !single
712 
713    test:
714       suffix: 3
715       nsize: 4
716       requires: mumps
717       args: -da_refine 3 -snes_monitor_short -pc_type redundant -dm_mat_type mpiaij -redundant_ksp_type preonly -redundant_pc_factor_mat_solver_type mumps -pc_redundant_number 2
718 
719    test:
720       suffix: 4
721       nsize: 12
722       requires: mumps
723       args: -da_refine 3 -snes_monitor_short -pc_type redundant -dm_mat_type mpiaij -redundant_ksp_type preonly -redundant_pc_factor_mat_solver_type mumps -pc_redundant_number 5
724       output_file: output/ex19_3.out
725 
726    test:
727       suffix: 6
728       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -snes_view -ksp_type fgmres -da_refine 1
729       requires: !single
730 
731    test:
732       suffix: 7
733       nsize: 3
734       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -snes_view -da_refine 1 -ksp_type fgmres
735 
736       requires: !single
737    test:
738       suffix: 8
739       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -pc_fieldsplit_block_size 2 -pc_fieldsplit_0_fields 0,1 -pc_fieldsplit_1_fields 0,1 -pc_fieldsplit_type multiplicative -snes_view -fieldsplit_pc_type lu -da_refine 1 -ksp_type fgmres
740       requires: !single
741 
742    test:
743       suffix: 9
744       nsize: 3
745       args: -snes_monitor_short -ksp_monitor_short -pc_type fieldsplit -pc_fieldsplit_type multiplicative -snes_view -da_refine 1 -ksp_type fgmres
746       requires: !single
747 
748    test:
749       suffix: aspin
750       nsize: 4
751       args: -da_refine 3 -da_overlap 2 -snes_monitor_short -snes_type aspin -grashof 4e4 -lidvelocity 100 -ksp_monitor_short
752       requires: !single
753 
754    test:
755       suffix: bcgsl
756       nsize: 2
757       args: -ksp_type bcgsl -ksp_monitor_short -da_refine 2 -ksp_bcgsl_ell 3 -snes_view
758       requires: !single
759 
760    test:
761       suffix: bcols1
762       nsize: 2
763       args: -da_refine 3 -snes_monitor_short -pc_type mg -ksp_type fgmres -pc_mg_type full -mat_fd_coloring_bcols 1
764       output_file: output/ex19_1.out
765       requires: !single
766 
767    test:
768       suffix: bjacobi
769       nsize: 4
770       args: -da_refine 4 -ksp_type fgmres -pc_type bjacobi -pc_bjacobi_blocks 2 -sub_ksp_type gmres -sub_ksp_max_it 2 -sub_pc_type bjacobi -sub_sub_ksp_type preonly -sub_sub_pc_type ilu -snes_monitor_short
771       requires: !single
772 
773    test:
774       suffix: cgne
775       args: -da_refine 2 -pc_type lu -ksp_type cgne -ksp_monitor_short -ksp_converged_reason -ksp_view -ksp_norm_type unpreconditioned
776       filter: grep -v HERMITIAN
777       requires: !single
778 
779    test:
780       suffix: cgs
781       args: -da_refine 1 -ksp_monitor_short -ksp_type cgs
782       requires: !single
783 
784    test:
785       suffix: composite_fieldsplit
786       args: -ksp_type fgmres -pc_type composite -pc_composite_type MULTIPLICATIVE -pc_composite_pcs fieldsplit,none -sub_0_pc_fieldsplit_block_size 4 -sub_0_pc_fieldsplit_type additive -sub_0_pc_fieldsplit_0_fields 0,1,2 -sub_0_pc_fieldsplit_1_fields 3 -snes_monitor_short -ksp_monitor_short
787       requires: !single
788 
789    test:
790       suffix: composite_fieldsplit_bjacobi
791       args: -ksp_type fgmres -pc_type composite -pc_composite_type MULTIPLICATIVE -pc_composite_pcs fieldsplit,bjacobi -sub_0_pc_fieldsplit_block_size 4 -sub_0_pc_fieldsplit_type additive -sub_0_pc_fieldsplit_0_fields 0,1,2 -sub_0_pc_fieldsplit_1_fields 3 -sub_1_pc_bjacobi_blocks 16 -sub_1_sub_pc_type lu -snes_monitor_short -ksp_monitor_short
792       requires: !single
793 
794    test:
795       suffix: composite_fieldsplit_bjacobi_2
796       nsize: 4
797       args: -ksp_type fgmres -pc_type composite -pc_composite_type MULTIPLICATIVE -pc_composite_pcs fieldsplit,bjacobi -sub_0_pc_fieldsplit_block_size 4 -sub_0_pc_fieldsplit_type additive -sub_0_pc_fieldsplit_0_fields 0,1,2 -sub_0_pc_fieldsplit_1_fields 3 -sub_1_pc_bjacobi_blocks 16 -sub_1_sub_pc_type lu -snes_monitor_short -ksp_monitor_short
798       requires: !single
799 
800    test:
801       suffix: composite_gs_newton
802       nsize: 2
803       args: -da_refine 3 -grashof 4e4 -lidvelocity 100 -snes_monitor_short -snes_type composite -snes_composite_type additiveoptimal -snes_composite_sneses ngs,newtonls -sub_0_snes_max_it 20 -sub_1_pc_type mg
804       requires: !single
805 
806    test:
807       suffix: cuda
808       requires: cuda !single
809       args: -dm_vec_type cuda -dm_mat_type aijcusparse -pc_type none -ksp_type fgmres -snes_monitor_short -snes_rtol 1.e-5
810 
811    test:
812       suffix: draw
813       args: -pc_type fieldsplit -snes_view draw -fieldsplit_x_velocity_pc_type mg -fieldsplit_x_velocity_pc_mg_galerkin pmat -fieldsplit_x_velocity_pc_mg_levels 2 -da_refine 1 -fieldsplit_x_velocity_mg_coarse_pc_type svd
814       requires: x !single
815 
816    test:
817       suffix: drawports
818       args: -snes_monitor_solution draw::draw_ports -da_refine 1
819       output_file: output/ex19_draw.out
820       requires: x !single
821 
822    test:
823       suffix: fas
824       args: -da_refine 4 -snes_monitor_short -snes_type fas -fas_levels_snes_type ngs -fas_levels_snes_ngs_sweeps 3 -fas_levels_snes_ngs_atol 0.0 -fas_levels_snes_ngs_stol 0.0 -grashof 4e4 -snes_fas_smoothup 6 -snes_fas_smoothdown 6 -lidvelocity 100
825       requires: !single
826 
827    test:
828       suffix: fas_full
829       args: -da_refine 4 -snes_monitor_short -snes_type fas -snes_fas_type full -snes_fas_full_downsweep -fas_levels_snes_type ngs -fas_levels_snes_ngs_sweeps 3 -fas_levels_snes_ngs_atol 0.0 -fas_levels_snes_ngs_stol 0.0 -grashof 4e4 -snes_fas_smoothup 6 -snes_fas_smoothdown 6 -lidvelocity 100
830       requires: !single
831 
832    test:
833       suffix: fdcoloring_ds
834       args: -da_refine 3 -snes_converged_reason -pc_type mg -mat_fd_type ds
835       output_file: output/ex19_2.out
836       requires: !single
837 
838    test:
839       suffix: fdcoloring_ds_baij
840       args: -da_refine 3 -snes_converged_reason -pc_type mg -mat_fd_type ds -dm_mat_type baij
841       output_file: output/ex19_2.out
842       requires: !single
843 
844    test:
845       suffix: fdcoloring_ds_bcols1
846       args: -da_refine 3 -snes_converged_reason -pc_type mg -mat_fd_type ds -mat_fd_coloring_bcols 1
847       output_file: output/ex19_2.out
848       requires: !single
849 
850    test:
851       suffix: fdcoloring_wp
852       args: -da_refine 3 -snes_monitor_short -pc_type mg
853       requires: !single
854 
855    test:
856       suffix: fdcoloring_wp_baij
857       args: -da_refine 3 -snes_monitor_short -pc_type mg -dm_mat_type baij
858       output_file: output/ex19_fdcoloring_wp.out
859       requires: !single
860 
861    test:
862       suffix: fdcoloring_wp_bcols1
863       args: -da_refine 3 -snes_monitor_short -pc_type mg -mat_fd_coloring_bcols 1
864       output_file: output/ex19_fdcoloring_wp.out
865       requires: !single
866 
867    test:
868       suffix: fieldsplit_2
869       args: -ksp_type fgmres -pc_type fieldsplit -pc_fieldsplit_block_size 4 -pc_fieldsplit_type additive -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -snes_monitor_short -ksp_monitor_short
870       requires: !single
871 
872    test:
873       suffix: fieldsplit_3
874       args: -ksp_type fgmres -pc_type fieldsplit -pc_fieldsplit_block_size 4 -pc_fieldsplit_type additive -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -fieldsplit_0_pc_type lu -fieldsplit_1_pc_type lu -snes_monitor_short -ksp_monitor_short
875       requires: !single
876 
877    test:
878       suffix: fieldsplit_4
879       args: -ksp_type fgmres -pc_type fieldsplit -pc_fieldsplit_block_size 4 -pc_fieldsplit_type SCHUR -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -fieldsplit_0_pc_type lu -fieldsplit_1_pc_type lu -snes_monitor_short -ksp_monitor_short
880       requires: !single
881 
882    # HYPRE PtAP broken with complex numbers
883    test:
884       suffix: fieldsplit_hypre
885       nsize: 2
886       requires: hypre mumps !complex !defined(PETSC_HAVE_HYPRE_DEVICE)
887       args: -pc_type fieldsplit -pc_fieldsplit_block_size 4 -pc_fieldsplit_type SCHUR -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -fieldsplit_0_pc_type lu -fieldsplit_0_pc_factor_mat_solver_type mumps -fieldsplit_1_pc_type hypre -fieldsplit_1_pc_hypre_type boomeramg -snes_monitor_short -ksp_monitor_short
888 
889    test:
890       suffix: fieldsplit_mumps
891       nsize: 2
892       requires: mumps
893       args: -pc_type fieldsplit -pc_fieldsplit_block_size 4 -pc_fieldsplit_type SCHUR -pc_fieldsplit_0_fields 0,1,2 -pc_fieldsplit_1_fields 3 -fieldsplit_0_pc_type lu -fieldsplit_1_pc_type lu -snes_monitor_short -ksp_monitor_short -fieldsplit_0_pc_factor_mat_solver_type mumps -fieldsplit_1_pc_factor_mat_solver_type mumps
894       output_file: output/ex19_fieldsplit_5.out
895 
896    test:
897       suffix: greedy_coloring
898       nsize: 2
899       args: -da_refine 3 -snes_monitor_short -snes_fd_color -snes_fd_color_use_mat -mat_coloring_type greedy -mat_coloring_weight_type lf -mat_coloring_view> ex19_greedy_coloring.tmp 2>&1
900       requires: !single
901 
902    # HYPRE PtAP broken with complex numbers
903    test:
904       suffix: hypre
905       nsize: 2
906       requires: hypre !complex !defined(PETSC_HAVE_HYPRE_DEVICE)
907       args: -da_refine 3 -snes_monitor_short -pc_type hypre -ksp_norm_type unpreconditioned
908 
909    # ibcgs is broken when using device vectors
910    test:
911       suffix: ibcgs
912       nsize: 2
913       args: -ksp_type ibcgs -ksp_monitor_short -da_refine 2 -snes_view
914       requires: !complex !single
915 
916    test:
917       suffix: kaczmarz
918       nsize: 2
919       args: -pc_type kaczmarz -ksp_monitor_short -snes_monitor_short -snes_view
920       requires: !single
921 
922    test:
923       suffix: klu
924       requires: suitesparse
925       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type klu
926       output_file: output/ex19_superlu.out
927 
928    test:
929       suffix: klu_2
930       requires: suitesparse
931       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type klu -pc_factor_mat_ordering_type nd
932       output_file: output/ex19_superlu.out
933 
934    test:
935       suffix: klu_3
936       requires: suitesparse
937       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type klu -mat_klu_use_btf 0
938       output_file: output/ex19_superlu.out
939 
940    test:
941       suffix: ml
942       nsize: 2
943       requires: ml
944       args: -da_refine 3 -snes_monitor_short -pc_type ml
945 
946    test:
947       suffix: ngmres_fas
948       args: -da_refine 4 -snes_monitor_short -snes_type ngmres -npc_fas_levels_snes_type ngs -npc_fas_levels_snes_ngs_sweeps 3 -npc_fas_levels_snes_ngs_atol 0.0 -npc_fas_levels_snes_ngs_stol 0.0 -npc_snes_type fas -npc_fas_levels_snes_type ngs -npc_snes_max_it 1 -npc_snes_fas_smoothup 6 -npc_snes_fas_smoothdown 6 -lidvelocity 100 -grashof 4e4
949       requires: !single
950 
951    test:
952       suffix: ngmres_fas_gssecant
953       args: -da_refine 3 -snes_monitor_short -snes_type ngmres -npc_snes_type fas -npc_fas_levels_snes_type ngs -npc_fas_levels_snes_max_it 6 -npc_fas_levels_snes_ngs_secant -npc_fas_levels_snes_ngs_max_it 1 -npc_fas_coarse_snes_max_it 1 -lidvelocity 100 -grashof 4e4
954       requires: !single
955 
956    test:
957       suffix: ngmres_fas_ms
958       nsize: 2
959       args: -snes_grid_sequence 2 -lidvelocity 200 -grashof 1e4 -snes_monitor_short -snes_view -snes_converged_reason -snes_type ngmres -npc_snes_type fas -npc_fas_coarse_snes_type newtonls -npc_fas_coarse_ksp_type preonly -npc_snes_max_it 1
960       requires: !single
961 
962    test:
963       suffix: ngmres_nasm
964       nsize: 4
965       args: -da_refine 4 -da_overlap 2 -snes_monitor_short -snes_type ngmres -snes_max_it 10 -npc_snes_type nasm -npc_snes_nasm_type basic -grashof 4e4 -lidvelocity 100
966       requires: !single
967 
968    test:
969       suffix: ngs
970       args: -snes_type ngs -snes_view -snes_monitor -snes_rtol 1e-4
971       requires: !single
972 
973    test:
974       suffix: ngs_fd
975       args: -snes_type ngs -snes_ngs_secant -snes_view -snes_monitor -snes_rtol 1e-4
976       requires: !single
977 
978    test:
979       suffix: parms
980       nsize: 2
981       requires: parms
982       args: -pc_type parms -ksp_monitor_short -snes_view
983 
984    test:
985       suffix: superlu
986       requires: superlu
987       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type superlu
988 
989    test:
990       suffix: superlu_sell
991       requires: superlu
992       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type superlu -dm_mat_type sell -pc_factor_mat_ordering_type natural
993       output_file: output/ex19_superlu.out
994 
995    test:
996       suffix: superlu_dist
997       requires: superlu_dist
998       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type superlu_dist
999       output_file: output/ex19_superlu.out
1000 
1001    test:
1002       suffix: superlu_dist_2
1003       nsize: 2
1004       requires: superlu_dist
1005       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type superlu_dist
1006       output_file: output/ex19_superlu.out
1007 
1008    test:
1009       suffix: superlu_dist_3d
1010       nsize: 4
1011       requires: superlu_dist !defined(PETSCTEST_VALGRIND)
1012       filter: grep -v iam | grep -v openMP
1013       args: -da_grid_x 20 -da_grid_y 20 -pc_type lu -pc_factor_mat_solver_type superlu_dist -mat_superlu_dist_3d -mat_superlu_dist_d 2 -snes_view -snes_monitor -ksp_monitor
1014 
1015    test:
1016       suffix: superlu_equil
1017       requires: superlu
1018       args: -da_grid_x 20 -da_grid_y 20 -{snes,ksp}_monitor_short -pc_type lu -pc_factor_mat_solver_type superlu -mat_superlu_equil
1019 
1020    test:
1021       suffix: superlu_equil_sell
1022       requires: superlu
1023       args: -da_grid_x 20 -da_grid_y 20 -{snes,ksp}_monitor_short -pc_type lu -pc_factor_mat_solver_type superlu -mat_superlu_equil -dm_mat_type sell -pc_factor_mat_ordering_type natural
1024       output_file: output/ex19_superlu_equil.out
1025 
1026    test:
1027       suffix: tcqmr
1028       args: -da_refine 1 -ksp_monitor_short -ksp_type tcqmr
1029       requires: !single
1030 
1031    test:
1032       suffix: tfqmr
1033       args: -da_refine 1 -ksp_monitor_short -ksp_type tfqmr
1034       requires: !single
1035 
1036    test:
1037       suffix: umfpack
1038       requires: suitesparse
1039       args: -da_refine 2 -pc_type lu -pc_factor_mat_solver_type umfpack -snes_view -snes_monitor_short -ksp_monitor_short -pc_factor_mat_ordering_type external
1040 
1041    test:
1042       suffix: tut_1
1043       nsize: 4
1044       requires: !single
1045       args: -da_refine 5 -snes_monitor -ksp_monitor -snes_view
1046 
1047    test:
1048       suffix: tut_2
1049       nsize: 4
1050       requires: !single
1051       args: -da_refine 5 -snes_monitor -ksp_monitor -snes_view -pc_type mg
1052 
1053    # HYPRE PtAP broken with complex numbers
1054    test:
1055       suffix: tut_3
1056       nsize: 4
1057       requires: hypre !single !complex !defined(PETSC_HAVE_HYPRE_DEVICE)
1058       args: -da_refine 5 -snes_monitor -ksp_monitor -snes_view -pc_type hypre
1059 
1060    test:
1061       suffix: tut_8
1062       nsize: 4
1063       requires: ml !single
1064       args: -da_refine 5 -snes_monitor -ksp_monitor -snes_view -pc_type ml
1065 
1066    test:
1067       suffix: tut_4
1068       nsize: 1
1069       requires: !single
1070       args: -da_refine 5 -log_view
1071       filter: head -n 2
1072       filter_output: head -n 2
1073 
1074    test:
1075       suffix: tut_5
1076       nsize: 1
1077       requires: !single
1078       args: -da_refine 5 -log_view -pc_type mg
1079       filter: head -n 2
1080       filter_output: head -n 2
1081 
1082    test:
1083       suffix: tut_6
1084       nsize: 4
1085       requires: !single
1086       args: -da_refine 5 -log_view
1087       filter: head -n 2
1088       filter_output: head -n 2
1089 
1090    test:
1091       suffix: tut_7
1092       nsize: 4
1093       requires: !single
1094       args: -da_refine 5 -log_view -pc_type mg
1095       filter: head -n 2
1096       filter_output: head -n 2
1097 
1098    test:
1099       suffix: cuda_1
1100       nsize: 1
1101       requires: cuda
1102       args: -snes_monitor -dm_mat_type seqaijcusparse -dm_vec_type seqcuda -pc_type gamg -pc_gamg_esteig_ksp_max_it 10 -ksp_monitor -mg_levels_ksp_max_it 3
1103 
1104    test:
1105       suffix: cuda_2
1106       nsize: 3
1107       requires: cuda !single
1108       args: -snes_monitor -dm_mat_type mpiaijcusparse -dm_vec_type mpicuda -pc_type gamg -pc_gamg_esteig_ksp_max_it 10 -ksp_monitor  -mg_levels_ksp_max_it 3
1109 
1110    test:
1111       suffix: cuda_dm_bind_below
1112       nsize: 2
1113       requires: cuda
1114       args: -dm_mat_type aijcusparse -dm_vec_type cuda -da_refine 3 -pc_type mg -mg_levels_ksp_type chebyshev -mg_levels_pc_type jacobi -log_view -pc_mg_log -dm_bind_below 10000
1115       filter: awk "/Level/ {print \$NF}"
1116 
1117    test:
1118       suffix: viennacl_dm_bind_below
1119       nsize: 2
1120       requires: viennacl
1121       args: -dm_mat_type aijviennacl -dm_vec_type viennacl -da_refine 3 -pc_type mg -mg_levels_ksp_type chebyshev -mg_levels_pc_type jacobi -log_view -pc_mg_log -dm_bind_below 10000
1122       filter: awk "/Level/ {print \$NF}"
1123 
1124    test:
1125       suffix: seqbaijmkl
1126       nsize: 1
1127       requires: defined(PETSC_HAVE_MKL_SPARSE_OPTIMIZE)
1128       args: -dm_mat_type baij -snes_monitor -ksp_monitor -snes_view
1129 
1130    test:
1131       suffix: mpibaijmkl
1132       nsize: 2
1133       requires:  defined(PETSC_HAVE_MKL_SPARSE_OPTIMIZE)
1134       args: -dm_mat_type baij -snes_monitor -ksp_monitor -snes_view
1135 
1136    test:
1137      suffix: cpardiso
1138      nsize: 4
1139      requires: mkl_cpardiso
1140      args: -pc_type lu -pc_factor_mat_solver_type mkl_cpardiso -ksp_monitor
1141 
1142    test:
1143      suffix: logviewmemory
1144      requires: defined(PETSC_USE_LOG) !defined(PETSCTEST_VALGRIND)
1145      args: -log_view -log_view_memory -da_refine 4
1146      filter: grep MatFDColorSetUp | wc -w | xargs  -I % sh -c "expr % \> 21"
1147 
1148    test:
1149      suffix: fs
1150      args: -pc_type fieldsplit -da_refine 3  -all_ksp_monitor -fieldsplit_y_velocity_pc_type lu  -fieldsplit_temperature_pc_type lu -fieldsplit_x_velocity_pc_type lu  -snes_view
1151 
1152    test:
1153      suffix: asm_matconvert
1154      args: -mat_type aij -pc_type asm -pc_asm_sub_mat_type dense -snes_view
1155 
1156    test:
1157       suffix: euclid
1158       nsize: 2
1159       requires: hypre !single !complex !defined(PETSC_HAVE_HYPRE_MIXEDINT) !defined(PETSC_HAVE_HYPRE_DEVICE)
1160       args: -da_refine 2 -ksp_monitor -snes_monitor -snes_view -pc_type hypre -pc_hypre_type euclid
1161 
1162    test:
1163       suffix: euclid_bj
1164       nsize: 2
1165       requires: hypre !single !complex !defined(PETSC_HAVE_HYPRE_MIXEDINT) !defined(PETSC_HAVE_HYPRE_DEVICE)
1166       args: -da_refine 2 -ksp_monitor -snes_monitor -snes_view -pc_type hypre -pc_hypre_type euclid -pc_hypre_euclid_bj
1167 
1168    test:
1169       suffix: euclid_droptolerance
1170       nsize: 1
1171       requires: hypre !single !complex !defined(PETSC_HAVE_HYPRE_MIXEDINT) !defined(PETSC_HAVE_HYPRE_DEVICE)
1172       args: -da_refine 2 -ksp_monitor -snes_monitor -snes_view -pc_type hypre -pc_hypre_type euclid -pc_hypre_euclid_droptolerance .1
1173 
1174    test:
1175       suffix: failure_size
1176       nsize: 1
1177       requires: !defined(PETSC_USE_64BIT_INDICES) !defined(PETSCTEST_VALGRIND)
1178       args: -da_refine 100 -petsc_ci_portable_error_output -error_output_stdout
1179       filter: egrep -v "(options_left|memory block|leaked context|is not freed before MPI_Finalize|Could be the program crashed)"
1180 
1181 TEST*/
1182