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: {
42: MatAssemblyBegin(jac,MAT_FINAL_ASSEMBLY);
43: MatAssemblyEnd(jac,MAT_FINAL_ASSEMBLY);
44: return 0;
45: }
47: PETSC_EXTERN PetscErrorCode MatAssemblyEnd_MFFD(Mat,MatAssemblyType);
48: PETSC_EXTERN PetscErrorCode MatMFFDSetBase_MFFD(Mat,Vec,Vec);
50: /*@
51: MatSNESMFGetSNES - returns the SNES associated with a matrix created with MatCreateSNESMF()
53: Not collective
55: Input Parameter:
56: . J - the matrix
58: Output Parameter:
59: . snes - the SNES object
61: Level: advanced
63: .seealso: MatCreateSNESMF()
64: @*/
65: PetscErrorCode MatSNESMFGetSNES(Mat J,SNES *snes)
66: {
67: MatMFFD j;
69: MatShellGetContext(J,&j);
70: *snes = (SNES)j->ctx;
71: return 0;
72: }
74: /*
75: MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
76: base from the SNES context
78: */
79: static PetscErrorCode MatAssemblyEnd_SNESMF(Mat J,MatAssemblyType mt)
80: {
81: MatMFFD j;
82: SNES snes;
83: Vec u,f;
84: DM dm;
85: DMSNES dms;
87: MatShellGetContext(J,&j);
88: snes = (SNES)j->ctx;
89: MatAssemblyEnd_MFFD(J,mt);
91: SNESGetSolution(snes,&u);
92: SNESGetDM(snes,&dm);
93: DMGetDMSNES(dm,&dms);
94: if ((j->func == (PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction) && !dms->ops->computemffunction) {
95: SNESGetFunction(snes,&f,NULL,NULL);
96: MatMFFDSetBase_MFFD(J,u,f);
97: } else {
98: /* f value known by SNES is not correct for other differencing function */
99: MatMFFDSetBase_MFFD(J,u,NULL);
100: }
101: return 0;
102: }
104: /*
105: MatAssemblyEnd_SNESMF_UseBase - Calls MatAssemblyEnd_MFFD() and then sets the
106: base from the SNES context. This version will cause the base to be used for differencing
107: even if the func is not SNESComputeFunction. See: MatSNESMFUseBase()
109: */
110: static PetscErrorCode MatAssemblyEnd_SNESMF_UseBase(Mat J,MatAssemblyType mt)
111: {
112: MatMFFD j;
113: SNES snes;
114: Vec u,f;
116: MatAssemblyEnd_MFFD(J,mt);
117: MatShellGetContext(J,&j);
118: snes = (SNES)j->ctx;
119: SNESGetSolution(snes,&u);
120: SNESGetFunction(snes,&f,NULL,NULL);
121: MatMFFDSetBase_MFFD(J,u,f);
122: return 0;
123: }
125: /*
126: This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
127: uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
128: */
129: static PetscErrorCode MatMFFDSetBase_SNESMF(Mat J,Vec U,Vec F)
130: {
131: MatMFFDSetBase_MFFD(J,U,F);
132: J->ops->assemblyend = MatAssemblyEnd_MFFD;
133: return 0;
134: }
136: static PetscErrorCode MatSNESMFSetReuseBase_SNESMF(Mat J,PetscBool use)
137: {
138: if (use) {
139: J->ops->assemblyend = MatAssemblyEnd_SNESMF_UseBase;
140: } else {
141: J->ops->assemblyend = MatAssemblyEnd_SNESMF;
142: }
143: return 0;
144: }
146: /*@
147: MatSNESMFSetReuseBase - Causes the base vector to be used for differencing even if the function provided to SNESSetFunction() is not the
148: same as that provided to MatMFFDSetFunction().
150: Logically Collective on Mat
152: Input Parameters:
153: + J - the MatMFFD matrix
154: - use - if true always reuse the base vector instead of recomputing f(u) even if the function in the MatSNESMF is
155: not SNESComputeFunction()
157: Notes:
158: Care must be taken when using this routine to insure that the function provided to MatMFFDSetFunction(), call it F_MF() is compatible with
159: 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
161: Developer Notes:
162: This was provided for the MOOSE team who desired to have a SNESSetFunction() function that could change configurations due
163: to contacts while the function provided to MatMFFDSetFunction() cannot. Except for the possibility of changing the configuration
164: both functions compute the same mathematical function so the differencing makes sense.
166: Level: advanced
168: .seealso: MatCreateSNESMF(), MatSNESMFGetReuseBase()
169: @*/
170: PetscErrorCode MatSNESMFSetReuseBase(Mat J,PetscBool use)
171: {
173: PetscTryMethod(J,"MatSNESMFSetReuseBase_C",(Mat,PetscBool),(J,use));
174: return 0;
175: }
177: static PetscErrorCode MatSNESMFGetReuseBase_SNESMF(Mat J,PetscBool *use)
178: {
179: if (J->ops->assemblyend == MatAssemblyEnd_SNESMF_UseBase) *use = PETSC_TRUE;
180: else *use = PETSC_FALSE;
181: return 0;
182: }
184: /*@
185: MatSNESMFGetReuseBase - Determines if the base vector is to be used for differencing even if the function provided to SNESSetFunction() is not the
186: same as that provided to MatMFFDSetFunction().
188: Logically Collective on Mat
190: Input Parameter:
191: . J - the MatMFFD matrix
193: Output Parameter:
194: . use - if true always reuse the base vector instead of recomputing f(u) even if the function in the MatSNESMF is
195: not SNESComputeFunction()
197: Notes:
198: Care must be taken when using this routine to insure that the function provided to MatMFFDSetFunction(), call it F_MF() is compatible with
199: 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
201: Developer Notes:
202: This was provided for the MOOSE team who desired to have a SNESSetFunction() function that could change configurations due
203: to contacts while the function provided to MatMFFDSetFunction() cannot. Except for the possibility of changing the configuration
204: both functions compute the same mathematical function so the differencing makes sense.
206: Level: advanced
208: .seealso: MatCreateSNESMF(), MatSNESMFSetReuseBase()
209: @*/
210: PetscErrorCode MatSNESMFGetReuseBase(Mat J,PetscBool *use)
211: {
213: PetscUseMethod(J,"MatSNESMFGetReuseBase_C",(Mat,PetscBool*),(J,use));
214: return 0;
215: }
217: /*@
218: MatCreateSNESMF - Creates a matrix-free matrix context for use with
219: a SNES solver. This matrix can be used as the Jacobian argument for
220: the routine SNESSetJacobian(). See MatCreateMFFD() for details on how
221: the finite difference computation is done.
223: Collective on SNES
225: Input Parameters:
226: . snes - the SNES context
228: Output Parameter:
229: . J - the matrix-free matrix
231: Level: advanced
233: Notes:
234: You can call SNESSetJacobian() with MatMFFDComputeJacobian() if you are using matrix and not a different
235: preconditioner matrix
237: If you wish to provide a different function to do differencing on to compute the matrix-free operator than
238: that provided to SNESSetFunction() then call MatMFFDSetFunction() with your function after this call.
240: The difference between this routine and MatCreateMFFD() is that this matrix
241: automatically gets the current base vector from the SNES object and not from an
242: explicit call to MatMFFDSetBase().
244: Warning:
245: If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
246: the x from the SNES object and MatMFFDSetBase() must from that point on be used to
247: change the base vector x.
249: Warning:
250: Using a different function for the differencing will not work if you are using non-linear left preconditioning.
252: .seealso: MatDestroy(), MatMFFDSetFunction(), MatMFFDSetFunctionError(), MatMFFDDSSetUmin()
253: MatMFFDSetHHistory(), MatMFFDResetHHistory(), MatCreateMFFD(),
254: MatMFFDGetH(), MatMFFDRegister(), MatMFFDComputeJacobian(), MatSNESMFSetReuseBase(), MatSNESMFGetReuseBase()
256: @*/
257: PetscErrorCode MatCreateSNESMF(SNES snes,Mat *J)
258: {
259: PetscInt n,N;
260: MatMFFD mf;
262: if (snes->vec_func) {
263: VecGetLocalSize(snes->vec_func,&n);
264: VecGetSize(snes->vec_func,&N);
265: } else if (snes->dm) {
266: Vec tmp;
267: DMGetGlobalVector(snes->dm,&tmp);
268: VecGetLocalSize(tmp,&n);
269: VecGetSize(tmp,&N);
270: DMRestoreGlobalVector(snes->dm,&tmp);
271: } else SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() or SNESSetDM() first");
272: MatCreateMFFD(PetscObjectComm((PetscObject)snes),n,n,N,N,J);
273: MatShellGetContext(*J,&mf);
274: mf->ctx = snes;
276: if (snes->npc && snes->npcside== PC_LEFT) {
277: MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunctionDefaultNPC,snes);
278: } else {
279: DM dm;
280: DMSNES dms;
282: SNESGetDM(snes,&dm);
283: DMGetDMSNES(dm,&dms);
284: MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))(dms->ops->computemffunction ? SNESComputeMFFunction : SNESComputeFunction),snes);
285: }
286: (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
288: PetscObjectComposeFunction((PetscObject)*J,"MatMFFDSetBase_C",MatMFFDSetBase_SNESMF);
289: PetscObjectComposeFunction((PetscObject)*J,"MatSNESMFSetReuseBase_C",MatSNESMFSetReuseBase_SNESMF);
290: PetscObjectComposeFunction((PetscObject)*J,"MatSNESMFGetReuseBase_C",MatSNESMFGetReuseBase_SNESMF);
291: return 0;
292: }