xref: /petsc/src/snes/mf/snesmfj.c (revision e0f5d30fb0a27649a1e2f909cb9b25b212638ba7)
1 
2 #include <private/snesimpl.h>  /*I  "petscsnes.h" I*/
3 #include <../src/mat/impls/mffd/mffdimpl.h>
4 #include <private/matimpl.h>
5 
6 #undef __FUNCT__
7 #define __FUNCT__ "MatMFFDComputeJacobian"
8 /*@C
9    MatMFFDComputeJacobian - Tells the matrix-free Jacobian object the new location at which
10        Jacobian matrix vector products will be computed at, i.e. J(x) * a. The x is obtained
11        from the SNES object (using SNESGetSolution()).
12 
13    Logically 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    Warning:
26       If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
27     the x from the SNES object and MatMFFDSetBase() must from that point on be used to
28     change the base vector x.
29 
30    Notes:
31      This can be passed into SNESSetJacobian() 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.
35 
36 .seealso: MatMFFDGetH(), MatCreateSNESMF(), MatCreateMFFD(), MATMFFD,
37           MatMFFDSetHHistory(), MatMFFDSetFunctionError(), MatCreateMFFD(), SNESSetJacobian()
38 
39 @*/
40 PetscErrorCode  MatMFFDComputeJacobian(SNES snes,Vec x,Mat *jac,Mat *B,MatStructure *flag,void *dummy)
41 {
42   PetscErrorCode ierr;
43   PetscFunctionBegin;
44   ierr = MatAssemblyBegin(*jac,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
45   ierr = MatAssemblyEnd(*jac,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr);
46   PetscFunctionReturn(0);
47 }
48 
49 PetscErrorCode MatAssemblyEnd_MFFD(Mat,MatAssemblyType);
50 PetscErrorCode MatMFFDSetBase_MFFD(Mat,Vec,Vec);
51 #undef __FUNCT__
52 #define __FUNCT__ "MatAssemblyEnd_SNESMF"
53 /*
54    MatAssemblyEnd_SNESMF - Calls MatAssemblyEnd_MFFD() and then sets the
55     base from the SNES context
56 
57 */
58 PetscErrorCode MatAssemblyEnd_SNESMF(Mat J,MatAssemblyType mt)
59 {
60   PetscErrorCode ierr;
61   MatMFFD        j = (MatMFFD)J->data;
62   SNES           snes = (SNES)j->funcctx;
63   Vec            u,f;
64 
65   PetscFunctionBegin;
66   ierr = MatAssemblyEnd_MFFD(J,mt);CHKERRQ(ierr);
67 
68   ierr = SNESGetSolution(snes,&u);CHKERRQ(ierr);
69   ierr = SNESGetFunction(snes,&f,PETSC_NULL,PETSC_NULL);CHKERRQ(ierr);
70   ierr = MatMFFDSetBase_MFFD(J,u,f);CHKERRQ(ierr);
71   PetscFunctionReturn(0);
72 }
73 
74 EXTERN_C_BEGIN
75 extern PetscErrorCode  MatMFFDSetBase_MFFD(Mat,Vec,Vec);
76 /*
77     This routine resets the MatAssemblyEnd() for the MatMFFD created from MatCreateSNESMF() so that it NO longer
78   uses the solution in the SNES object to update the base. See the warning in MatCreateSNESMF().
79 */
80 #undef __FUNCT__
81 #define __FUNCT__ "MatMFFDSetBase_SNESMF"
82 PetscErrorCode  MatMFFDSetBase_SNESMF(Mat J,Vec U,Vec F)
83 {
84   PetscErrorCode ierr;
85 
86   PetscFunctionBegin;
87   ierr = MatMFFDSetBase_MFFD(J,U,F);CHKERRQ(ierr);
88   J->ops->assemblyend = MatAssemblyEnd_MFFD;
89   PetscFunctionReturn(0);
90 }
91 EXTERN_C_END
92 
93 #undef __FUNCT__
94 #define __FUNCT__ "MatCreateSNESMF"
95 /*@
96    MatCreateSNESMF - Creates a matrix-free matrix context for use with
97    a SNES solver.  This matrix can be used as the Jacobian argument for
98    the routine SNESSetJacobian(). See MatCreateMFFD() for details on how
99    the finite difference computation is done.
100 
101    Collective on SNES and Vec
102 
103    Input Parameters:
104 .  snes - the SNES context
105 
106    Output Parameter:
107 .  J - the matrix-free matrix
108 
109    Level: advanced
110 
111    Warning:
112       If MatMFFDSetBase() is ever called on jac then this routine will NO longer get
113     the x from the SNES object and MatMFFDSetBase() must from that point on be used to
114     change the base vector x.
115 
116    Notes: The difference between this routine and MatCreateMFFD() is that this matrix
117      automatically gets the current base vector from the SNES object and not from an
118      explicit call to MatMFFDSetBase().
119 
120 .seealso: MatDestroy(), MatMFFDSetFunctionError(), MatMFFDDSSetUmin()
121           MatMFFDSetHHistory(), MatMFFDResetHHistory(), MatCreateMFFD(),
122           MatMFFDGetH(), MatMFFDRegisterDynamic), MatMFFDComputeJacobian()
123 
124 @*/
125 PetscErrorCode  MatCreateSNESMF(SNES snes,Mat *J)
126 {
127   PetscErrorCode ierr;
128   PetscInt       n,N;
129 
130   PetscFunctionBegin;
131   if (snes->vec_func) {
132     ierr = VecGetLocalSize(snes->vec_func,&n);CHKERRQ(ierr);
133     ierr = VecGetSize(snes->vec_func,&N);CHKERRQ(ierr);
134   } else if (snes->dm) {
135     Vec tmp;
136     ierr = DMGetGlobalVector(snes->dm,&tmp);CHKERRQ(ierr);
137     ierr = VecGetLocalSize(tmp,&n);CHKERRQ(ierr);
138     ierr = VecGetSize(tmp,&N);CHKERRQ(ierr);
139     ierr = DMRestoreGlobalVector(snes->dm,&tmp);CHKERRQ(ierr);
140   } else SETERRQ(((PetscObject)snes)->comm,PETSC_ERR_ARG_WRONGSTATE,"Must call SNESSetFunction() or SNESSetDM() first");
141   ierr = MatCreateMFFD(((PetscObject)snes)->comm,n,n,N,N,J);CHKERRQ(ierr);
142   ierr = MatMFFDSetFunction(*J,(PetscErrorCode (*)(void*,Vec,Vec))SNESComputeFunction,snes);CHKERRQ(ierr);
143   (*J)->ops->assemblyend = MatAssemblyEnd_SNESMF;
144   ierr = PetscObjectComposeFunctionDynamic((PetscObject)*J,"MatMFFDSetBase_C","MatMFFDSetBase_SNESMF",MatMFFDSetBase_SNESMF);CHKERRQ(ierr);
145   PetscFunctionReturn(0);
146 }
147 
148 
149 
150 
151 
152 
153