Actual source code: snesmfj.c
1: #define PETSCSNES_DLL
3: #include private/snesimpl.h
4: #include private/matimpl.h
5: #include ../src/mat/impls/mffd/mffdimpl.h
9: /*@C
10: MatMFFDComputeJacobian - Tells the matrix-free Jacobian object the new location at which
11: Jacobian matrix vector products will be computed at, i.e. J(x) * a. The x is obtained
12: from the SNES object (using SNESGetSolution()).
14: Collective on SNES
16: Input Parameters:
17: + snes - the nonlinear solver context
18: . x - the point at which the Jacobian vector products will be performed
19: . jac - the matrix-free Jacobian object
20: . B - either the same as jac or another matrix type (ignored)
21: . flag - not relevent for matrix-free form
22: - dummy - the user context (ignored)
24: Level: developer
26: Warning:
27: If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
28: the x from the SNES object and MatMFFDSetBase() must from that point on be used to
29: change the base vector x.
31: Notes:
32: This can be passed into SNESSetJacobian() when using a completely matrix-free solver,
33: that is the B matrix is also the same matrix operator. This is used when you select
34: -snes_mf but rarely used directly by users. (All this routine does is call MatAssemblyBegin/End() on
35: the Mat jac.
37: This is not callable or usable in Fortran
39: .seealso: MatMFFDGetH(), MatCreateSNESMF(), MatCreateMFFD(), MATMFFD,
40: MatMFFDSetHHistory(),
41: MatMFFDKSPMonitor(), MatMFFDSetFunctionError(), MatMFFDCreate(), SNESSetJacobian()
43: @*/
44: PetscErrorCode MatMFFDComputeJacobian(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
45: {
48: MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);
49: MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);
50: return(0);
51: }
53: PetscErrorCode MatAssemblyEnd_MFFD(Mat,MatAssemblyType);
56: /*
57: MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
58: base from the SNES context
60: */
61: PetscErrorCode MatAssemblyEnd_SNESMF(Mat J,MatAssemblyType mt)
62: {
64: MatMFFD j = (MatMFFD)J->data;
65: SNES snes = (SNES)j->funcctx;
66: Vec u,f;
69: MatAssemblyEnd_MFFD(J,mt);
71: SNESGetSolution(snes,&u);
72: SNESGetFunction(snes,&f,PETSC_NULL,PETSC_NULL);
73: MatMFFDSetBase(J,u,f);
74: return(0);
75: }
79: /*
80: This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
81: uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
82: */
85: PetscErrorCode MatMFFDSetBase_SNESMF(Mat J,Vec U,Vec F)
86: {
90: MatMFFDSetBase_MFFD(J,U,F);
91: J->ops->assemblyend = MatAssemblyEnd_MFFD;
92: return(0);
93: }
98: /*@
99: MatCreateSNESMF - Creates a matrix-free matrix context for use with
100: a SNES solver. This matrix can be used as the Jacobian argument for
101: the routine SNESSetJacobian(). See MatCreateMFFD() for details on how
102: the finite difference computation is done.
104: Collective on SNES and Vec
106: Input Parameters:
107: . snes - the SNES context
109: Output Parameter:
110: . J - the matrix-free matrix
112: Level: advanced
114: Warning:
115: If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
116: the x from the SNES object and MatMFFDSetBase() must from that point on be used to
117: change the base vector x.
119: Notes: The difference between this routine and MatCreateMFFD() is that this matrix
120: automatically gets the current base vector from the SNES object and not from an
121: explicit call to MatMFFDSetBase().
123: .seealso: MatDestroy(), MatMFFDSetFunctionError(), MatMFFDDefaultSetUmin()
124: MatMFFDSetHHistory(), MatMFFDResetHHistory(), MatCreateMFFD(),
125: MatMFFDGetH(),MatMFFDKSPMonitor(), MatMFFDRegisterDynamic), MatMFFDComputeJacobian()
126:
127: @*/
128: PetscErrorCode MatCreateSNESMF(SNES snes,Mat *J)
129: {
131: PetscInt n,N;
134: if (!snes->vec_func) SETERRQ(PETSC_ERR_ARG_WRONGSTATE,"SNESSetFunction() must be called first");
135:
136: VecGetLocalSize(snes->vec_func,&n);
137: VecGetSize(snes->vec_func,&N);
138: MatCreateMFFD(((PetscObject)snes)->comm,n,n,N,N,J);
139: MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction,snes);
140: (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
141: PetscObjectComposeFunctionDynamic((PetscObject)*J,"MatMFFDSetBase_C","MatMFFDSetBase_SNESMF",MatMFFDSetBase_SNESMF);
142: return(0);
143: }