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