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