Actual source code: snesmfj.c
2: #include <petsc/private/snesimpl.h>
3: #include <petscdm.h>
4: #include <../src/mat/impls/mffd/mffdimpl.h>
5: #include <petsc/private/matimpl.h>
7: /*@C
8: MatMFFDComputeJacobian - Tells the matrix-free Jacobian object the new location at which
9: Jacobian matrix vector products will be computed at, i.e. J(x) * a. The x is obtained
10: from the `SNES` object (using `SNESGetSolution()`).
12: Collective
14: Input Parameters:
15: + snes - the nonlinear solver context
16: . x - the point at which the Jacobian vector products will be performed
17: . jac - the matrix-free Jacobian object of `MatType` `MATMFFD`, likely obtained with `MatCreateSNESMF()`
18: . B - either the same as jac or another matrix type (ignored)
19: . flag - not relevant for matrix-free form
20: - dummy - the user context (ignored)
22: Options Database Key:
23: . -snes_mf - use the matrix created with `MatSNESMFCreate()` to setup the Jacobian for each new solution in the Newton process
25: Level: developer
27: Notes:
28: If `MatMFFDSetBase()` is ever called on jac then this routine will NO longer get
29: the x from the `SNES` object and `MatMFFDSetBase()` must from that point on be used to
30: change the base vector x.
32: This can be passed into `SNESSetJacobian()` as the Jacobian evaluation function argument
33: when using a completely matrix-free solver,
34: that is the B matrix is also the same matrix operator. This is used when you select
35: -snes_mf but rarely used directly by users. (All this routine does is call `MatAssemblyBegin/End()` on
36: the `Mat` jac.)
38: .seealso: `MatMFFDGetH()`, `MatCreateSNESMF()`, `MatMFFDSetBase()`, `MatCreateMFFD()`, `MATMFFD`,
39: `MatMFFDSetHHistory()`, `MatMFFDSetFunctionError()`, `MatCreateMFFD()`, `SNESSetJacobian()`
40: @*/
41: PetscErrorCode MatMFFDComputeJacobian(SNES snes, Vec x, Mat jac, Mat B, void *dummy)
42: {
43: PetscFunctionBegin;
44: PetscCall(MatAssemblyBegin(jac, MAT_FINAL_ASSEMBLY));
45: PetscCall(MatAssemblyEnd(jac, MAT_FINAL_ASSEMBLY));
46: PetscFunctionReturn(PETSC_SUCCESS);
47: }
49: PETSC_EXTERN PetscErrorCode MatAssemblyEnd_MFFD(Mat, MatAssemblyType);
50: PETSC_EXTERN PetscErrorCode MatMFFDSetBase_MFFD(Mat, Vec, Vec);
52: /*@
53: MatSNESMFGetSNES - returns the `SNES` associated with a matrix created with `MatCreateSNESMF()`
55: Not Collective
57: Input Parameter:
58: . J - the matrix
60: Output Parameter:
61: . snes - the `SNES` object
63: Level: advanced
65: .seealso: `Mat`, `SNES`, `MatCreateSNESMF()`
66: @*/
67: PetscErrorCode MatSNESMFGetSNES(Mat J, SNES *snes)
68: {
69: MatMFFD j;
71: PetscFunctionBegin;
72: PetscCall(MatShellGetContext(J, &j));
73: *snes = (SNES)j->ctx;
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
77: /*
78: MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
79: base from the SNES context
81: */
82: static PetscErrorCode MatAssemblyEnd_SNESMF(Mat J, MatAssemblyType mt)
83: {
84: MatMFFD j;
85: SNES snes;
86: Vec u, f;
87: DM dm;
88: DMSNES dms;
90: PetscFunctionBegin;
91: PetscCall(MatShellGetContext(J, &j));
92: snes = (SNES)j->ctx;
93: PetscCall(MatAssemblyEnd_MFFD(J, mt));
95: PetscCall(SNESGetSolution(snes, &u));
96: PetscCall(SNESGetDM(snes, &dm));
97: PetscCall(DMGetDMSNES(dm, &dms));
98: if ((j->func == (PetscErrorCode(*)(void *, Vec, Vec))SNESComputeFunction) && !dms->ops->computemffunction) {
99: PetscCall(SNESGetFunction(snes, &f, NULL, NULL));
100: PetscCall(MatMFFDSetBase_MFFD(J, u, f));
101: } else {
102: /* f value known by SNES is not correct for other differencing function */
103: PetscCall(MatMFFDSetBase_MFFD(J, u, NULL));
104: }
105: PetscFunctionReturn(PETSC_SUCCESS);
106: }
108: /*
109: MatAssemblyEnd_SNESMF_UseBase - Calls MatAssemblyEnd_MFFD() and then sets the
110: base from the SNES context. This version will cause the base to be used for differencing
111: even if the func is not SNESComputeFunction. See: MatSNESMFUseBase()
113: */
114: static PetscErrorCode MatAssemblyEnd_SNESMF_UseBase(Mat J, MatAssemblyType mt)
115: {
116: MatMFFD j;
117: SNES snes;
118: Vec u, f;
120: PetscFunctionBegin;
121: PetscCall(MatAssemblyEnd_MFFD(J, mt));
122: PetscCall(MatShellGetContext(J, &j));
123: snes = (SNES)j->ctx;
124: PetscCall(SNESGetSolution(snes, &u));
125: PetscCall(SNESGetFunction(snes, &f, NULL, NULL));
126: PetscCall(MatMFFDSetBase_MFFD(J, u, f));
127: PetscFunctionReturn(PETSC_SUCCESS);
128: }
130: /*
131: This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
132: uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
133: */
134: static PetscErrorCode MatMFFDSetBase_SNESMF(Mat J, Vec U, Vec F)
135: {
136: PetscFunctionBegin;
137: PetscCall(MatMFFDSetBase_MFFD(J, U, F));
138: J->ops->assemblyend = MatAssemblyEnd_MFFD;
139: PetscFunctionReturn(PETSC_SUCCESS);
140: }
142: static PetscErrorCode MatSNESMFSetReuseBase_SNESMF(Mat J, PetscBool use)
143: {
144: PetscFunctionBegin;
145: if (use) {
146: J->ops->assemblyend = MatAssemblyEnd_SNESMF_UseBase;
147: } else {
148: J->ops->assemblyend = MatAssemblyEnd_SNESMF;
149: }
150: PetscFunctionReturn(PETSC_SUCCESS);
151: }
153: /*@
154: MatSNESMFSetReuseBase - Causes the base vector to be used for differencing even if the function provided to `SNESSetFunction()` is not the
155: same as that provided to `MatMFFDSetFunction()`.
157: Logically Collective
159: Input Parameters:
160: + J - the `MATMFFD` matrix
161: - use - if true always reuse the base vector instead of recomputing f(u) even if the function in the `MATMFFD` is
162: not `SNESComputeFunction()`
164: Level: advanced
166: Note:
167: Care must be taken when using this routine to insure that the function provided to `MatMFFDSetFunction()`, call it F_MF() is compatible with
168: with that provided to `SNESSetFunction()`, call it F_SNES(). That is, (F_MF(u + h*d) - F_SNES(u))/h has to approximate the derivative
170: Developer Note:
171: This was provided for the MOOSE team who desired to have a `SNESSetFunction()` function that could change configurations (similar to variable
172: switching) to contacts while the function provided to `MatMFFDSetFunction()` cannot. Except for the possibility of changing the configuration
173: both functions compute the same mathematical function so the differencing makes sense.
175: .seealso: `MATMFFD`, `MatMFFDSetFunction()`, `SNESSetFunction()`, `MatCreateSNESMF()`, `MatSNESMFGetReuseBase()`
176: @*/
177: PetscErrorCode MatSNESMFSetReuseBase(Mat J, PetscBool use)
178: {
179: PetscFunctionBegin;
181: PetscTryMethod(J, "MatSNESMFSetReuseBase_C", (Mat, PetscBool), (J, use));
182: PetscFunctionReturn(PETSC_SUCCESS);
183: }
185: static PetscErrorCode MatSNESMFGetReuseBase_SNESMF(Mat J, PetscBool *use)
186: {
187: PetscFunctionBegin;
188: if (J->ops->assemblyend == MatAssemblyEnd_SNESMF_UseBase) *use = PETSC_TRUE;
189: else *use = PETSC_FALSE;
190: PetscFunctionReturn(PETSC_SUCCESS);
191: }
193: /*@
194: MatSNESMFGetReuseBase - Determines if the base vector is to be used for differencing even if the function provided to `SNESSetFunction()` is not the
195: same as that provided to `MatMFFDSetFunction()`.
197: Logically Collective
199: Input Parameter:
200: . J - the `MATMFFD` matrix
202: Output Parameter:
203: . use - if true always reuse the base vector instead of recomputing f(u) even if the function in the `MATMFFD` is
204: not `SNESComputeFunction()`
206: Level: advanced
208: Note:
209: See `MatSNESMFSetReuseBase()`
211: .seealso: `MatSNESMFSetReuseBase()`, `MatCreateSNESMF()`, `MatSNESMFSetReuseBase()`
212: @*/
213: PetscErrorCode MatSNESMFGetReuseBase(Mat J, PetscBool *use)
214: {
215: PetscFunctionBegin;
217: PetscUseMethod(J, "MatSNESMFGetReuseBase_C", (Mat, PetscBool *), (J, use));
218: PetscFunctionReturn(PETSC_SUCCESS);
219: }
221: /*@
222: MatCreateSNESMF - Creates a matrix-free matrix context for use with
223: a `SNES` solver. This matrix can be used as the Jacobian argument for
224: the routine `SNESSetJacobian()`. See `MatCreateMFFD()` for details on how
225: the finite difference computation is done.
227: Collective
229: Input Parameters:
230: . snes - the `SNES` context
232: Output Parameter:
233: . J - the matrix-free matrix which is of type `MATMFFD`
235: Level: advanced
237: Notes:
238: You can call `SNESSetJacobian()` with `MatMFFDComputeJacobian()` if you are using matrix and not a different
239: preconditioner matrix
241: If you wish to provide a different function to do differencing on to compute the matrix-free operator than
242: that provided to `SNESSetFunction()` then call `MatMFFDSetFunction()` with your function after this call.
244: The difference between this routine and `MatCreateMFFD()` is that this matrix
245: automatically gets the current base vector from the `SNES` object and not from an
246: explicit call to `MatMFFDSetBase()`.
248: If `MatMFFDSetBase()` is ever called on jac then this routine will NO longer get
249: the x from the `SNES` object and `MatMFFDSetBase()` must from that point on be used to
250: change the base vector x.
252: Using a different function for the differencing will not work if you are using non-linear left preconditioning.
254: .seealso: `MATMFFD, `MatDestroy()`, `MatMFFDSetFunction()`, `MatMFFDSetFunctionError()`, `MatMFFDDSSetUmin()`
255: `MatMFFDSetHHistory()`, `MatMFFDResetHHistory()`, `MatCreateMFFD()`,
256: `MatMFFDGetH()`, `MatMFFDRegister()`, `MatMFFDComputeJacobian()`, `MatSNESMFSetReuseBase()`, `MatSNESMFGetReuseBase()`
257: @*/
258: PetscErrorCode MatCreateSNESMF(SNES snes, Mat *J)
259: {
260: PetscInt n, N;
261: MatMFFD mf;
263: PetscFunctionBegin;
264: if (snes->vec_func) {
265: PetscCall(VecGetLocalSize(snes->vec_func, &n));
266: PetscCall(VecGetSize(snes->vec_func, &N));
267: } else if (snes->dm) {
268: Vec tmp;
269: PetscCall(DMGetGlobalVector(snes->dm, &tmp));
270: PetscCall(VecGetLocalSize(tmp, &n));
271: PetscCall(VecGetSize(tmp, &N));
272: PetscCall(DMRestoreGlobalVector(snes->dm, &tmp));
273: } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetFunction() or SNESSetDM() first");
274: PetscCall(MatCreateMFFD(PetscObjectComm((PetscObject)snes), n, n, N, N, J));
275: PetscCall(MatShellGetContext(*J, &mf));
276: mf->ctx = snes;
278: if (snes->npc && snes->npcside == PC_LEFT) {
279: PetscCall(MatMFFDSetFunction(*J, (PetscErrorCode(*)(void *, Vec, Vec))SNESComputeFunctionDefaultNPC, snes));
280: } else {
281: DM dm;
282: DMSNES dms;
284: PetscCall(SNESGetDM(snes, &dm));
285: PetscCall(DMGetDMSNES(dm, &dms));
286: PetscCall(MatMFFDSetFunction(*J, (PetscErrorCode(*)(void *, Vec, Vec))(dms->ops->computemffunction ? SNESComputeMFFunction : SNESComputeFunction), snes));
287: }
288: (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
290: PetscCall(PetscObjectComposeFunction((PetscObject)*J, "MatMFFDSetBase_C", MatMFFDSetBase_SNESMF));
291: PetscCall(PetscObjectComposeFunction((PetscObject)*J, "MatSNESMFSetReuseBase_C", MatSNESMFSetReuseBase_SNESMF));
292: PetscCall(PetscObjectComposeFunction((PetscObject)*J, "MatSNESMFGetReuseBase_C", MatSNESMFGetReuseBase_SNESMF));
293: PetscFunctionReturn(PETSC_SUCCESS);
294: }