xref: /petsc/src/snes/interface/snesj2.c (revision a69119a591a03a9d906b29c0a4e9802e4d7c9795)
1 
2 #include <petsc/private/snesimpl.h> /*I  "petscsnes.h"  I*/
3 #include <petscdm.h>                /*I  "petscdm.h"    I*/
4 
5 /*
6    MatFDColoringSetFunction() takes a function with four arguments, we want to use SNESComputeFunction()
7    since it logs function computation information.
8 */
9 static PetscErrorCode SNESComputeFunctionCtx(SNES snes, Vec x, Vec f, void *ctx) {
10   return SNESComputeFunction(snes, x, f);
11 }
12 static PetscErrorCode SNESComputeMFFunctionCtx(SNES snes, Vec x, Vec f, void *ctx) {
13   return SNESComputeMFFunction(snes, x, f);
14 }
15 
16 /*@C
17     SNESComputeJacobianDefaultColor - Computes the Jacobian using
18     finite differences and coloring to exploit matrix sparsity.
19 
20     Collective on SNES
21 
22     Input Parameters:
23 +   snes - nonlinear solver object
24 .   x1 - location at which to evaluate Jacobian
25 -   ctx - MatFDColoring context or NULL
26 
27     Output Parameters:
28 +   J - Jacobian matrix (not altered in this routine)
29 -   B - newly computed Jacobian matrix to use with preconditioner (generally the same as J)
30 
31     Level: intermediate
32 
33    Options Database Key:
34 +  -snes_fd_color_use_mat - use a matrix coloring from the explicit matrix nonzero pattern instead of from the DM providing the matrix
35 .  -snes_fd_color - Activates SNESComputeJacobianDefaultColor() in SNESSetFromOptions()
36 .  -mat_fd_coloring_err <err> - Sets <err> (square root of relative error in the function)
37 .  -mat_fd_coloring_umin <umin> - Sets umin, the minimum allowable u-value magnitude
38 .  -mat_fd_type - Either wp or ds (see MATMFFD_WP or MATMFFD_DS)
39 .  -snes_mf_operator - Use matrix free application of Jacobian
40 -  -snes_mf - Use matrix free Jacobian with no explicit Jacobian representation
41 
42     Notes:
43         If the coloring is not provided through the context, this will first try to get the
44         coloring from the DM.  If the DM type has no coloring routine, then it will try to
45         get the coloring from the matrix.  This requires that the matrix have nonzero entries
46         precomputed.
47 
48        SNES supports three approaches for computing (approximate) Jacobians: user provided via SNESSetJacobian(), matrix free via SNESSetUseMatrixFree(),
49        and computing explicitly with finite differences and coloring using MatFDColoring. It is also possible to use automatic differentiation
50        and the MatFDColoring object, see src/ts/tutorials/autodiff/ex16adj_tl.cxx
51 
52 .seealso: `SNESSetJacobian()`, `SNESTestJacobian()`, `SNESComputeJacobianDefault()`, `SNESSetUseMatrixFree()`,
53           `MatFDColoringCreate()`, `MatFDColoringSetFunction()`
54 
55 @*/
56 PetscErrorCode SNESComputeJacobianDefaultColor(SNES snes, Vec x1, Mat J, Mat B, void *ctx) {
57   MatFDColoring color = (MatFDColoring)ctx;
58   DM            dm;
59   MatColoring   mc;
60   ISColoring    iscoloring;
61   PetscBool     hascolor;
62   PetscBool     solvec, matcolor = PETSC_FALSE;
63   DMSNES        dms;
64 
65   PetscFunctionBegin;
66   if (color) PetscValidHeaderSpecific(color, MAT_FDCOLORING_CLASSID, 5);
67   if (!color) PetscCall(PetscObjectQuery((PetscObject)B, "SNESMatFDColoring", (PetscObject *)&color));
68 
69   if (!color) {
70     PetscCall(SNESGetDM(snes, &dm));
71     PetscCall(DMHasColoring(dm, &hascolor));
72     matcolor = PETSC_FALSE;
73     PetscCall(PetscOptionsGetBool(((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_fd_color_use_mat", &matcolor, NULL));
74     if (hascolor && !matcolor) {
75       PetscCall(DMCreateColoring(dm, IS_COLORING_GLOBAL, &iscoloring));
76     } else {
77       PetscCall(MatColoringCreate(B, &mc));
78       PetscCall(MatColoringSetDistance(mc, 2));
79       PetscCall(MatColoringSetType(mc, MATCOLORINGSL));
80       PetscCall(MatColoringSetFromOptions(mc));
81       PetscCall(MatColoringApply(mc, &iscoloring));
82       PetscCall(MatColoringDestroy(&mc));
83     }
84     PetscCall(MatFDColoringCreate(B, iscoloring, &color));
85     PetscCall(DMGetDMSNES(dm, &dms));
86     if (dms->ops->computemffunction) {
87       PetscCall(MatFDColoringSetFunction(color, (PetscErrorCode(*)(void))SNESComputeMFFunctionCtx, NULL));
88     } else {
89       PetscCall(MatFDColoringSetFunction(color, (PetscErrorCode(*)(void))SNESComputeFunctionCtx, NULL));
90     }
91     PetscCall(MatFDColoringSetFromOptions(color));
92     PetscCall(MatFDColoringSetUp(B, iscoloring, color));
93     PetscCall(ISColoringDestroy(&iscoloring));
94     PetscCall(PetscObjectCompose((PetscObject)B, "SNESMatFDColoring", (PetscObject)color));
95     PetscCall(PetscObjectDereference((PetscObject)color));
96   }
97 
98   /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
99   PetscCall(VecEqual(x1, snes->vec_sol, &solvec));
100   if (!snes->vec_rhs && solvec) {
101     Vec F;
102     PetscCall(SNESGetFunction(snes, &F, NULL, NULL));
103     PetscCall(MatFDColoringSetF(color, F));
104   }
105   PetscCall(MatFDColoringApply(B, color, x1, snes));
106   if (J != B) {
107     PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
108     PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
109   }
110   PetscFunctionReturn(0);
111 }
112