Actual source code: snesj2.c
1: #include <petsc/private/snesimpl.h>
2: #include <petscdm.h>
4: /*
5: MatFDColoringSetFunction() takes a function with four arguments, we want to use SNESComputeFunction()
6: since it logs function computation information.
7: */
8: static PetscErrorCode SNESComputeFunctionCtx(void *snes, Vec x, Vec f, void *ctx)
9: {
10: return SNESComputeFunction((SNES)snes, x, f);
11: }
12: static PetscErrorCode SNESComputeMFFunctionCtx(void *snes, Vec x, Vec f, void *ctx)
13: {
14: return SNESComputeMFFunction((SNES)snes, x, f);
15: }
17: /*@C
18: SNESComputeJacobianDefaultColor - Computes the Jacobian using
19: finite differences and coloring to exploit matrix sparsity.
21: Collective
23: Input Parameters:
24: + snes - nonlinear solver object
25: . x1 - location at which to evaluate Jacobian
26: - ctx - `MatFDColoring` context or `NULL`
28: Output Parameters:
29: + J - Jacobian matrix (not altered in this routine)
30: - B - newly computed Jacobian matrix to use with preconditioner (generally the same as `J`)
32: Options Database Keys:
33: + -snes_fd_color_use_mat - use a matrix coloring from the explicit matrix nonzero pattern instead of from the `DM` providing the matrix
34: . -snes_fd_color - Activates `SNESComputeJacobianDefaultColor()` in `SNESSetFromOptions()`
35: . -mat_fd_coloring_err <err> - Sets <err> (square root of relative error in the function)
36: . -mat_fd_coloring_umin <umin> - Sets umin, the minimum allowable u-value magnitude
37: . -mat_fd_type - Either wp or ds (see `MATMFFD_WP` or `MATMFFD_DS`)
38: . -snes_mf_operator - Use matrix-free application of Jacobian
39: - -snes_mf - Use matrix-free Jacobian with no explicit Jacobian representation
41: Notes:
43: Level: intermediate
45: If the coloring is not provided through the context, this will first try to get the
46: coloring from the `DM`. If the `DM` has no coloring routine, then it will try to
47: get the coloring from the matrix. This requires that the matrix have its nonzero locations already provided.
49: `SNES` supports three approaches for computing (approximate) Jacobians: user provided via `SNESSetJacobian()`, matrix-free via `SNESSetUseMatrixFree()`,
50: and computing explicitly with finite differences and coloring using `MatFDColoring`. It is also possible to use automatic differentiation
51: and the `MatFDColoring` object, see src/ts/tutorials/autodiff/ex16adj_tl.cxx
53: This function can be provided to `SNESSetJacobian()` along with an appropriate sparse matrix to hold the Jacobian
55: Developer Note:
56: The function has a poorly chosen name since it does not mention the use of finite differences
58: .seealso: [](ch_snes), `SNES`, `SNESSetJacobian()`, `SNESTestJacobian()`, `SNESComputeJacobianDefault()`, `SNESSetUseMatrixFree()`,
59: `MatFDColoringCreate()`, `MatFDColoringSetFunction()`
60: @*/
61: PetscErrorCode SNESComputeJacobianDefaultColor(SNES snes, Vec x1, Mat J, Mat B, void *ctx)
62: {
63: MatFDColoring color = (MatFDColoring)ctx;
64: DM dm;
65: MatColoring mc;
66: ISColoring iscoloring;
67: PetscBool hascolor;
68: PetscBool solvec, matcolor = PETSC_FALSE;
69: DMSNES dms;
71: PetscFunctionBegin;
73: if (!color) PetscCall(PetscObjectQuery((PetscObject)B, "SNESMatFDColoring", (PetscObject *)&color));
75: if (!color) {
76: PetscCall(SNESGetDM(snes, &dm));
77: PetscCall(DMHasColoring(dm, &hascolor));
78: matcolor = PETSC_FALSE;
79: PetscCall(PetscOptionsGetBool(((PetscObject)snes)->options, ((PetscObject)snes)->prefix, "-snes_fd_color_use_mat", &matcolor, NULL));
80: if (hascolor && !matcolor) {
81: PetscCall(DMCreateColoring(dm, IS_COLORING_GLOBAL, &iscoloring));
82: } else {
83: PetscCall(MatColoringCreate(B, &mc));
84: PetscCall(MatColoringSetDistance(mc, 2));
85: PetscCall(MatColoringSetType(mc, MATCOLORINGSL));
86: PetscCall(MatColoringSetFromOptions(mc));
87: PetscCall(MatColoringApply(mc, &iscoloring));
88: PetscCall(MatColoringDestroy(&mc));
89: }
90: PetscCall(MatFDColoringCreate(B, iscoloring, &color));
91: PetscCall(DMGetDMSNES(dm, &dms));
92: if (dms->ops->computemffunction) {
93: PetscCall(MatFDColoringSetFunction(color, (PetscErrorCode(*)(void))SNESComputeMFFunctionCtx, NULL));
94: } else {
95: PetscCall(MatFDColoringSetFunction(color, (PetscErrorCode(*)(void))SNESComputeFunctionCtx, NULL));
96: }
97: PetscCall(MatFDColoringSetFromOptions(color));
98: PetscCall(MatFDColoringSetUp(B, iscoloring, color));
99: PetscCall(ISColoringDestroy(&iscoloring));
100: PetscCall(PetscObjectCompose((PetscObject)B, "SNESMatFDColoring", (PetscObject)color));
101: PetscCall(PetscObjectDereference((PetscObject)color));
102: }
104: /* F is only usable if there is no RHS on the SNES and the full solution corresponds to x1 */
105: PetscCall(VecEqual(x1, snes->vec_sol, &solvec));
106: if (!snes->vec_rhs && solvec) {
107: Vec F;
108: PetscCall(SNESGetFunction(snes, &F, NULL, NULL));
109: PetscCall(MatFDColoringSetF(color, F));
110: }
111: PetscCall(MatFDColoringApply(B, color, x1, snes));
112: if (J != B) {
113: PetscCall(MatAssemblyBegin(J, MAT_FINAL_ASSEMBLY));
114: PetscCall(MatAssemblyEnd(J, MAT_FINAL_ASSEMBLY));
115: }
116: PetscFunctionReturn(PETSC_SUCCESS);
117: }
119: /*@C
120: SNESPruneJacobianColor - Remove nondiagonal zeros in the Jacobian matrix and update the `MatMFFD` coloring information based on the new nonzero structure
122: Collective
124: Input Parameters:
125: + snes - the `SNES` context
126: . J - Jacobian matrix (not altered in this routine)
127: - B - newly computed Jacobian matrix to use with preconditioner (generally the same as `J`)
129: Level: intermediate
131: Notes:
132: This function improves the `MatMFFD` coloring performance when the Jacobian matrix is overallocated or contains
133: many constant zeros entries, which is typically the case when the matrix is generated by a `DM`
134: and multiple fields are involved.
136: Users need to make sure that the Jacobian matrix is properly filled to reflect the sparsity
137: structure. For `MatMFFD` coloring, the values of nonzero entries are not important. So one can
138: usually call `SNESComputeJacobian()` with randomized input vectors to generate a dummy Jacobian.
139: `SNESComputeJacobian()` should be called before `SNESSolve()` but after `SNESSetUp()`.
141: .seealso: [](ch_snes), `SNESComputeJacobianDefaultColor()`, `MatEliminateZeros()`, `MatFDColoringCreate()`, `MatFDColoringSetFunction()`
142: @*/
143: PetscErrorCode SNESPruneJacobianColor(SNES snes, Mat J, Mat B)
144: {
145: DM dm;
146: DMSNES dms;
147: MatColoring mc;
148: ISColoring iscoloring;
149: MatFDColoring matfdcoloring;
151: PetscFunctionBegin;
152: /* Generate new coloring after eliminating zeros in the matrix */
153: PetscCall(MatEliminateZeros(B, PETSC_TRUE));
154: PetscCall(MatColoringCreate(B, &mc));
155: PetscCall(MatColoringSetDistance(mc, 2));
156: PetscCall(MatColoringSetType(mc, MATCOLORINGSL));
157: PetscCall(MatColoringSetFromOptions(mc));
158: PetscCall(MatColoringApply(mc, &iscoloring));
159: PetscCall(MatColoringDestroy(&mc));
160: /* Replace the old coloring with the new one */
161: PetscCall(MatFDColoringCreate(B, iscoloring, &matfdcoloring));
162: PetscCall(SNESGetDM(snes, &dm));
163: PetscCall(DMGetDMSNES(dm, &dms));
164: // Comment out the following branch to bypass the coverage test. You can uncomment it when needed.
165: //if (dms->ops->computemffunction) {
166: // PetscCall(MatFDColoringSetFunction(matfdcoloring, (PetscErrorCode(*)(void))SNESComputeMFFunctionCtx, NULL));
167: //} else {
168: PetscCall(MatFDColoringSetFunction(matfdcoloring, (PetscErrorCode(*)(void))SNESComputeFunctionCtx, NULL));
169: //}
170: PetscCall(MatFDColoringSetFromOptions(matfdcoloring));
171: PetscCall(MatFDColoringSetUp(B, iscoloring, matfdcoloring));
172: PetscCall(PetscObjectCompose((PetscObject)B, "SNESMatFDColoring", (PetscObject)matfdcoloring));
173: PetscCall(PetscObjectDereference((PetscObject)matfdcoloring));
174: PetscCall(ISColoringDestroy(&iscoloring));
175: PetscFunctionReturn(PETSC_SUCCESS);
176: }