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: Option 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: MatAssemblyBegin(jac, MAT_FINAL_ASSEMBLY);
44: MatAssemblyEnd(jac, MAT_FINAL_ASSEMBLY);
45: return 0;
46: }
48: PETSC_EXTERN PetscErrorCode MatAssemblyEnd_MFFD(Mat, MatAssemblyType);
49: PETSC_EXTERN PetscErrorCode MatMFFDSetBase_MFFD(Mat, Vec, Vec);
51: /*@
52: MatSNESMFGetSNES - returns the `SNES` associated with a matrix created with `MatCreateSNESMF()`
54: Not collective
56: Input Parameter:
57: . J - the matrix
59: Output Parameter:
60: . snes - the `SNES` object
62: Level: advanced
64: .seealso: `Mat`, `SNES`, `MatCreateSNESMF()`
65: @*/
66: PetscErrorCode MatSNESMFGetSNES(Mat J, SNES *snes)
67: {
68: MatMFFD j;
70: MatShellGetContext(J, &j);
71: *snes = (SNES)j->ctx;
72: return 0;
73: }
75: /*
76: MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
77: base from the SNES context
79: */
80: static PetscErrorCode MatAssemblyEnd_SNESMF(Mat J, MatAssemblyType mt)
81: {
82: MatMFFD j;
83: SNES snes;
84: Vec u, f;
85: DM dm;
86: DMSNES dms;
88: MatShellGetContext(J, &j);
89: snes = (SNES)j->ctx;
90: MatAssemblyEnd_MFFD(J, mt);
92: SNESGetSolution(snes, &u);
93: SNESGetDM(snes, &dm);
94: DMGetDMSNES(dm, &dms);
95: if ((j->func == (PetscErrorCode(*)(void *, Vec, Vec))SNESComputeFunction) && !dms->ops->computemffunction) {
96: SNESGetFunction(snes, &f, NULL, NULL);
97: MatMFFDSetBase_MFFD(J, u, f);
98: } else {
99: /* f value known by SNES is not correct for other differencing function */
100: MatMFFDSetBase_MFFD(J, u, NULL);
101: }
102: return 0;
103: }
105: /*
106: MatAssemblyEnd_SNESMF_UseBase - Calls MatAssemblyEnd_MFFD() and then sets the
107: base from the SNES context. This version will cause the base to be used for differencing
108: even if the func is not SNESComputeFunction. See: MatSNESMFUseBase()
110: */
111: static PetscErrorCode MatAssemblyEnd_SNESMF_UseBase(Mat J, MatAssemblyType mt)
112: {
113: MatMFFD j;
114: SNES snes;
115: Vec u, f;
117: MatAssemblyEnd_MFFD(J, mt);
118: MatShellGetContext(J, &j);
119: snes = (SNES)j->ctx;
120: SNESGetSolution(snes, &u);
121: SNESGetFunction(snes, &f, NULL, NULL);
122: MatMFFDSetBase_MFFD(J, u, f);
123: return 0;
124: }
126: /*
127: This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
128: uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
129: */
130: static PetscErrorCode MatMFFDSetBase_SNESMF(Mat J, Vec U, Vec F)
131: {
132: MatMFFDSetBase_MFFD(J, U, F);
133: J->ops->assemblyend = MatAssemblyEnd_MFFD;
134: return 0;
135: }
137: static PetscErrorCode MatSNESMFSetReuseBase_SNESMF(Mat J, PetscBool use)
138: {
139: if (use) {
140: J->ops->assemblyend = MatAssemblyEnd_SNESMF_UseBase;
141: } else {
142: J->ops->assemblyend = MatAssemblyEnd_SNESMF;
143: }
144: return 0;
145: }
147: /*@
148: MatSNESMFSetReuseBase - Causes the base vector to be used for differencing even if the function provided to `SNESSetFunction()` is not the
149: same as that provided to `MatMFFDSetFunction()`.
151: Logically Collective
153: Input Parameters:
154: + J - the `MATMFFD` matrix
155: - use - if true always reuse the base vector instead of recomputing f(u) even if the function in the `MATMFFD` is
156: not `SNESComputeFunction()`
158: Note:
159: Care must be taken when using this routine to insure that the function provided to `MatMFFDSetFunction()`, call it F_MF() is compatible with
160: 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
162: Developer Note:
163: This was provided for the MOOSE team who desired to have a `SNESSetFunction()` function that could change configurations (similar to variable
164: switching) to contacts while the function provided to `MatMFFDSetFunction()` cannot. Except for the possibility of changing the configuration
165: both functions compute the same mathematical function so the differencing makes sense.
167: Level: advanced
169: .seealso: `MATMFFD`, `MatMFFDSetFunction()`, `SNESSetFunction()`, `MatCreateSNESMF()`, `MatSNESMFGetReuseBase()`
170: @*/
171: PetscErrorCode MatSNESMFSetReuseBase(Mat J, PetscBool use)
172: {
174: PetscTryMethod(J, "MatSNESMFSetReuseBase_C", (Mat, PetscBool), (J, use));
175: return 0;
176: }
178: static PetscErrorCode MatSNESMFGetReuseBase_SNESMF(Mat J, PetscBool *use)
179: {
180: if (J->ops->assemblyend == MatAssemblyEnd_SNESMF_UseBase) *use = PETSC_TRUE;
181: else *use = PETSC_FALSE;
182: return 0;
183: }
185: /*@
186: MatSNESMFGetReuseBase - Determines if the base vector is to be used for differencing even if the function provided to `SNESSetFunction()` is not the
187: same as that provided to `MatMFFDSetFunction()`.
189: Logically Collective
191: Input Parameter:
192: . J - the `MATMFFD` matrix
194: Output Parameter:
195: . use - if true always reuse the base vector instead of recomputing f(u) even if the function in the `MATMFFD` is
196: not `SNESComputeFunction()`
198: Level: advanced
200: Note:
201: See `MatSNESMFSetReuseBase()`
203: .seealso: `MatSNESMFSetReuseBase()`, `MatCreateSNESMF()`, `MatSNESMFSetReuseBase()`
204: @*/
205: PetscErrorCode MatSNESMFGetReuseBase(Mat J, PetscBool *use)
206: {
208: PetscUseMethod(J, "MatSNESMFGetReuseBase_C", (Mat, PetscBool *), (J, use));
209: return 0;
210: }
212: /*@
213: MatCreateSNESMF - Creates a matrix-free matrix context for use with
214: a `SNES` solver. This matrix can be used as the Jacobian argument for
215: the routine `SNESSetJacobian()`. See `MatCreateMFFD()` for details on how
216: the finite difference computation is done.
218: Collective
220: Input Parameters:
221: . snes - the `SNES` context
223: Output Parameter:
224: . J - the matrix-free matrix which is of type `MATMFFD`
226: Level: advanced
228: Notes:
229: You can call `SNESSetJacobian()` with `MatMFFDComputeJacobian()` if you are using matrix and not a different
230: preconditioner matrix
232: If you wish to provide a different function to do differencing on to compute the matrix-free operator than
233: that provided to `SNESSetFunction()` then call `MatMFFDSetFunction()` with your function after this call.
235: The difference between this routine and `MatCreateMFFD()` is that this matrix
236: automatically gets the current base vector from the `SNES` object and not from an
237: explicit call to `MatMFFDSetBase()`.
239: If `MatMFFDSetBase()` is ever called on jac then this routine will NO longer get
240: the x from the `SNES` object and `MatMFFDSetBase()` must from that point on be used to
241: change the base vector x.
243: Using a different function for the differencing will not work if you are using non-linear left preconditioning.
245: .seealso: `MATMFFD, `MatDestroy()`, `MatMFFDSetFunction()`, `MatMFFDSetFunctionError()`, `MatMFFDDSSetUmin()`
246: `MatMFFDSetHHistory()`, `MatMFFDResetHHistory()`, `MatCreateMFFD()`,
247: `MatMFFDGetH()`, `MatMFFDRegister()`, `MatMFFDComputeJacobian()`, `MatSNESMFSetReuseBase()`, `MatSNESMFGetReuseBase()`
248: @*/
249: PetscErrorCode MatCreateSNESMF(SNES snes, Mat *J)
250: {
251: PetscInt n, N;
252: MatMFFD mf;
254: if (snes->vec_func) {
255: VecGetLocalSize(snes->vec_func, &n);
256: VecGetSize(snes->vec_func, &N);
257: } else if (snes->dm) {
258: Vec tmp;
259: DMGetGlobalVector(snes->dm, &tmp);
260: VecGetLocalSize(tmp, &n);
261: VecGetSize(tmp, &N);
262: DMRestoreGlobalVector(snes->dm, &tmp);
263: } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must call SNESSetFunction() or SNESSetDM() first");
264: MatCreateMFFD(PetscObjectComm((PetscObject)snes), n, n, N, N, J);
265: MatShellGetContext(*J, &mf);
266: mf->ctx = snes;
268: if (snes->npc && snes->npcside == PC_LEFT) {
269: MatMFFDSetFunction(*J, (PetscErrorCode(*)(void *, Vec, Vec))SNESComputeFunctionDefaultNPC, snes);
270: } else {
271: DM dm;
272: DMSNES dms;
274: SNESGetDM(snes, &dm);
275: DMGetDMSNES(dm, &dms);
276: MatMFFDSetFunction(*J, (PetscErrorCode(*)(void *, Vec, Vec))(dms->ops->computemffunction ? SNESComputeMFFunction : SNESComputeFunction), snes);
277: }
278: (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
280: PetscObjectComposeFunction((PetscObject)*J, "MatMFFDSetBase_C", MatMFFDSetBase_SNESMF);
281: PetscObjectComposeFunction((PetscObject)*J, "MatSNESMFSetReuseBase_C", MatSNESMFSetReuseBase_SNESMF);
282: PetscObjectComposeFunction((PetscObject)*J, "MatSNESMFGetReuseBase_C", MatSNESMFGetReuseBase_SNESMF);
283: return 0;
284: }