xref: /petsc/src/snes/mf/snesmfj.c (revision 3ec795f1650ef19b9b64b98f5b2fa31503b03712)
1 #define PETSCSNES_DLL
2 
3 #include "include/private/snesimpl.h"  /*I  "petscsnes.h" I*/
4 #include "include/private/matimpl.h"
5 #include "src/mat/impls/mffd/mffdimpl.h"
6 
7 #undef __FUNCT__
8 #define __FUNCT__ "MatMFFDComputeJacobian"
9 /*@
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.
12 
13    Collective on SNES
14 
15    Input Parameters:
16 +   snes - the nonlinear solver context
17 .   x - the point at which the Jacobian vector products will be performed
18 .   jac - the matrix-free Jacobian object
19 .   B - either the same as jac or another matrix type (ignored)
20 .   flag - not relevent for matrix-free form
21 -   dummy - the user context (ignored)
22 
23    Level: developer
24 
25    Notes:
26      This can be passed into SNESSetJacobian() when using a completely matrix-free solver,
27      that is the B matrix is also the same matrix operator. This is used when you select
28      -mat_mffd but rarely used directly by users.
29 
30 .seealso: MatMFFDGetH(), MatCreateSNESMF(), MatCreateMFFD(), MATMFFD,
31           MatMFFDSetHHistory(),
32           MatMFFDKSPMonitor(), MatMFFDSetFunctionError(), MatMFFDCreate(), SNESSetJacobian()
33 
34 @*/
35 PetscErrorCode PETSCSNES_DLLEXPORT MatMFFDComputeJacobian(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
36 {
37   PetscErrorCode ierr;
38   PetscFunctionBegin;
39   ierr = MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
40   ierr = MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
41   PetscFunctionReturn(0);
42 }
43 
44 PetscErrorCode MatAssemblyEnd_MFFD(Mat,MatAssemblyType);
45 #undef __FUNCT__
46 #define __FUNCT__ "MatAssemblyEnd_SNESMF"
47 /*
48    MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
49     base from the SNES context
50 
51 */
52 PetscErrorCode MatAssemblyEnd_SNESMF(Mat J,MatAssemblyType mt)
53 {
54   PetscErrorCode ierr;
55   MatMFFD        j = (MatMFFD)J->data;
56   SNES           snes = (SNES)j->funcctx;
57 
58   PetscFunctionBegin;
59   ierr = MatAssemblyEnd_MFFD(J,mt);CHKERRQ(ierr);
60 
61   ierr = SNESGetSolution(snes,&j->current_u);CHKERRQ(ierr);
62   ierr = SNESGetFunction(snes,&j->current_f,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
63   if (!j->w) {
64     ierr = VecDuplicate(j->current_u, &j->w);CHKERRQ(ierr);
65   }
66   PetscFunctionReturn(0);
67 }
68 
69 #undef __FUNCT__
70 #define __FUNCT__ "MatCreateSNESMF"
71 /*@
72    MatCreateSNESMF - Creates a matrix-free matrix context for use with
73    a SNES solver.  This matrix can be used as the Jacobian argument for
74    the routine SNESSetJacobian().
75 
76    Collective on SNES and Vec
77 
78    Input Parameters:
79 +  snes - the SNES context
80 -  x - vector where SNES solution is to be stored.
81 
82    Output Parameter:
83 .  J - the matrix-free matrix
84 
85    Level: advanced
86 
87    Notes:
88    The matrix-free matrix context merely contains the function pointers
89    and work space for performing finite difference approximations of
90    Jacobian-vector products, F'(u)*a,
91 
92    The default code uses the following approach to compute h
93 
94 .vb
95      F'(u)*a = [F(u+h*a) - F(u)]/h where
96      h = error_rel*u'a/||a||^2                        if  |u'a| > umin*||a||_{1}
97        = error_rel*umin*sign(u'a)*||a||_{1}/||a||^2   otherwise
98  where
99      error_rel = square root of relative error in function evaluation
100      umin = minimum iterate parameter
101 .ve
102    (see MATMFFD_WP or MATMFFD_DS)
103 
104    The user can set the error_rel via MatMFFDSetFunctionError() and
105    umin via MatMFFDDefaultSetUmin(); see the nonlinear solvers chapter
106    of the users manual for details.
107 
108    The user should call MatDestroy() when finished with the matrix-free
109    matrix context.
110 
111    Options Database Keys:
112 +  -mat_mffd_err <error_rel> - Sets error_rel
113 +  -mat_mffd_type - wp or ds (see MATMFFD_WP or MATMFFD_DS)
114 .  -mat_mffd_unim <umin> - Sets umin (for default PETSc routine that computes h only)
115 -  -mat_mffd_ksp_monitor - KSP monitor routine that prints differencing h
116 
117 .keywords: SNES, default, matrix-free, create, matrix
118 
119 .seealso: MatDestroy(), MatMFFDSetFunctionError(), MatMFFDDefaultSetUmin()
120           MatMFFDSetHHistory(), MatMFFDResetHHistory(), MatCreateMF(),
121           MatMFFDGetH(),MatMFFDKSPMonitor(), MatMFFDRegisterDynamic), MatMFFDComputeJacobian()
122 
123 @*/
124 PetscErrorCode PETSCSNES_DLLEXPORT MatCreateSNESMF(SNES snes,Vec x,Mat *J)
125 {
126   PetscErrorCode ierr;
127 
128   PetscFunctionBegin;
129   ierr = MatCreateMFFD(x,J);CHKERRQ(ierr);
130   ierr = MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*, _p_Vec*, _p_Vec*))SNESComputeFunction,snes);CHKERRQ(ierr);
131   (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
132 
133   PetscFunctionReturn(0);
134 }
135 
136 
137 
138 
139 
140 
141