1 #include <petscdmda.h> /*I "petscdmda.h" I*/ 2 #include <petsc/private/dmimpl.h> 3 #include <petsc/private/snesimpl.h> /*I "petscsnes.h" I*/ 4 5 /* This structure holds the user-provided DMDA callbacks */ 6 typedef struct { 7 PetscErrorCode (*residuallocal)(DMDALocalInfo*,void*,void*,void*); 8 PetscErrorCode (*jacobianlocal)(DMDALocalInfo*,void*,Mat,Mat,void*); 9 PetscErrorCode (*objectivelocal)(DMDALocalInfo*,void*,PetscReal*,void*); 10 void *residuallocalctx; 11 void *jacobianlocalctx; 12 void *objectivelocalctx; 13 InsertMode residuallocalimode; 14 15 /* For Picard iteration defined locally */ 16 PetscErrorCode (*rhsplocal)(DMDALocalInfo*,void*,void*,void*); 17 PetscErrorCode (*jacobianplocal)(DMDALocalInfo*,void*,Mat,Mat,void*); 18 void *picardlocalctx; 19 } DMSNES_DA; 20 21 static PetscErrorCode DMSNESDestroy_DMDA(DMSNES sdm) 22 { 23 PetscErrorCode ierr; 24 25 PetscFunctionBegin; 26 ierr = PetscFree(sdm->data);CHKERRQ(ierr); 27 PetscFunctionReturn(0); 28 } 29 30 static PetscErrorCode DMSNESDuplicate_DMDA(DMSNES oldsdm,DMSNES sdm) 31 { 32 PetscErrorCode ierr; 33 34 PetscFunctionBegin; 35 ierr = PetscNewLog(sdm,(DMSNES_DA**)&sdm->data);CHKERRQ(ierr); 36 if (oldsdm->data) { 37 ierr = PetscMemcpy(sdm->data,oldsdm->data,sizeof(DMSNES_DA));CHKERRQ(ierr); 38 } 39 PetscFunctionReturn(0); 40 } 41 42 43 static PetscErrorCode DMDASNESGetContext(DM dm,DMSNES sdm,DMSNES_DA **dmdasnes) 44 { 45 PetscErrorCode ierr; 46 47 PetscFunctionBegin; 48 *dmdasnes = NULL; 49 if (!sdm->data) { 50 ierr = PetscNewLog(dm,(DMSNES_DA**)&sdm->data);CHKERRQ(ierr); 51 sdm->ops->destroy = DMSNESDestroy_DMDA; 52 sdm->ops->duplicate = DMSNESDuplicate_DMDA; 53 } 54 *dmdasnes = (DMSNES_DA*)sdm->data; 55 PetscFunctionReturn(0); 56 } 57 58 static PetscErrorCode SNESComputeFunction_DMDA(SNES snes,Vec X,Vec F,void *ctx) 59 { 60 PetscErrorCode ierr; 61 DM dm; 62 DMSNES_DA *dmdasnes = (DMSNES_DA*)ctx; 63 DMDALocalInfo info; 64 Vec Xloc; 65 void *x,*f; 66 67 PetscFunctionBegin; 68 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 69 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 70 PetscValidHeaderSpecific(F,VEC_CLASSID,3); 71 if (!dmdasnes->residuallocal) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Corrupt context"); 72 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 73 ierr = DMGetLocalVector(dm,&Xloc);CHKERRQ(ierr); 74 ierr = DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 75 ierr = DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 76 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 77 ierr = DMDAVecGetArray(dm,Xloc,&x);CHKERRQ(ierr); 78 switch (dmdasnes->residuallocalimode) { 79 case INSERT_VALUES: { 80 ierr = DMDAVecGetArray(dm,F,&f);CHKERRQ(ierr); 81 ierr = PetscLogEventBegin(SNES_FunctionEval,snes,X,F,0);CHKERRQ(ierr); 82 CHKMEMQ; 83 ierr = (*dmdasnes->residuallocal)(&info,x,f,dmdasnes->residuallocalctx);CHKERRQ(ierr); 84 CHKMEMQ; 85 ierr = PetscLogEventEnd(SNES_FunctionEval,snes,X,F,0);CHKERRQ(ierr); 86 ierr = DMDAVecRestoreArray(dm,F,&f);CHKERRQ(ierr); 87 } break; 88 case ADD_VALUES: { 89 Vec Floc; 90 ierr = DMGetLocalVector(dm,&Floc);CHKERRQ(ierr); 91 ierr = VecZeroEntries(Floc);CHKERRQ(ierr); 92 ierr = DMDAVecGetArray(dm,Floc,&f);CHKERRQ(ierr); 93 ierr = PetscLogEventBegin(SNES_FunctionEval,snes,X,F,0);CHKERRQ(ierr); 94 CHKMEMQ; 95 ierr = (*dmdasnes->residuallocal)(&info,x,f,dmdasnes->residuallocalctx);CHKERRQ(ierr); 96 CHKMEMQ; 97 ierr = PetscLogEventEnd(SNES_FunctionEval,snes,X,F,0);CHKERRQ(ierr); 98 ierr = DMDAVecRestoreArray(dm,Floc,&f);CHKERRQ(ierr); 99 ierr = VecZeroEntries(F);CHKERRQ(ierr); 100 ierr = DMLocalToGlobalBegin(dm,Floc,ADD_VALUES,F);CHKERRQ(ierr); 101 ierr = DMLocalToGlobalEnd(dm,Floc,ADD_VALUES,F);CHKERRQ(ierr); 102 ierr = DMRestoreLocalVector(dm,&Floc);CHKERRQ(ierr); 103 } break; 104 default: SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_INCOMP,"Cannot use imode=%d",(int)dmdasnes->residuallocalimode); 105 } 106 ierr = DMDAVecRestoreArray(dm,Xloc,&x);CHKERRQ(ierr); 107 ierr = DMRestoreLocalVector(dm,&Xloc);CHKERRQ(ierr); 108 if (snes->domainerror) { 109 ierr = VecSetInf(F);CHKERRQ(ierr); 110 } 111 PetscFunctionReturn(0); 112 } 113 114 static PetscErrorCode SNESComputeObjective_DMDA(SNES snes,Vec X,PetscReal *ob,void *ctx) 115 { 116 PetscErrorCode ierr; 117 DM dm; 118 DMSNES_DA *dmdasnes = (DMSNES_DA*)ctx; 119 DMDALocalInfo info; 120 Vec Xloc; 121 void *x; 122 123 PetscFunctionBegin; 124 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 125 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 126 PetscValidPointer(ob,3); 127 if (!dmdasnes->objectivelocal) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Corrupt context"); 128 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 129 ierr = DMGetLocalVector(dm,&Xloc);CHKERRQ(ierr); 130 ierr = DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 131 ierr = DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 132 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 133 ierr = DMDAVecGetArray(dm,Xloc,&x);CHKERRQ(ierr); 134 CHKMEMQ; 135 ierr = (*dmdasnes->objectivelocal)(&info,x,ob,dmdasnes->objectivelocalctx);CHKERRQ(ierr); 136 CHKMEMQ; 137 ierr = DMDAVecRestoreArray(dm,Xloc,&x);CHKERRQ(ierr); 138 ierr = DMRestoreLocalVector(dm,&Xloc);CHKERRQ(ierr); 139 PetscFunctionReturn(0); 140 } 141 142 143 /* Routine is called by example, hence must be labeled PETSC_EXTERN */ 144 PETSC_EXTERN PetscErrorCode SNESComputeJacobian_DMDA(SNES snes,Vec X,Mat A,Mat B,void *ctx) 145 { 146 PetscErrorCode ierr; 147 DM dm; 148 DMSNES_DA *dmdasnes = (DMSNES_DA*)ctx; 149 DMDALocalInfo info; 150 Vec Xloc; 151 void *x; 152 153 PetscFunctionBegin; 154 if (!dmdasnes->residuallocal) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Corrupt context"); 155 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 156 157 if (dmdasnes->jacobianlocal) { 158 ierr = DMGetLocalVector(dm,&Xloc);CHKERRQ(ierr); 159 ierr = DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 160 ierr = DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 161 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 162 ierr = DMDAVecGetArray(dm,Xloc,&x);CHKERRQ(ierr); 163 CHKMEMQ; 164 ierr = (*dmdasnes->jacobianlocal)(&info,x,A,B,dmdasnes->jacobianlocalctx);CHKERRQ(ierr); 165 CHKMEMQ; 166 ierr = DMDAVecRestoreArray(dm,Xloc,&x);CHKERRQ(ierr); 167 ierr = DMRestoreLocalVector(dm,&Xloc);CHKERRQ(ierr); 168 } else { 169 MatFDColoring fdcoloring; 170 ierr = PetscObjectQuery((PetscObject)dm,"DMDASNES_FDCOLORING",(PetscObject*)&fdcoloring);CHKERRQ(ierr); 171 if (!fdcoloring) { 172 ISColoring coloring; 173 174 ierr = DMCreateColoring(dm,dm->coloringtype,&coloring);CHKERRQ(ierr); 175 ierr = MatFDColoringCreate(B,coloring,&fdcoloring);CHKERRQ(ierr); 176 switch (dm->coloringtype) { 177 case IS_COLORING_GLOBAL: 178 ierr = MatFDColoringSetFunction(fdcoloring,(PetscErrorCode (*)(void))SNESComputeFunction_DMDA,dmdasnes);CHKERRQ(ierr); 179 break; 180 default: SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_SUP,"No support for coloring type '%s'",ISColoringTypes[dm->coloringtype]); 181 } 182 ierr = PetscObjectSetOptionsPrefix((PetscObject)fdcoloring,((PetscObject)dm)->prefix);CHKERRQ(ierr); 183 ierr = MatFDColoringSetFromOptions(fdcoloring);CHKERRQ(ierr); 184 ierr = MatFDColoringSetUp(B,coloring,fdcoloring);CHKERRQ(ierr); 185 ierr = ISColoringDestroy(&coloring);CHKERRQ(ierr); 186 ierr = PetscObjectCompose((PetscObject)dm,"DMDASNES_FDCOLORING",(PetscObject)fdcoloring);CHKERRQ(ierr); 187 ierr = PetscObjectDereference((PetscObject)fdcoloring);CHKERRQ(ierr); 188 189 /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call 190 * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the 191 * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually 192 * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be 193 * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed. 194 */ 195 ierr = PetscObjectDereference((PetscObject)dm);CHKERRQ(ierr); 196 } 197 ierr = MatFDColoringApply(B,fdcoloring,X,snes);CHKERRQ(ierr); 198 } 199 /* This will be redundant if the user called both, but it's too common to forget. */ 200 if (A != B) { 201 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 202 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 203 } 204 PetscFunctionReturn(0); 205 } 206 207 /*@C 208 DMDASNESSetFunctionLocal - set a local residual evaluation function 209 210 Logically Collective 211 212 Input Arguments: 213 + dm - DM to associate callback with 214 . imode - INSERT_VALUES if local function computes owned part, ADD_VALUES if it contributes to ghosted part 215 . func - local residual evaluation 216 - ctx - optional context for local residual evaluation 217 218 Calling sequence: 219 For PetscErrorCode (*func)(DMDALocalInfo *info,void *x, void *f, void *ctx), 220 + info - DMDALocalInfo defining the subdomain to evaluate the residual on 221 . x - dimensional pointer to state at which to evaluate residual (e.g. PetscScalar *x or **x or ***x) 222 . f - dimensional pointer to residual, write the residual here (e.g. PetscScalar *f or **f or ***f) 223 - ctx - optional context passed above 224 225 Level: beginner 226 227 .seealso: DMDASNESSetJacobianLocal(), DMSNESSetFunction(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 228 @*/ 229 PetscErrorCode DMDASNESSetFunctionLocal(DM dm,InsertMode imode,PetscErrorCode (*func)(DMDALocalInfo*,void*,void*,void*),void *ctx) 230 { 231 PetscErrorCode ierr; 232 DMSNES sdm; 233 DMSNES_DA *dmdasnes; 234 235 PetscFunctionBegin; 236 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 237 ierr = DMGetDMSNESWrite(dm,&sdm);CHKERRQ(ierr); 238 ierr = DMDASNESGetContext(dm,sdm,&dmdasnes);CHKERRQ(ierr); 239 240 dmdasnes->residuallocalimode = imode; 241 dmdasnes->residuallocal = func; 242 dmdasnes->residuallocalctx = ctx; 243 244 ierr = DMSNESSetFunction(dm,SNESComputeFunction_DMDA,dmdasnes);CHKERRQ(ierr); 245 if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */ 246 ierr = DMSNESSetJacobian(dm,SNESComputeJacobian_DMDA,dmdasnes);CHKERRQ(ierr); 247 } 248 PetscFunctionReturn(0); 249 } 250 251 /*@C 252 DMDASNESSetJacobianLocal - set a local Jacobian evaluation function 253 254 Logically Collective 255 256 Input Arguments: 257 + dm - DM to associate callback with 258 . func - local Jacobian evaluation 259 - ctx - optional context for local Jacobian evaluation 260 261 Calling sequence: 262 For PetscErrorCode (*func)(DMDALocalInfo *info,void *x,Mat J,Mat M,void *ctx), 263 + info - DMDALocalInfo defining the subdomain to evaluate the Jacobian at 264 . x - dimensional pointer to state at which to evaluate Jacobian (e.g. PetscScalar *x or **x or ***x) 265 . J - Mat object for the Jacobian 266 . M - Mat object for the Jacobian preconditioner matrix 267 - ctx - optional context passed above 268 269 Level: beginner 270 271 .seealso: DMDASNESSetFunctionLocal(), DMSNESSetJacobian(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 272 @*/ 273 PetscErrorCode DMDASNESSetJacobianLocal(DM dm,PetscErrorCode (*func)(DMDALocalInfo*,void*,Mat,Mat,void*),void *ctx) 274 { 275 PetscErrorCode ierr; 276 DMSNES sdm; 277 DMSNES_DA *dmdasnes; 278 279 PetscFunctionBegin; 280 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 281 ierr = DMGetDMSNESWrite(dm,&sdm);CHKERRQ(ierr); 282 ierr = DMDASNESGetContext(dm,sdm,&dmdasnes);CHKERRQ(ierr); 283 284 dmdasnes->jacobianlocal = func; 285 dmdasnes->jacobianlocalctx = ctx; 286 287 ierr = DMSNESSetJacobian(dm,SNESComputeJacobian_DMDA,dmdasnes);CHKERRQ(ierr); 288 PetscFunctionReturn(0); 289 } 290 291 292 /*@C 293 DMDASNESSetObjectiveLocal - set a local residual evaluation function 294 295 Logically Collective 296 297 Input Arguments: 298 + dm - DM to associate callback with 299 . func - local objective evaluation 300 - ctx - optional context for local residual evaluation 301 302 Calling sequence for func: 303 + info - DMDALocalInfo defining the subdomain to evaluate the residual on 304 . x - dimensional pointer to state at which to evaluate residual 305 . ob - eventual objective value 306 - ctx - optional context passed above 307 308 Level: beginner 309 310 .seealso: DMSNESSetFunction(), DMDASNESSetJacobianLocal(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 311 @*/ 312 PetscErrorCode DMDASNESSetObjectiveLocal(DM dm,DMDASNESObjective func,void *ctx) 313 { 314 PetscErrorCode ierr; 315 DMSNES sdm; 316 DMSNES_DA *dmdasnes; 317 318 PetscFunctionBegin; 319 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 320 ierr = DMGetDMSNESWrite(dm,&sdm);CHKERRQ(ierr); 321 ierr = DMDASNESGetContext(dm,sdm,&dmdasnes);CHKERRQ(ierr); 322 323 dmdasnes->objectivelocal = func; 324 dmdasnes->objectivelocalctx = ctx; 325 326 ierr = DMSNESSetObjective(dm,SNESComputeObjective_DMDA,dmdasnes);CHKERRQ(ierr); 327 PetscFunctionReturn(0); 328 } 329 330 static PetscErrorCode SNESComputePicard_DMDA(SNES snes,Vec X,Vec F,void *ctx) 331 { 332 PetscErrorCode ierr; 333 DM dm; 334 DMSNES_DA *dmdasnes = (DMSNES_DA*)ctx; 335 DMDALocalInfo info; 336 Vec Xloc; 337 void *x,*f; 338 339 PetscFunctionBegin; 340 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 341 PetscValidHeaderSpecific(X,VEC_CLASSID,2); 342 PetscValidHeaderSpecific(F,VEC_CLASSID,3); 343 if (!dmdasnes->rhsplocal) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Corrupt context"); 344 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 345 ierr = DMGetLocalVector(dm,&Xloc);CHKERRQ(ierr); 346 ierr = DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 347 ierr = DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 348 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 349 ierr = DMDAVecGetArray(dm,Xloc,&x);CHKERRQ(ierr); 350 switch (dmdasnes->residuallocalimode) { 351 case INSERT_VALUES: { 352 ierr = DMDAVecGetArray(dm,F,&f);CHKERRQ(ierr); 353 CHKMEMQ; 354 ierr = (*dmdasnes->rhsplocal)(&info,x,f,dmdasnes->picardlocalctx);CHKERRQ(ierr); 355 CHKMEMQ; 356 ierr = DMDAVecRestoreArray(dm,F,&f);CHKERRQ(ierr); 357 } break; 358 case ADD_VALUES: { 359 Vec Floc; 360 ierr = DMGetLocalVector(dm,&Floc);CHKERRQ(ierr); 361 ierr = VecZeroEntries(Floc);CHKERRQ(ierr); 362 ierr = DMDAVecGetArray(dm,Floc,&f);CHKERRQ(ierr); 363 CHKMEMQ; 364 ierr = (*dmdasnes->rhsplocal)(&info,x,f,dmdasnes->picardlocalctx);CHKERRQ(ierr); 365 CHKMEMQ; 366 ierr = DMDAVecRestoreArray(dm,Floc,&f);CHKERRQ(ierr); 367 ierr = VecZeroEntries(F);CHKERRQ(ierr); 368 ierr = DMLocalToGlobalBegin(dm,Floc,ADD_VALUES,F);CHKERRQ(ierr); 369 ierr = DMLocalToGlobalEnd(dm,Floc,ADD_VALUES,F);CHKERRQ(ierr); 370 ierr = DMRestoreLocalVector(dm,&Floc);CHKERRQ(ierr); 371 } break; 372 default: SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_INCOMP,"Cannot use imode=%d",(int)dmdasnes->residuallocalimode); 373 } 374 ierr = DMDAVecRestoreArray(dm,Xloc,&x);CHKERRQ(ierr); 375 ierr = DMRestoreLocalVector(dm,&Xloc);CHKERRQ(ierr); 376 PetscFunctionReturn(0); 377 } 378 379 static PetscErrorCode SNESComputePicardJacobian_DMDA(SNES snes,Vec X,Mat A,Mat B,void *ctx) 380 { 381 PetscErrorCode ierr; 382 DM dm; 383 DMSNES_DA *dmdasnes = (DMSNES_DA*)ctx; 384 DMDALocalInfo info; 385 Vec Xloc; 386 void *x; 387 388 PetscFunctionBegin; 389 if (!dmdasnes->jacobianplocal) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_PLIB,"Corrupt context"); 390 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 391 392 ierr = DMGetLocalVector(dm,&Xloc);CHKERRQ(ierr); 393 ierr = DMGlobalToLocalBegin(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 394 ierr = DMGlobalToLocalEnd(dm,X,INSERT_VALUES,Xloc);CHKERRQ(ierr); 395 ierr = DMDAGetLocalInfo(dm,&info);CHKERRQ(ierr); 396 ierr = DMDAVecGetArray(dm,Xloc,&x);CHKERRQ(ierr); 397 CHKMEMQ; 398 ierr = (*dmdasnes->jacobianplocal)(&info,x,A,B,dmdasnes->picardlocalctx);CHKERRQ(ierr); 399 CHKMEMQ; 400 ierr = DMDAVecRestoreArray(dm,Xloc,&x);CHKERRQ(ierr); 401 ierr = DMRestoreLocalVector(dm,&Xloc);CHKERRQ(ierr); 402 if (A != B) { 403 ierr = MatAssemblyBegin(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 404 ierr = MatAssemblyEnd(A,MAT_FINAL_ASSEMBLY);CHKERRQ(ierr); 405 } 406 PetscFunctionReturn(0); 407 } 408 409 /*@C 410 DMDASNESSetPicardLocal - set a local right hand side and matrix evaluation function for Picard iteration 411 412 Logically Collective 413 414 Input Arguments: 415 + dm - DM to associate callback with 416 . imode - INSERT_VALUES if local function computes owned part, ADD_VALUES if it contributes to ghosted part 417 . func - local residual evaluation 418 - ctx - optional context for local residual evaluation 419 420 Calling sequence for func: 421 + info - DMDALocalInfo defining the subdomain to evaluate the residual on 422 . x - dimensional pointer to state at which to evaluate residual 423 . f - dimensional pointer to residual, write the residual here 424 - ctx - optional context passed above 425 426 Notes: The user must use 427 extern PetscErrorCode SNESPicardComputeFunction(SNES,Vec,Vec,void*); 428 extern PetscErrorCode SNESPicardComputeJacobian(SNES,Vec,Mat,Mat,MatStructure*,void*); 429 ierr = SNESSetFunction(snes,NULL,SNESPicardComputeFunction,&user);CHKERRQ(ierr); 430 in their code before calling this routine. 431 432 433 Level: beginner 434 435 .seealso: DMSNESSetFunction(), DMDASNESSetJacobian(), DMDACreate1d(), DMDACreate2d(), DMDACreate3d() 436 @*/ 437 PetscErrorCode DMDASNESSetPicardLocal(DM dm,InsertMode imode,PetscErrorCode (*func)(DMDALocalInfo*,void*,void*,void*), 438 PetscErrorCode (*jac)(DMDALocalInfo*,void*,Mat,Mat,void*),void *ctx) 439 { 440 PetscErrorCode ierr; 441 DMSNES sdm; 442 DMSNES_DA *dmdasnes; 443 444 PetscFunctionBegin; 445 PetscValidHeaderSpecific(dm,DM_CLASSID,1); 446 ierr = DMGetDMSNESWrite(dm,&sdm);CHKERRQ(ierr); 447 ierr = DMDASNESGetContext(dm,sdm,&dmdasnes);CHKERRQ(ierr); 448 449 dmdasnes->residuallocalimode = imode; 450 dmdasnes->rhsplocal = func; 451 dmdasnes->jacobianplocal = jac; 452 dmdasnes->picardlocalctx = ctx; 453 454 ierr = DMSNESSetPicard(dm,SNESComputePicard_DMDA,SNESComputePicardJacobian_DMDA,dmdasnes);CHKERRQ(ierr); 455 PetscFunctionReturn(0); 456 } 457