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: Logically Collective on SNES
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
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: Level: developer
24: Warning:
25: If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
26: the x from the SNES object and MatMFFDSetBase() must from that point on be used to
27: change the base vector x.
29: Notes:
30: This can be passed into SNESSetJacobian() as the Jacobian evaluation function argument
31: when using a completely matrix-free solver,
32: that is the B matrix is also the same matrix operator. This is used when you select
33: -snes_mf but rarely used directly by users. (All this routine does is call MatAssemblyBegin/End() on
34: the Mat jac.)
36: .seealso: MatMFFDGetH(), MatCreateSNESMF(), MatCreateMFFD(), MATMFFD,
37: MatMFFDSetHHistory(), MatMFFDSetFunctionError(), MatCreateMFFD(), SNESSetJacobian()
39: @*/
40: PetscErrorCode MatMFFDComputeJacobian(SNES snes,Vec x,Mat jac,Mat B,void *dummy)
41: {
45: MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
46: MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
47: return(0);
48: }
50: PETSC_EXTERN PetscErrorCode MatAssemblyEnd_MFFD(Mat,MatAssemblyType);
51: PETSC_EXTERN PetscErrorCode MatMFFDSetBase_MFFD(Mat,Vec,Vec);
53: /*@
54: MatSNESMFGetSNES - returns the SNES associated with a matrix created with MatCreateSNESMF()
56: Not collective
58: Input Parameter:
59: . J - the matrix
61: Output Parameter:
62: . snes - the SNES object
64: Level: advanced
66: .seealso: MatCreateSNESMF()
67: @*/
68: PetscErrorCode MatSNESMFGetSNES(Mat J,SNES *snes)
69: {
70: MatMFFD j;
74: MatShellGetContext(J,&j);
75: *snes = (SNES)j->ctx;
76: return(0);
77: }
79: /*
80: MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
81: base from the SNES context
83: */
84: static PetscErrorCode MatAssemblyEnd_SNESMF(Mat J,MatAssemblyType mt)
85: {
87: MatMFFD j;
88: SNES snes;
89: Vec u,f;
90: DM dm;
91: DMSNES dms;
94: MatShellGetContext(J,&j);
95: snes = (SNES)j->ctx;
96: MatAssemblyEnd_MFFD(J,mt);
98: SNESGetSolution(snes,&u);
99: SNESGetDM(snes,&dm);
100: DMGetDMSNES(dm,&dms);
101: if ((j->func == (PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction) && !dms->ops->computemffunction) {
102: SNESGetFunction(snes,&f,NULL,NULL);
103: MatMFFDSetBase_MFFD(J,u,f);
104: } else {
105: /* f value known by SNES is not correct for other differencing function */
106: MatMFFDSetBase_MFFD(J,u,NULL);
107: }
108: return(0);
109: }
111: /*
112: MatAssemblyEnd_SNESMF_UseBase - Calls MatAssemblyEnd_MFFD() and then sets the
113: base from the SNES context. This version will cause the base to be used for differencing
114: even if the func is not SNESComputeFunction. See: MatSNESMFUseBase()
116: */
117: static PetscErrorCode MatAssemblyEnd_SNESMF_UseBase(Mat J,MatAssemblyType mt)
118: {
120: MatMFFD j;
121: SNES snes;
122: Vec u,f;
125: MatAssemblyEnd_MFFD(J,mt);
126: MatShellGetContext(J,&j);
127: snes = (SNES)j->ctx;
128: SNESGetSolution(snes,&u);
129: SNESGetFunction(snes,&f,NULL,NULL);
130: MatMFFDSetBase_MFFD(J,u,f);
131: return(0);
132: }
134: /*
135: This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
136: uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
137: */
138: static PetscErrorCode MatMFFDSetBase_SNESMF(Mat J,Vec U,Vec F)
139: {
143: MatMFFDSetBase_MFFD(J,U,F);
144: J->ops->assemblyend = MatAssemblyEnd_MFFD;
145: return(0);
146: }
148: static PetscErrorCode MatSNESMFSetReuseBase_SNESMF(Mat J,PetscBool use)
149: {
151: if (use) {
152: J->ops->assemblyend = MatAssemblyEnd_SNESMF_UseBase;
153: } else {
154: J->ops->assemblyend = MatAssemblyEnd_SNESMF;
155: }
156: return(0);
157: }
159: /*@
160: MatSNESMFSetReuseBase - Causes the base vector to be used for differencing even if the function provided to SNESSetFunction() is not the
161: same as that provided to MatMFFDSetFunction().
163: Logically Collective on Mat
165: Input Parameters:
166: + J - the MatMFFD matrix
167: - use - if true always reuse the base vector instead of recomputing f(u) even if the function in the MatSNESMF is
168: not SNESComputeFunction()
170: Notes:
171: Care must be taken when using this routine to insure that the function provided to MatMFFDSetFunction(), call it F_MF() is compatible with
172: 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
174: Developer Notes:
175: This was provided for the MOOSE team who desired to have a SNESSetFunction() function that could change configurations due
176: to contacts while the function provided to MatMFFDSetFunction() cannot. Except for the possibility of changing the configuration
177: both functions compute the same mathematical function so the differencing makes sense.
179: Level: advanced
181: .seealso: MatCreateSNESMF(), MatSNESMFGetReuseBase()
182: @*/
183: PetscErrorCode MatSNESMFSetReuseBase(Mat J,PetscBool use)
184: {
189: PetscTryMethod(J,"MatSNESMFSetReuseBase_C",(Mat,PetscBool),(J,use));
190: return(0);
191: }
193: static PetscErrorCode MatSNESMFGetReuseBase_SNESMF(Mat J,PetscBool *use)
194: {
196: if (J->ops->assemblyend == MatAssemblyEnd_SNESMF_UseBase) *use = PETSC_TRUE;
197: else *use = PETSC_FALSE;
198: return(0);
199: }
201: /*@
202: MatSNESMFGetReuseBase - Determines if the base vector is to be used for differencing even if the function provided to SNESSetFunction() is not the
203: same as that provided to MatMFFDSetFunction().
205: Logically Collective on Mat
207: Input Parameter:
208: . J - the MatMFFD matrix
210: Output Parameter:
211: . use - if true always reuse the base vector instead of recomputing f(u) even if the function in the MatSNESMF is
212: not SNESComputeFunction()
214: Notes:
215: Care must be taken when using this routine to insure that the function provided to MatMFFDSetFunction(), call it F_MF() is compatible with
216: 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
218: Developer Notes:
219: This was provided for the MOOSE team who desired to have a SNESSetFunction() function that could change configurations due
220: to contacts while the function provided to MatMFFDSetFunction() cannot. Except for the possibility of changing the configuration
221: both functions compute the same mathematical function so the differencing makes sense.
223: Level: advanced
225: .seealso: MatCreateSNESMF(), MatSNESMFSetReuseBase()
226: @*/
227: PetscErrorCode MatSNESMFGetReuseBase(Mat J,PetscBool *use)
228: {
233: PetscUseMethod(J,"MatSNESMFGetReuseBase_C",(Mat,PetscBool*),(J,use));
234: return(0);
235: }
237: /*@
238: MatCreateSNESMF - Creates a matrix-free matrix context for use with
239: a SNES solver. This matrix can be used as the Jacobian argument for
240: the routine SNESSetJacobian(). See MatCreateMFFD() for details on how
241: the finite difference computation is done.
243: Collective on SNES
245: Input Parameters:
246: . snes - the SNES context
248: Output Parameter:
249: . J - the matrix-free matrix
251: Level: advanced
253: Notes:
254: You can call SNESSetJacobian() with MatMFFDComputeJacobian() if you are using matrix and not a different
255: preconditioner matrix
257: If you wish to provide a different function to do differencing on to compute the matrix-free operator than
258: that provided to SNESSetFunction() then call MatMFFDSetFunction() with your function after this call.
260: The difference between this routine and MatCreateMFFD() is that this matrix
261: automatically gets the current base vector from the SNES object and not from an
262: explicit call to MatMFFDSetBase().
264: Warning:
265: If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
266: the x from the SNES object and MatMFFDSetBase() must from that point on be used to
267: change the base vector x.
269: Warning:
270: Using a different function for the differencing will not work if you are using non-linear left preconditioning.
272: .seealso: MatDestroy(), MatMFFDSetFunction(), MatMFFDSetFunctionError(), MatMFFDDSSetUmin()
273: MatMFFDSetHHistory(), MatMFFDResetHHistory(), MatCreateMFFD(),
274: MatMFFDGetH(), MatMFFDRegister(), MatMFFDComputeJacobian(), MatSNESMFSetReuseBase(), MatSNESMFGetReuseBase()
276: @*/
277: PetscErrorCode MatCreateSNESMF(SNES snes,Mat *J)
278: {
280: PetscInt n,N;
281: MatMFFD mf;
284: if (snes->vec_func) {
285: VecGetLocalSize(snes->vec_func,&n);
286: VecGetSize(snes->vec_func,&N);
287: } else if (snes->dm) {
288: Vec tmp;
289: DMGetGlobalVector(snes->dm,&tmp);
290: VecGetLocalSize(tmp,&n);
291: VecGetSize(tmp,&N);
292: DMRestoreGlobalVector(snes->dm,&tmp);
293: } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() or SNESSetDM() first");
294: MatCreateMFFD(PetscObjectComm((PetscObject)snes),n,n,N,N,J);
295: MatShellGetContext(*J,&mf);
296: mf->ctx = snes;
298: if (snes->npc && snes->npcside== PC_LEFT) {
299: MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunctionDefaultNPC,snes);
300: } else {
301: DM dm;
302: DMSNES dms;
304: SNESGetDM(snes,&dm);
305: DMGetDMSNES(dm,&dms);
306: MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))(dms->ops->computemffunction ? SNESComputeMFFunction : SNESComputeFunction),snes);
307: }
308: (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
310: PetscObjectComposeFunction((PetscObject)*J,"MatMFFDSetBase_C",MatMFFDSetBase_SNESMF);
311: PetscObjectComposeFunction((PetscObject)*J,"MatSNESMFSetReuseBase_C",MatSNESMFSetReuseBase_SNESMF);
312: PetscObjectComposeFunction((PetscObject)*J,"MatSNESMFGetReuseBase_C",MatSNESMFGetReuseBase_SNESMF);
313: return(0);
314: }