1 static char help[] = "Bratu nonlinear PDE in 3d.\n\
2 We solve the Bratu (SFI - solid fuel ignition) problem in a 3D rectangular\n\
3 domain, using distributed arrays (DMDAs) to partition the parallel grid.\n\
4 The command line options include:\n\
5 -par <parameter>, where <parameter> indicates the problem's nonlinearity\n\
6 problem SFI: <parameter> = Bratu parameter (0 <= par <= 6.81)\n\n";
7
8 /* ------------------------------------------------------------------------
9
10 Solid Fuel Ignition (SFI) problem. This problem is modeled by
11 the partial differential equation
12
13 -Laplacian u - lambda*exp(u) = 0, 0 < x,y < 1,
14
15 with boundary conditions
16
17 u = 0 for x = 0, x = 1, y = 0, y = 1, z = 0, z = 1
18
19 A finite difference approximation with the usual 7-point stencil
20 is used to discretize the boundary value problem to obtain a nonlinear
21 system of equations.
22
23 ------------------------------------------------------------------------- */
24
25 /*
26 Include "petscdmda.h" so that we can use distributed arrays (DMDAs).
27 Include "petscsnes.h" so that we can use SNES solvers. Note that this
28 file automatically includes:
29 petscsys.h - base PETSc routines petscvec.h - vectors
30 petscmat.h - matrices
31 petscis.h - index sets petscksp.h - Krylov subspace methods
32 petscviewer.h - viewers petscpc.h - preconditioners
33 petscksp.h - linear solvers
34 */
35 #include <petscdm.h>
36 #include <petscdmda.h>
37 #include <petscsnes.h>
38
39 /*
40 User-defined application context - contains data needed by the
41 application-provided call-back routines, FormJacobian() and
42 FormFunction().
43 */
44 typedef struct {
45 PetscReal param; /* test problem parameter */
46 DM da; /* distributed array data structure */
47 } AppCtx;
48
49 /*
50 User-defined routines
51 */
52 extern PetscErrorCode FormFunctionLocal(SNES, Vec, Vec, void *);
53 extern PetscErrorCode FormFunction(SNES, Vec, Vec, void *);
54 extern PetscErrorCode FormInitialGuess(AppCtx *, Vec);
55 extern PetscErrorCode FormJacobian(SNES, Vec, Mat, Mat, void *);
56
main(int argc,char ** argv)57 int main(int argc, char **argv)
58 {
59 SNES snes; /* nonlinear solver */
60 Vec x, r; /* solution, residual vectors */
61 Mat J = NULL; /* Jacobian matrix */
62 AppCtx user; /* user-defined work context */
63 PetscInt its; /* iterations for convergence */
64 MatFDColoring matfdcoloring = NULL;
65 PetscBool matrix_free = PETSC_FALSE, coloring = PETSC_FALSE, coloring_ds = PETSC_FALSE, local_coloring = PETSC_FALSE;
66 PetscReal bratu_lambda_max = 6.81, bratu_lambda_min = 0., fnorm;
67
68 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
69 Initialize program
70 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
71
72 PetscFunctionBeginUser;
73 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
74
75 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
76 Initialize problem parameters
77 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
78 user.param = 6.0;
79 PetscCall(PetscOptionsGetReal(NULL, NULL, "-par", &user.param, NULL));
80 PetscCheck(user.param < bratu_lambda_max && user.param > bratu_lambda_min, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lambda is out of range");
81
82 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
83 Create nonlinear solver context
84 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
85 PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
86
87 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
88 Create distributed array (DMDA) to manage parallel grid and vectors
89 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
90 PetscCall(DMDACreate3d(PETSC_COMM_WORLD, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DM_BOUNDARY_NONE, DMDA_STENCIL_STAR, 4, 4, 4, PETSC_DECIDE, PETSC_DECIDE, PETSC_DECIDE, 1, 1, NULL, NULL, NULL, &user.da));
91 PetscCall(DMSetFromOptions(user.da));
92 PetscCall(DMSetUp(user.da));
93 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
94 Extract global vectors from DMDA; then duplicate for remaining
95 vectors that are the same types
96 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
97 PetscCall(DMCreateGlobalVector(user.da, &x));
98 PetscCall(VecDuplicate(x, &r));
99
100 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
101 Set function evaluation routine and vector
102 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
103 PetscCall(SNESSetFunction(snes, r, FormFunction, (void *)&user));
104
105 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
106 Create matrix data structure; set Jacobian evaluation routine
107
108 Set Jacobian matrix data structure and default Jacobian evaluation
109 routine. User can override with:
110 -snes_mf : matrix-free Newton-Krylov method with no preconditioning
111 (unless user explicitly sets preconditioner)
112 -snes_mf_operator : form matrix used to construct the preconditioner as set by the user,
113 but use matrix-free approx for Jacobian-vector
114 products within Newton-Krylov method
115 -fdcoloring : using finite differences with coloring to compute the Jacobian
116
117 Note one can use -matfd_coloring wp or ds the only reason for the -fdcoloring_ds option
118 below is to test the call to MatFDColoringSetType().
119 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
120 PetscCall(PetscOptionsGetBool(NULL, NULL, "-snes_mf", &matrix_free, NULL));
121 PetscCall(PetscOptionsGetBool(NULL, NULL, "-fdcoloring", &coloring, NULL));
122 PetscCall(PetscOptionsGetBool(NULL, NULL, "-fdcoloring_ds", &coloring_ds, NULL));
123 PetscCall(PetscOptionsGetBool(NULL, NULL, "-fdcoloring_local", &local_coloring, NULL));
124 if (!matrix_free) {
125 PetscCall(DMSetMatType(user.da, MATAIJ));
126 PetscCall(DMCreateMatrix(user.da, &J));
127 if (coloring) {
128 ISColoring iscoloring;
129 if (!local_coloring) {
130 PetscCall(DMCreateColoring(user.da, IS_COLORING_GLOBAL, &iscoloring));
131 PetscCall(MatFDColoringCreate(J, iscoloring, &matfdcoloring));
132 PetscCall(MatFDColoringSetFunction(matfdcoloring, (MatFDColoringFn *)FormFunction, &user));
133 } else {
134 PetscCall(DMCreateColoring(user.da, IS_COLORING_LOCAL, &iscoloring));
135 PetscCall(MatFDColoringCreate(J, iscoloring, &matfdcoloring));
136 PetscCall(MatFDColoringUseDM(J, matfdcoloring));
137 PetscCall(MatFDColoringSetFunction(matfdcoloring, (MatFDColoringFn *)FormFunctionLocal, &user));
138 }
139 if (coloring_ds) PetscCall(MatFDColoringSetType(matfdcoloring, MATMFFD_DS));
140 PetscCall(MatFDColoringSetFromOptions(matfdcoloring));
141 PetscCall(MatFDColoringSetUp(J, iscoloring, matfdcoloring));
142 PetscCall(SNESSetJacobian(snes, J, J, SNESComputeJacobianDefaultColor, matfdcoloring));
143 PetscCall(ISColoringDestroy(&iscoloring));
144 } else {
145 PetscCall(SNESSetJacobian(snes, J, J, FormJacobian, &user));
146 }
147 }
148
149 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
150 Customize nonlinear solver; set runtime options
151 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
152 PetscCall(SNESSetDM(snes, user.da));
153 PetscCall(SNESSetFromOptions(snes));
154
155 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
156 Evaluate initial guess
157 Note: The user should initialize the vector, x, with the initial guess
158 for the nonlinear solver prior to calling SNESSolve(). In particular,
159 to employ an initial guess of zero, the user should explicitly set
160 this vector to zero by calling VecSet().
161 */
162 PetscCall(FormInitialGuess(&user, x));
163
164 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
165 Solve nonlinear system
166 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
167 PetscCall(SNESSolve(snes, NULL, x));
168 PetscCall(SNESGetIterationNumber(snes, &its));
169
170 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
171 Explicitly check norm of the residual of the solution
172 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
173 PetscCall(FormFunction(snes, x, r, (void *)&user));
174 PetscCall(VecNorm(r, NORM_2, &fnorm));
175 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of SNES iterations = %" PetscInt_FMT " fnorm %g\n", its, (double)fnorm));
176
177 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
178 Free work space. All PETSc objects should be destroyed when they
179 are no longer needed.
180 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
181
182 PetscCall(MatDestroy(&J));
183 PetscCall(VecDestroy(&x));
184 PetscCall(VecDestroy(&r));
185 PetscCall(SNESDestroy(&snes));
186 PetscCall(DMDestroy(&user.da));
187 PetscCall(MatFDColoringDestroy(&matfdcoloring));
188 PetscCall(PetscFinalize());
189 return 0;
190 }
191 /* ------------------------------------------------------------------- */
192 /*
193 FormInitialGuess - Forms initial approximation.
194
195 Input Parameters:
196 user - user-defined application context
197 X - vector
198
199 Output Parameter:
200 X - vector
201 */
FormInitialGuess(AppCtx * user,Vec X)202 PetscErrorCode FormInitialGuess(AppCtx *user, Vec X)
203 {
204 PetscInt i, j, k, Mx, My, Mz, xs, ys, zs, xm, ym, zm;
205 PetscReal lambda, temp1, hx, hy, hz, tempk, tempj;
206 PetscScalar ***x;
207
208 PetscFunctionBeginUser;
209 PetscCall(DMDAGetInfo(user->da, PETSC_IGNORE, &Mx, &My, &Mz, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE));
210
211 lambda = user->param;
212 hx = 1.0 / (PetscReal)(Mx - 1);
213 hy = 1.0 / (PetscReal)(My - 1);
214 hz = 1.0 / (PetscReal)(Mz - 1);
215 temp1 = lambda / (lambda + 1.0);
216
217 /*
218 Get a pointer to vector data.
219 - For default PETSc vectors, VecGetArray() returns a pointer to
220 the data array. Otherwise, the routine is implementation dependent.
221 - You MUST call VecRestoreArray() when you no longer need access to
222 the array.
223 */
224 PetscCall(DMDAVecGetArray(user->da, X, &x));
225
226 /*
227 Get local grid boundaries (for 3-dimensional DMDA):
228 xs, ys, zs - starting grid indices (no ghost points)
229 xm, ym, zm - widths of local grid (no ghost points)
230
231 */
232 PetscCall(DMDAGetCorners(user->da, &xs, &ys, &zs, &xm, &ym, &zm));
233
234 /*
235 Compute initial guess over the locally owned part of the grid
236 */
237 for (k = zs; k < zs + zm; k++) {
238 tempk = (PetscReal)(PetscMin(k, Mz - k - 1)) * hz;
239 for (j = ys; j < ys + ym; j++) {
240 tempj = PetscMin((PetscReal)(PetscMin(j, My - j - 1)) * hy, tempk);
241 for (i = xs; i < xs + xm; i++) {
242 if (i == 0 || j == 0 || k == 0 || i == Mx - 1 || j == My - 1 || k == Mz - 1) {
243 /* boundary conditions are all zero Dirichlet */
244 x[k][j][i] = 0.0;
245 } else {
246 x[k][j][i] = temp1 * PetscSqrtReal(PetscMin((PetscReal)(PetscMin(i, Mx - i - 1)) * hx, tempj));
247 }
248 }
249 }
250 }
251
252 /*
253 Restore vector
254 */
255 PetscCall(DMDAVecRestoreArray(user->da, X, &x));
256 PetscFunctionReturn(PETSC_SUCCESS);
257 }
258 /* ------------------------------------------------------------------- */
259 /*
260 FormFunctionLocal - Evaluates nonlinear function, F(x) on a ghosted patch
261
262 Input Parameters:
263 . snes - the SNES context
264 . localX - input vector, this contains the ghosted region needed
265 . ptr - optional user-defined context, as set by SNESSetFunction()
266
267 Output Parameter:
268 . F - function vector, this does not contain a ghosted region
269 */
FormFunctionLocal(SNES snes,Vec localX,Vec F,void * ptr)270 PetscErrorCode FormFunctionLocal(SNES snes, Vec localX, Vec F, void *ptr)
271 {
272 AppCtx *user = (AppCtx *)ptr;
273 PetscInt i, j, k, Mx, My, Mz, xs, ys, zs, xm, ym, zm;
274 PetscReal two = 2.0, lambda, hx, hy, hz, hxhzdhy, hyhzdhx, hxhydhz, sc;
275 PetscScalar u_north, u_south, u_east, u_west, u_up, u_down, u, u_xx, u_yy, u_zz, ***x, ***f;
276 DM da;
277
278 PetscFunctionBeginUser;
279 PetscCall(SNESGetDM(snes, &da));
280 PetscCall(DMDAGetInfo(da, PETSC_IGNORE, &Mx, &My, &Mz, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE));
281
282 lambda = user->param;
283 hx = 1.0 / (PetscReal)(Mx - 1);
284 hy = 1.0 / (PetscReal)(My - 1);
285 hz = 1.0 / (PetscReal)(Mz - 1);
286 sc = hx * hy * hz * lambda;
287 hxhzdhy = hx * hz / hy;
288 hyhzdhx = hy * hz / hx;
289 hxhydhz = hx * hy / hz;
290
291 /*
292 Get pointers to vector data
293 */
294 PetscCall(DMDAVecGetArrayRead(da, localX, &x));
295 PetscCall(DMDAVecGetArray(da, F, &f));
296
297 /*
298 Get local grid boundaries
299 */
300 PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
301
302 /*
303 Compute function over the locally owned part of the grid
304 */
305 for (k = zs; k < zs + zm; k++) {
306 for (j = ys; j < ys + ym; j++) {
307 for (i = xs; i < xs + xm; i++) {
308 if (i == 0 || j == 0 || k == 0 || i == Mx - 1 || j == My - 1 || k == Mz - 1) {
309 f[k][j][i] = x[k][j][i];
310 } else {
311 u = x[k][j][i];
312 u_east = x[k][j][i + 1];
313 u_west = x[k][j][i - 1];
314 u_north = x[k][j + 1][i];
315 u_south = x[k][j - 1][i];
316 u_up = x[k + 1][j][i];
317 u_down = x[k - 1][j][i];
318 u_xx = (-u_east + two * u - u_west) * hyhzdhx;
319 u_yy = (-u_north + two * u - u_south) * hxhzdhy;
320 u_zz = (-u_up + two * u - u_down) * hxhydhz;
321 f[k][j][i] = u_xx + u_yy + u_zz - sc * PetscExpScalar(u);
322 }
323 }
324 }
325 }
326
327 /*
328 Restore vectors
329 */
330 PetscCall(DMDAVecRestoreArrayRead(da, localX, &x));
331 PetscCall(DMDAVecRestoreArray(da, F, &f));
332 PetscCall(PetscLogFlops(11.0 * ym * xm));
333 PetscFunctionReturn(PETSC_SUCCESS);
334 }
335 /* ------------------------------------------------------------------- */
336 /*
337 FormFunction - Evaluates nonlinear function, F(x) on the entire domain
338
339 Input Parameters:
340 . snes - the SNES context
341 . X - input vector
342 . ptr - optional user-defined context, as set by SNESSetFunction()
343
344 Output Parameter:
345 . F - function vector
346 */
FormFunction(SNES snes,Vec X,Vec F,void * ptr)347 PetscErrorCode FormFunction(SNES snes, Vec X, Vec F, void *ptr)
348 {
349 Vec localX;
350 DM da;
351
352 PetscFunctionBeginUser;
353 PetscCall(SNESGetDM(snes, &da));
354 PetscCall(DMGetLocalVector(da, &localX));
355
356 /*
357 Scatter ghost points to local vector,using the 2-step process
358 DMGlobalToLocalBegin(),DMGlobalToLocalEnd().
359 By placing code between these two statements, computations can be
360 done while messages are in transition.
361 */
362 PetscCall(DMGlobalToLocalBegin(da, X, INSERT_VALUES, localX));
363 PetscCall(DMGlobalToLocalEnd(da, X, INSERT_VALUES, localX));
364
365 PetscCall(FormFunctionLocal(snes, localX, F, ptr));
366 PetscCall(DMRestoreLocalVector(da, &localX));
367 PetscFunctionReturn(PETSC_SUCCESS);
368 }
369 /* ------------------------------------------------------------------- */
370 /*
371 FormJacobian - Evaluates Jacobian matrix.
372
373 Input Parameters:
374 . snes - the SNES context
375 . x - input vector
376 . ptr - optional user-defined context, as set by SNESSetJacobian()
377
378 Output Parameters:
379 . A - Jacobian matrix
380 . B - optionally different matrix used to construct the preconditioner
381
382 */
FormJacobian(SNES snes,Vec X,Mat J,Mat jac,void * ptr)383 PetscErrorCode FormJacobian(SNES snes, Vec X, Mat J, Mat jac, void *ptr)
384 {
385 AppCtx *user = (AppCtx *)ptr; /* user-defined application context */
386 Vec localX;
387 PetscInt i, j, k, Mx, My, Mz;
388 MatStencil col[7], row;
389 PetscInt xs, ys, zs, xm, ym, zm;
390 PetscScalar lambda, v[7], hx, hy, hz, hxhzdhy, hyhzdhx, hxhydhz, sc, ***x;
391 DM da;
392
393 PetscFunctionBeginUser;
394 PetscCall(SNESGetDM(snes, &da));
395 PetscCall(DMGetLocalVector(da, &localX));
396 PetscCall(DMDAGetInfo(da, PETSC_IGNORE, &Mx, &My, &Mz, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE, PETSC_IGNORE));
397
398 lambda = user->param;
399 hx = 1.0 / (PetscReal)(Mx - 1);
400 hy = 1.0 / (PetscReal)(My - 1);
401 hz = 1.0 / (PetscReal)(Mz - 1);
402 sc = hx * hy * hz * lambda;
403 hxhzdhy = hx * hz / hy;
404 hyhzdhx = hy * hz / hx;
405 hxhydhz = hx * hy / hz;
406
407 /*
408 Scatter ghost points to local vector, using the 2-step process
409 DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
410 By placing code between these two statements, computations can be
411 done while messages are in transition.
412 */
413 PetscCall(DMGlobalToLocalBegin(da, X, INSERT_VALUES, localX));
414 PetscCall(DMGlobalToLocalEnd(da, X, INSERT_VALUES, localX));
415
416 /*
417 Get pointer to vector data
418 */
419 PetscCall(DMDAVecGetArrayRead(da, localX, &x));
420
421 /*
422 Get local grid boundaries
423 */
424 PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
425
426 /*
427 Compute entries for the locally owned part of the Jacobian.
428 - Currently, all PETSc parallel matrix formats are partitioned by
429 contiguous chunks of rows across the processors.
430 - Each processor needs to insert only elements that it owns
431 locally (but any non-local elements will be sent to the
432 appropriate processor during matrix assembly).
433 - Here, we set all entries for a particular row at once.
434 - We can set matrix entries either using either
435 MatSetValuesLocal() or MatSetValues(), as discussed above.
436 */
437 for (k = zs; k < zs + zm; k++) {
438 for (j = ys; j < ys + ym; j++) {
439 for (i = xs; i < xs + xm; i++) {
440 row.k = k;
441 row.j = j;
442 row.i = i;
443 /* boundary points */
444 if (i == 0 || j == 0 || k == 0 || i == Mx - 1 || j == My - 1 || k == Mz - 1) {
445 v[0] = 1.0;
446 PetscCall(MatSetValuesStencil(jac, 1, &row, 1, &row, v, INSERT_VALUES));
447 } else {
448 /* interior grid points */
449 v[0] = -hxhydhz;
450 col[0].k = k - 1;
451 col[0].j = j;
452 col[0].i = i;
453 v[1] = -hxhzdhy;
454 col[1].k = k;
455 col[1].j = j - 1;
456 col[1].i = i;
457 v[2] = -hyhzdhx;
458 col[2].k = k;
459 col[2].j = j;
460 col[2].i = i - 1;
461 v[3] = 2.0 * (hyhzdhx + hxhzdhy + hxhydhz) - sc * PetscExpScalar(x[k][j][i]);
462 col[3].k = row.k;
463 col[3].j = row.j;
464 col[3].i = row.i;
465 v[4] = -hyhzdhx;
466 col[4].k = k;
467 col[4].j = j;
468 col[4].i = i + 1;
469 v[5] = -hxhzdhy;
470 col[5].k = k;
471 col[5].j = j + 1;
472 col[5].i = i;
473 v[6] = -hxhydhz;
474 col[6].k = k + 1;
475 col[6].j = j;
476 col[6].i = i;
477 PetscCall(MatSetValuesStencil(jac, 1, &row, 7, col, v, INSERT_VALUES));
478 }
479 }
480 }
481 }
482 PetscCall(DMDAVecRestoreArrayRead(da, localX, &x));
483 PetscCall(DMRestoreLocalVector(da, &localX));
484
485 /*
486 Assemble matrix, using the 2-step process:
487 MatAssemblyBegin(), MatAssemblyEnd().
488 */
489 PetscCall(MatAssemblyBegin(jac, MAT_FINAL_ASSEMBLY));
490 PetscCall(MatAssemblyEnd(jac, MAT_FINAL_ASSEMBLY));
491
492 /*
493 Normally since the matrix has already been assembled above; this
494 would do nothing. But in the matrix-free mode -snes_mf_operator
495 this tells the "matrix-free" matrix that a new linear system solve
496 is about to be done.
497 */
498
499 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
500 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
501
502 /*
503 Tell the matrix we will never add a new nonzero location to the
504 matrix. If we do, it will generate an error.
505 */
506 PetscCall(MatSetOption(jac, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_TRUE));
507 PetscFunctionReturn(PETSC_SUCCESS);
508 }
509
510 /*TEST
511
512 test:
513 nsize: 4
514 args: -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
515
516 test:
517 suffix: 2
518 nsize: 4
519 args: -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
520
521 test:
522 suffix: 3
523 nsize: 4
524 args: -fdcoloring -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
525
526 test:
527 suffix: 3_ds
528 nsize: 4
529 args: -fdcoloring -fdcoloring_ds -snes_monitor_short -ksp_gmres_cgs_refinement_type refine_always
530
531 test:
532 suffix: 4
533 nsize: 4
534 args: -fdcoloring_local -fdcoloring -ksp_monitor_short -da_refine 1
535 requires: !single
536
537 test:
538 suffix: 5
539 nsize: 4
540 args: -fdcoloring_local -fdcoloring -ksp_monitor_short -da_refine 1 -snes_type newtontrdc
541 requires: !single
542
543 test:
544 suffix: 6
545 nsize: 4
546 args: -fdcoloring_local -fdcoloring -da_refine 1 -snes_type newtontr -snes_tr_fallback_type dogleg
547 requires: !single
548
549 TEST*/
550