1 static char help[] = "Solves the nonlinear system, the Bratu (SFI - solid fuel ignition) problem in a 2D rectangular domain.\n\
2 This example also illustrates the use of matrix coloring. Runtime options include:\n\
3 -par <parameter>, where <parameter> indicates the problem's nonlinearity\n\
4 problem SFI: <parameter> = Bratu parameter (0 <= par <= 6.81)\n\
5 -mx <xg>, where <xg> = number of grid points in the x-direction\n\
6 -my <yg>, where <yg> = number of grid points in the y-direction\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.
18
19 A finite difference approximation with the usual 5-point stencil
20 is used to discretize the boundary value problem to obtain a nonlinear
21 system of equations.
22
23 The parallel version of this code is snes/tutorials/ex5.c
24
25 */
26
27 /*
28 Include "petscsnes.h" so that we can use SNES solvers. Note that
29 this file automatically includes:
30 petscsys.h - base PETSc routines petscvec.h - vectors
31 petscmat.h - matrices
32 petscis.h - index sets petscksp.h - Krylov subspace methods
33 petscviewer.h - viewers petscpc.h - preconditioners
34 petscksp.h - linear solvers
35 */
36
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 PetscInt mx; /* Discretization in x-direction */
47 PetscInt my; /* Discretization in y-direction */
48 } AppCtx;
49
50 /*
51 User-defined routines
52 */
53 static PetscErrorCode FormJacobian(SNES, Vec, Mat, Mat, void *);
54 static PetscErrorCode FormFunction(SNES, Vec, Vec, void *);
55 static PetscErrorCode FormInitialGuess(AppCtx *, Vec);
56 static PetscErrorCode ConvergenceTest(KSP, PetscInt, PetscReal, KSPConvergedReason *, void *);
57 static PetscErrorCode ConvergenceDestroy(PetscCtxRt);
58 static PetscErrorCode postcheck(SNES, Vec, Vec, Vec, PetscBool *, PetscBool *, void *);
59 static PetscErrorCode monitor_change_deltamax(SNES, PetscInt, PetscReal, void *);
60
main(int argc,char ** argv)61 int main(int argc, char **argv)
62 {
63 SNES snes; /* nonlinear solver context */
64 Vec x, r; /* solution, residual vectors */
65 Mat J; /* Jacobian matrix */
66 AppCtx user; /* user-defined application context */
67 PetscInt i, its, N, hist_its[50];
68 PetscMPIInt size;
69 PetscReal bratu_lambda_max = 6.81, bratu_lambda_min = 0., history[50];
70 MatFDColoring fdcoloring;
71 PetscBool matrix_free = PETSC_FALSE, flg, fd_coloring = PETSC_FALSE, use_convergence_test = PETSC_FALSE, pc = PETSC_FALSE, prunejacobian = PETSC_FALSE, null_appctx = PETSC_TRUE, test_tr_deltamax = PETSC_FALSE;
72 KSP ksp;
73 PetscInt *testarray;
74
75 PetscFunctionBeginUser;
76 PetscCall(PetscInitialize(&argc, &argv, NULL, help));
77 PetscCallMPI(MPI_Comm_size(PETSC_COMM_WORLD, &size));
78 PetscCheck(size == 1, PETSC_COMM_WORLD, PETSC_ERR_WRONG_MPI_SIZE, "This is a uniprocessor example only!");
79
80 /*
81 Initialize problem parameters
82 */
83 user.mx = 4;
84 user.my = 4;
85 user.param = 6.0;
86 PetscCall(PetscOptionsGetInt(NULL, NULL, "-mx", &user.mx, NULL));
87 PetscCall(PetscOptionsGetInt(NULL, NULL, "-my", &user.my, NULL));
88 PetscCall(PetscOptionsGetReal(NULL, NULL, "-par", &user.param, NULL));
89 PetscCall(PetscOptionsGetBool(NULL, NULL, "-pc", &pc, NULL));
90 PetscCheck(user.param < bratu_lambda_max && user.param > bratu_lambda_min, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Lambda is out of range");
91 N = user.mx * user.my;
92 PetscCall(PetscOptionsGetBool(NULL, NULL, "-use_convergence_test", &use_convergence_test, NULL));
93 PetscCall(PetscOptionsGetBool(NULL, NULL, "-prune_jacobian", &prunejacobian, NULL));
94 PetscCall(PetscOptionsGetBool(NULL, NULL, "-null_appctx", &null_appctx, NULL));
95 PetscCall(PetscOptionsGetBool(NULL, NULL, "-test_tr_deltamax", &test_tr_deltamax, NULL));
96
97 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
98 Create nonlinear solver context
99 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
100
101 PetscCall(SNESCreate(PETSC_COMM_WORLD, &snes));
102 if (pc) {
103 PetscCall(SNESSetType(snes, SNESNEWTONTR));
104 PetscCall(SNESNewtonTRSetPostCheck(snes, postcheck, NULL));
105 }
106 if (test_tr_deltamax) {
107 PetscCall(SNESSetType(snes, SNESNEWTONTR));
108 PetscCall(SNESMonitorSet(snes, monitor_change_deltamax, NULL, NULL));
109 }
110
111 /* Test application context handling from Python */
112 if (!null_appctx) PetscCall(SNESSetApplicationContext(snes, (void *)&user));
113
114 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
115 Create vector data structures; set function evaluation routine
116 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
117
118 PetscCall(VecCreate(PETSC_COMM_WORLD, &x));
119 PetscCall(VecSetSizes(x, PETSC_DECIDE, N));
120 PetscCall(VecSetFromOptions(x));
121 PetscCall(VecDuplicate(x, &r));
122
123 /*
124 Set function evaluation routine and vector. Whenever the nonlinear
125 solver needs to evaluate the nonlinear function, it will call this
126 routine.
127 - Note that the final routine argument is the user-defined
128 context that provides application-specific data for the
129 function evaluation routine.
130 */
131 PetscCall(SNESSetFunction(snes, r, FormFunction, (void *)&user));
132
133 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
134 Create matrix data structure; set Jacobian evaluation routine
135 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
136
137 /*
138 Create matrix. Here we only approximately preallocate storage space
139 for the Jacobian. See the users manual for a discussion of better
140 techniques for preallocating matrix memory.
141 */
142 PetscCall(PetscOptionsGetBool(NULL, NULL, "-snes_mf", &matrix_free, NULL));
143 if (!matrix_free) {
144 PetscBool matrix_free_operator = PETSC_FALSE;
145 PetscCall(PetscOptionsGetBool(NULL, NULL, "-snes_mf_operator", &matrix_free_operator, NULL));
146 if (matrix_free_operator) matrix_free = PETSC_FALSE;
147 }
148 if (!matrix_free) PetscCall(MatCreateSeqAIJ(PETSC_COMM_WORLD, N, N, 5, NULL, &J));
149
150 /*
151 This option will cause the Jacobian to be computed via finite differences
152 efficiently using a coloring of the columns of the matrix.
153 */
154 PetscCall(PetscOptionsGetBool(NULL, NULL, "-snes_fd_coloring", &fd_coloring, NULL));
155 PetscCheck(!matrix_free || !fd_coloring, PETSC_COMM_WORLD, PETSC_ERR_ARG_INCOMP, "Use only one of -snes_mf, -snes_fd_coloring options! You can do -snes_mf_operator -snes_fd_coloring");
156
157 if (fd_coloring) {
158 ISColoring iscoloring;
159 MatColoring mc;
160 if (prunejacobian) {
161 /* Initialize x with random nonzero values so that the nonzeros in the Jacobian
162 can better reflect the sparsity structure of the Jacobian. */
163 PetscRandom rctx;
164 PetscCall(PetscRandomCreate(PETSC_COMM_WORLD, &rctx));
165 PetscCall(PetscRandomSetInterval(rctx, 1.0, 2.0));
166 PetscCall(VecSetRandom(x, rctx));
167 PetscCall(PetscRandomDestroy(&rctx));
168 }
169
170 /*
171 This initializes the nonzero structure of the Jacobian. This is artificial
172 because clearly if we had a routine to compute the Jacobian we won't need
173 to use finite differences.
174 */
175 PetscCall(FormJacobian(snes, x, J, J, &user));
176
177 /*
178 Color the matrix, i.e. determine groups of columns that share no common
179 rows. These columns in the Jacobian can all be computed simultaneously.
180 */
181 PetscCall(MatColoringCreate(J, &mc));
182 PetscCall(MatColoringSetType(mc, MATCOLORINGSL));
183 PetscCall(MatColoringSetFromOptions(mc));
184 PetscCall(MatColoringApply(mc, &iscoloring));
185 PetscCall(MatColoringDestroy(&mc));
186 /*
187 Create the data structure that SNESComputeJacobianDefaultColor() uses
188 to compute the actual Jacobians via finite differences.
189 */
190 PetscCall(MatFDColoringCreate(J, iscoloring, &fdcoloring));
191 PetscCall(MatFDColoringSetFunction(fdcoloring, (MatFDColoringFn *)FormFunction, &user));
192 PetscCall(MatFDColoringSetFromOptions(fdcoloring));
193 PetscCall(MatFDColoringSetUp(J, iscoloring, fdcoloring));
194 /*
195 Tell SNES to use the routine SNESComputeJacobianDefaultColor()
196 to compute Jacobians.
197 */
198 PetscCall(SNESSetJacobian(snes, J, J, SNESComputeJacobianDefaultColor, fdcoloring));
199 PetscCall(ISColoringDestroy(&iscoloring));
200 if (prunejacobian) PetscCall(SNESPruneJacobianColor(snes, J, J));
201 }
202 /*
203 Set Jacobian matrix data structure and default Jacobian evaluation
204 routine. Whenever the nonlinear solver needs to compute the
205 Jacobian matrix, it will call this routine.
206 - Note that the final routine argument is the user-defined
207 context that provides application-specific data for the
208 Jacobian evaluation routine.
209 - The user can override with:
210 -snes_fd : default finite differencing approximation of Jacobian
211 -snes_mf : matrix-free Newton-Krylov method with no preconditioning
212 (unless user explicitly sets preconditioner)
213 -snes_mf_operator : form matrix used to construct the preconditioner as set by the user,
214 but use matrix-free approx for Jacobian-vector
215 products within Newton-Krylov method
216 */
217 else if (!matrix_free) {
218 PetscCall(SNESSetJacobian(snes, J, J, FormJacobian, (void *)&user));
219 }
220
221 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
222 Customize nonlinear solver; set runtime options
223 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
224
225 /*
226 Set runtime options (e.g., -snes_monitor -snes_rtol <rtol> -ksp_type <type>)
227 */
228 PetscCall(SNESSetFromOptions(snes));
229
230 /*
231 Set array that saves the function norms. This array is intended
232 when the user wants to save the convergence history for later use
233 rather than just to view the function norms via -snes_monitor.
234 */
235 PetscCall(SNESSetConvergenceHistory(snes, history, hist_its, 50, PETSC_TRUE));
236
237 /*
238 Add a user provided convergence test; this is to test that SNESNEWTONTR properly calls the
239 user provided test before the specialized test. The convergence context is just an array to
240 test that it gets properly freed at the end
241 */
242 if (use_convergence_test) {
243 PetscCall(SNESGetKSP(snes, &ksp));
244 PetscCall(PetscMalloc1(5, &testarray));
245 PetscCall(KSPSetConvergenceTest(ksp, ConvergenceTest, testarray, ConvergenceDestroy));
246 }
247
248 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
249 Evaluate initial guess; then solve nonlinear system
250 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
251 /*
252 Note: The user should initialize the vector, x, with the initial guess
253 for the nonlinear solver prior to calling SNESSolve(). In particular,
254 to employ an initial guess of zero, the user should explicitly set
255 this vector to zero by calling VecSet().
256 */
257 PetscCall(FormInitialGuess(&user, x));
258 PetscCall(SNESSolve(snes, NULL, x));
259 PetscCall(SNESGetIterationNumber(snes, &its));
260 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Number of SNES iterations = %" PetscInt_FMT "\n", its));
261
262 /*
263 Print the convergence history. This is intended just to demonstrate
264 use of the data attained via SNESSetConvergenceHistory().
265 */
266 PetscCall(PetscOptionsHasName(NULL, NULL, "-print_history", &flg));
267 if (flg) {
268 for (i = 0; i < its + 1; i++) PetscCall(PetscPrintf(PETSC_COMM_WORLD, "iteration %" PetscInt_FMT ": Linear iterations %" PetscInt_FMT " Function norm = %g\n", i, hist_its[i], (double)history[i]));
269 }
270
271 /* Test NewtonTR API */
272 PetscCall(SNESNewtonTRSetTolerances(snes, 1.0, 2.0, 3.0));
273 PetscCall(SNESNewtonTRSetUpdateParameters(snes, 4.0, 5.0, 6.0, 7.0, 8.0));
274 PetscCall(PetscObjectTypeCompare((PetscObject)snes, SNESNEWTONTR, &flg));
275 if (flg) {
276 PetscReal tmp[5];
277
278 PetscCall(SNESNewtonTRGetTolerances(snes, &tmp[0], &tmp[1], &tmp[2]));
279 PetscCheck(tmp[0] == 1.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
280 PetscCheck(tmp[1] == 2.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
281 PetscCheck(tmp[2] == 3.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
282 PetscCall(SNESNewtonTRGetUpdateParameters(snes, &tmp[0], &tmp[1], &tmp[2], &tmp[3], &tmp[4]));
283 PetscCheck(tmp[0] == 4.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
284 PetscCheck(tmp[1] == 5.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
285 PetscCheck(tmp[2] == 6.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
286 PetscCheck(tmp[3] == 7.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
287 PetscCheck(tmp[4] == 8.0, PETSC_COMM_WORLD, PETSC_ERR_PLIB, "Wrong value");
288 }
289
290 /* - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
291 Free work space. All PETSc objects should be destroyed when they
292 are no longer needed.
293 - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - */
294
295 if (!matrix_free) PetscCall(MatDestroy(&J));
296 if (fd_coloring) PetscCall(MatFDColoringDestroy(&fdcoloring));
297 PetscCall(VecDestroy(&x));
298 PetscCall(VecDestroy(&r));
299 PetscCall(SNESDestroy(&snes));
300 PetscCall(PetscFinalize());
301 return 0;
302 }
303
304 /*
305 FormInitialGuess - Forms initial approximation.
306
307 Input Parameters:
308 user - user-defined application context
309 X - vector
310
311 Output Parameter:
312 X - vector
313 */
FormInitialGuess(AppCtx * user,Vec X)314 PetscErrorCode FormInitialGuess(AppCtx *user, Vec X)
315 {
316 PetscInt i, j, row, mx, my;
317 PetscReal lambda, temp1, temp, hx, hy;
318 PetscScalar *x;
319
320 PetscFunctionBeginUser;
321 mx = user->mx;
322 my = user->my;
323 lambda = user->param;
324
325 hx = 1.0 / (PetscReal)(mx - 1);
326 hy = 1.0 / (PetscReal)(my - 1);
327
328 /*
329 Get a pointer to vector data.
330 - For default PETSc vectors, VecGetArray() returns a pointer to
331 the data array. Otherwise, the routine is implementation dependent.
332 - You MUST call VecRestoreArray() when you no longer need access to
333 the array.
334 */
335 PetscCall(VecGetArray(X, &x));
336 temp1 = lambda / (lambda + 1.0);
337 for (j = 0; j < my; j++) {
338 temp = (PetscReal)(PetscMin(j, my - j - 1)) * hy;
339 for (i = 0; i < mx; i++) {
340 row = i + j * mx;
341 if (i == 0 || j == 0 || i == mx - 1 || j == my - 1) {
342 x[row] = 0.0;
343 continue;
344 }
345 x[row] = temp1 * PetscSqrtReal(PetscMin((PetscReal)(PetscMin(i, mx - i - 1)) * hx, temp));
346 }
347 }
348
349 /*
350 Restore vector
351 */
352 PetscCall(VecRestoreArray(X, &x));
353 PetscFunctionReturn(PETSC_SUCCESS);
354 }
355
356 /*
357 FormFunction - Evaluates nonlinear function, F(x).
358
359 Input Parameters:
360 . snes - the SNES context
361 . X - input vector
362 . ptr - optional user-defined context, as set by SNESSetFunction()
363
364 Output Parameter:
365 . F - function vector
366 */
FormFunction(SNES snes,Vec X,Vec F,void * ptr)367 PetscErrorCode FormFunction(SNES snes, Vec X, Vec F, void *ptr)
368 {
369 AppCtx *user = (AppCtx *)ptr;
370 PetscInt i, j, row, mx, my;
371 PetscReal two = 2.0, one = 1.0, lambda, hx, hy, hxdhy, hydhx;
372 PetscScalar ut, ub, ul, ur, u, uxx, uyy, sc, *f;
373 const PetscScalar *x;
374
375 PetscFunctionBeginUser;
376 mx = user->mx;
377 my = user->my;
378 lambda = user->param;
379 hx = one / (PetscReal)(mx - 1);
380 hy = one / (PetscReal)(my - 1);
381 sc = hx * hy;
382 hxdhy = hx / hy;
383 hydhx = hy / hx;
384
385 /*
386 Get pointers to vector data
387 */
388 PetscCall(VecGetArrayRead(X, &x));
389 PetscCall(VecGetArray(F, &f));
390
391 /*
392 Compute function
393 */
394 for (j = 0; j < my; j++) {
395 for (i = 0; i < mx; i++) {
396 row = i + j * mx;
397 if (i == 0 || j == 0 || i == mx - 1 || j == my - 1) {
398 f[row] = x[row];
399 continue;
400 }
401 u = x[row];
402 ub = x[row - mx];
403 ul = x[row - 1];
404 ut = x[row + mx];
405 ur = x[row + 1];
406 uxx = (-ur + two * u - ul) * hydhx;
407 uyy = (-ut + two * u - ub) * hxdhy;
408 f[row] = uxx + uyy - sc * lambda * PetscExpScalar(u);
409 }
410 }
411
412 /*
413 Restore vectors
414 */
415 PetscCall(VecRestoreArrayRead(X, &x));
416 PetscCall(VecRestoreArray(F, &f));
417 PetscFunctionReturn(PETSC_SUCCESS);
418 }
419
420 /*
421 FormJacobian - Evaluates Jacobian matrix.
422
423 Input Parameters:
424 . snes - the SNES context
425 . x - input vector
426 . ptr - optional user-defined context, as set by SNESSetJacobian()
427
428 Output Parameters:
429 . A - Jacobian matrix
430 . B - optionally different matrix used to construct the preconditioner
431
432 */
FormJacobian(SNES snes,Vec X,Mat J,Mat jac,void * ptr)433 PetscErrorCode FormJacobian(SNES snes, Vec X, Mat J, Mat jac, void *ptr)
434 {
435 AppCtx *user = (AppCtx *)ptr; /* user-defined application context */
436 PetscInt i, j, row, mx, my, col[5];
437 PetscScalar two = 2.0, one = 1.0, lambda, v[5], sc;
438 const PetscScalar *x;
439 PetscReal hx, hy, hxdhy, hydhx;
440
441 PetscFunctionBeginUser;
442 mx = user->mx;
443 my = user->my;
444 lambda = user->param;
445 hx = 1.0 / (PetscReal)(mx - 1);
446 hy = 1.0 / (PetscReal)(my - 1);
447 sc = hx * hy;
448 hxdhy = hx / hy;
449 hydhx = hy / hx;
450
451 /*
452 Get pointer to vector data
453 */
454 PetscCall(VecGetArrayRead(X, &x));
455
456 /*
457 Compute entries of the Jacobian
458 */
459 for (j = 0; j < my; j++) {
460 for (i = 0; i < mx; i++) {
461 row = i + j * mx;
462 if (i == 0 || j == 0 || i == mx - 1 || j == my - 1) {
463 PetscCall(MatSetValues(jac, 1, &row, 1, &row, &one, INSERT_VALUES));
464 continue;
465 }
466 v[0] = -hxdhy;
467 col[0] = row - mx;
468 v[1] = -hydhx;
469 col[1] = row - 1;
470 v[2] = two * (hydhx + hxdhy) - sc * lambda * PetscExpScalar(x[row]);
471 col[2] = row;
472 v[3] = -hydhx;
473 col[3] = row + 1;
474 v[4] = -hxdhy;
475 col[4] = row + mx;
476 PetscCall(MatSetValues(jac, 1, &row, 5, col, v, INSERT_VALUES));
477 }
478 }
479
480 /*
481 Restore vector
482 */
483 PetscCall(VecRestoreArrayRead(X, &x));
484
485 /*
486 Assemble matrix
487 */
488 PetscCall(MatAssemblyBegin(jac, MAT_FINAL_ASSEMBLY));
489 PetscCall(MatAssemblyEnd(jac, MAT_FINAL_ASSEMBLY));
490
491 if (jac != J) {
492 PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
493 PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
494 }
495 PetscFunctionReturn(PETSC_SUCCESS);
496 }
497
ConvergenceTest(KSP ksp,PetscInt it,PetscReal nrm,KSPConvergedReason * reason,PetscCtx ctx)498 PetscErrorCode ConvergenceTest(KSP ksp, PetscInt it, PetscReal nrm, KSPConvergedReason *reason, PetscCtx ctx)
499 {
500 PetscFunctionBeginUser;
501 *reason = KSP_CONVERGED_ITERATING;
502 if (it > 1) {
503 *reason = KSP_CONVERGED_ITS;
504 PetscCall(PetscInfo(NULL, "User provided convergence test returning after 2 iterations\n"));
505 }
506 PetscFunctionReturn(PETSC_SUCCESS);
507 }
508
ConvergenceDestroy(PetscCtxRt ctx)509 PetscErrorCode ConvergenceDestroy(PetscCtxRt ctx)
510 {
511 PetscFunctionBeginUser;
512 PetscCall(PetscInfo(NULL, "User provided convergence destroy called\n"));
513 PetscCall(PetscFree(*(void **)ctx));
514 PetscFunctionReturn(PETSC_SUCCESS);
515 }
516
postcheck(SNES snes,Vec x,Vec y,Vec w,PetscBool * changed_y,PetscBool * changed_w,PetscCtx ctx)517 PetscErrorCode postcheck(SNES snes, Vec x, Vec y, Vec w, PetscBool *changed_y, PetscBool *changed_w, PetscCtx ctx)
518 {
519 PetscReal norm;
520 Vec tmp;
521
522 PetscFunctionBeginUser;
523 PetscCall(VecDuplicate(x, &tmp));
524 PetscCall(VecWAXPY(tmp, -1.0, x, w));
525 PetscCall(VecNorm(tmp, NORM_2, &norm));
526 PetscCall(VecDestroy(&tmp));
527 PetscCall(PetscPrintf(PETSC_COMM_WORLD, "Norm of search step %g\n", (double)norm));
528 PetscFunctionReturn(PETSC_SUCCESS);
529 }
530
monitor_change_deltamax(SNES snes,PetscInt it,PetscReal fnorm,PetscCtx ctx)531 PetscErrorCode monitor_change_deltamax(SNES snes, PetscInt it, PetscReal fnorm, PetscCtx ctx)
532 {
533 PetscFunctionBeginUser;
534 if (it == 0) PetscCall(SNESNewtonTRSetTolerances(snes, PETSC_CURRENT, 0.01, PETSC_CURRENT));
535 PetscFunctionReturn(PETSC_SUCCESS);
536 }
537
538 /*TEST
539
540 build:
541 requires: !single
542
543 test:
544 args: -ksp_gmres_cgs_refinement_type refine_always
545
546 test:
547 suffix: 2
548 args: -snes_monitor_short -snes_type newtontr -ksp_gmres_cgs_refinement_type refine_always
549
550 test:
551 suffix: 2_trdeltamax_change
552 args: -snes_monitor_short -snes_type newtontr -ksp_gmres_cgs_refinement_type refine_always -pc -test_tr_deltamax
553
554 test:
555 suffix: 2a
556 filter: grep -i KSPConvergedDefault > /dev/null && echo "Found KSPConvergedDefault"
557 args: -snes_monitor_short -snes_type newtontr -ksp_gmres_cgs_refinement_type refine_always -info
558 requires: defined(PETSC_USE_INFO)
559
560 test:
561 suffix: 2b
562 filter: grep -i "User provided convergence test" > /dev/null && echo "Found User provided convergence test"
563 args: -snes_monitor_short -snes_type newtontr -ksp_gmres_cgs_refinement_type refine_always -use_convergence_test -info
564 requires: defined(PETSC_USE_INFO)
565
566 test:
567 suffix: 2c
568 args: -snes_converged_reason -snes_type newtontr -snes_tr_qn {{same different}separate output} -pc_type mat -snes_view -snes_tr_qn_mat_type lmvmdfp -snes_tr_norm_type infinity
569
570 test:
571 suffix: 3
572 args: -snes_monitor_short -mat_coloring_type sl -snes_fd_coloring -mx 8 -my 11 -ksp_gmres_cgs_refinement_type refine_always
573
574 test:
575 suffix: 4
576 args: -pc -par 6.807 -snes_monitor -snes_converged_reason
577
578 test:
579 suffix: 5
580 args: -snes_monitor_short -mat_coloring_type sl -snes_fd_coloring -mx 8 -my 11 -ksp_gmres_cgs_refinement_type refine_always -prune_jacobian
581 output_file: output/ex1_3.out
582
583 test:
584 suffix: 6
585 args: -snes_monitor draw:image:testfile -viewer_view
586
587 test:
588 suffix: python
589 requires: petsc4py
590 args: -python -snes_type python -snes_python_type ex1.py:MySNES -snes_view -null_appctx {{0 1}separate output}
591 localrunfiles: ex1.py
592
593 TEST*/
594