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