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