Actual source code: dmdasnes.c

      1: #include <petscdmda.h>
      2: #include <petsc/private/dmimpl.h>
      3: #include <petsc/private/snesimpl.h>

      5: /* This structure holds the user-provided DMDA callbacks */
      6: typedef struct {
      7:   /* array versions for vector data */
      8:   DMDASNESFunctionFn  *residuallocal;
      9:   DMDASNESJacobianFn  *jacobianlocal;
     10:   DMDASNESObjectiveFn *objectivelocal;

     12:   /* Vec version for vector data */
     13:   DMDASNESFunctionVecFn  *residuallocalvec;
     14:   DMDASNESJacobianVecFn  *jacobianlocalvec;
     15:   DMDASNESObjectiveVecFn *objectivelocalvec;

     17:   /* user contexts */
     18:   void      *residuallocalctx;
     19:   void      *jacobianlocalctx;
     20:   void      *objectivelocalctx;
     21:   InsertMode residuallocalimode;

     23:   /*   For Picard iteration defined locally */
     24:   PetscErrorCode (*rhsplocal)(DMDALocalInfo *, void *, void *, void *);
     25:   PetscErrorCode (*jacobianplocal)(DMDALocalInfo *, void *, Mat, Mat, void *);
     26:   void *picardlocalctx;
     27: } DMSNES_DA;

     29: static PetscErrorCode DMSNESDestroy_DMDA(DMSNES sdm)
     30: {
     31:   PetscFunctionBegin;
     32:   PetscCall(PetscFree(sdm->data));
     33:   PetscFunctionReturn(PETSC_SUCCESS);
     34: }

     36: static PetscErrorCode DMSNESDuplicate_DMDA(DMSNES oldsdm, DMSNES sdm)
     37: {
     38:   PetscFunctionBegin;
     39:   PetscCall(PetscNew((DMSNES_DA **)&sdm->data));
     40:   if (oldsdm->data) PetscCall(PetscMemcpy(sdm->data, oldsdm->data, sizeof(DMSNES_DA)));
     41:   PetscFunctionReturn(PETSC_SUCCESS);
     42: }

     44: static PetscErrorCode DMDASNESGetContext(DM dm, DMSNES sdm, DMSNES_DA **dmdasnes)
     45: {
     46:   PetscFunctionBegin;
     47:   *dmdasnes = NULL;
     48:   if (!sdm->data) {
     49:     PetscCall(PetscNew((DMSNES_DA **)&sdm->data));
     50:     sdm->ops->destroy   = DMSNESDestroy_DMDA;
     51:     sdm->ops->duplicate = DMSNESDuplicate_DMDA;
     52:   }
     53:   *dmdasnes = (DMSNES_DA *)sdm->data;
     54:   PetscFunctionReturn(PETSC_SUCCESS);
     55: }

     57: static PetscErrorCode SNESComputeFunction_DMDA(SNES snes, Vec X, Vec F, void *ctx)
     58: {
     59:   DM            dm;
     60:   DMSNES_DA    *dmdasnes = (DMSNES_DA *)ctx;
     61:   DMDALocalInfo info;
     62:   Vec           Xloc;
     63:   void         *x, *f, *rctx;

     65:   PetscFunctionBegin;
     69:   PetscCheck(dmdasnes->residuallocal || dmdasnes->residuallocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
     70:   PetscCall(SNESGetDM(snes, &dm));
     71:   PetscCall(DMGetLocalVector(dm, &Xloc));
     72:   PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
     73:   PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
     74:   PetscCall(DMDAGetLocalInfo(dm, &info));
     75:   rctx = dmdasnes->residuallocalctx ? dmdasnes->residuallocalctx : snes->user;
     76:   switch (dmdasnes->residuallocalimode) {
     77:   case INSERT_VALUES: {
     78:     PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0));
     79:     if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, F, rctx));
     80:     else {
     81:       PetscCall(DMDAVecGetArray(dm, Xloc, &x));
     82:       PetscCall(DMDAVecGetArray(dm, F, &f));
     83:       PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, rctx));
     84:       PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
     85:       PetscCall(DMDAVecRestoreArray(dm, F, &f));
     86:     }
     87:     PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0));
     88:   } break;
     89:   case ADD_VALUES: {
     90:     Vec Floc;
     91:     PetscCall(DMGetLocalVector(dm, &Floc));
     92:     PetscCall(VecZeroEntries(Floc));
     93:     PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0));
     94:     if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, Floc, rctx));
     95:     else {
     96:       PetscCall(DMDAVecGetArray(dm, Xloc, &x));
     97:       PetscCall(DMDAVecGetArray(dm, Floc, &f));
     98:       PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, rctx));
     99:       PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
    100:       PetscCall(DMDAVecRestoreArray(dm, Floc, &f));
    101:     }
    102:     PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0));
    103:     PetscCall(VecZeroEntries(F));
    104:     PetscCall(DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F));
    105:     PetscCall(DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F));
    106:     PetscCall(DMRestoreLocalVector(dm, &Floc));
    107:   } break;
    108:   default:
    109:     SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
    110:   }
    111:   PetscCall(DMRestoreLocalVector(dm, &Xloc));
    112:   PetscCall(VecFlag(F, snes->domainerror));
    113:   PetscFunctionReturn(PETSC_SUCCESS);
    114: }

    116: static PetscErrorCode SNESComputeObjective_DMDA(SNES snes, Vec X, PetscReal *ob, void *ctx)
    117: {
    118:   DM            dm;
    119:   DMSNES_DA    *dmdasnes = (DMSNES_DA *)ctx;
    120:   DMDALocalInfo info;
    121:   Vec           Xloc;
    122:   void         *x, *octx;

    124:   PetscFunctionBegin;
    127:   PetscAssertPointer(ob, 3);
    128:   PetscCheck(dmdasnes->objectivelocal || dmdasnes->objectivelocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
    129:   PetscCall(SNESGetDM(snes, &dm));
    130:   PetscCall(DMGetLocalVector(dm, &Xloc));
    131:   PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
    132:   PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
    133:   PetscCall(DMDAGetLocalInfo(dm, &info));
    134:   octx = dmdasnes->objectivelocalctx ? dmdasnes->objectivelocalctx : snes->user;
    135:   if (dmdasnes->objectivelocalvec) PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocalvec)(&info, Xloc, ob, octx));
    136:   else {
    137:     PetscCall(DMDAVecGetArray(dm, Xloc, &x));
    138:     PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocal)(&info, x, ob, octx));
    139:     PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
    140:   }
    141:   PetscCall(DMRestoreLocalVector(dm, &Xloc));
    142:   PetscCallMPI(MPIU_Allreduce(MPI_IN_PLACE, ob, 1, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)snes)));
    143:   PetscFunctionReturn(PETSC_SUCCESS);
    144: }

    146: /* Routine is called by example, hence must be labeled PETSC_EXTERN */
    147: PETSC_EXTERN PetscErrorCode SNESComputeJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
    148: {
    149:   DM            dm;
    150:   DMSNES_DA    *dmdasnes = (DMSNES_DA *)ctx;
    151:   DMDALocalInfo info;
    152:   Vec           Xloc;
    153:   void         *x, *jctx;

    155:   PetscFunctionBegin;
    156:   PetscCheck(dmdasnes->residuallocal || dmdasnes->residuallocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
    157:   PetscCall(SNESGetDM(snes, &dm));
    158:   jctx = dmdasnes->jacobianlocalctx ? dmdasnes->jacobianlocalctx : snes->user;
    159:   if (dmdasnes->jacobianlocal || dmdasnes->jacobianlocalvec) {
    160:     PetscCall(DMGetLocalVector(dm, &Xloc));
    161:     PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
    162:     PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
    163:     PetscCall(DMDAGetLocalInfo(dm, &info));
    164:     if (dmdasnes->jacobianlocalvec) PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocalvec)(&info, Xloc, A, B, jctx));
    165:     else {
    166:       PetscCall(DMDAVecGetArray(dm, Xloc, &x));
    167:       PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocal)(&info, x, A, B, jctx));
    168:       PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
    169:     }
    170:     PetscCall(DMRestoreLocalVector(dm, &Xloc));
    171:   } else {
    172:     MatFDColoring fdcoloring;
    173:     PetscCall(PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring));
    174:     if (!fdcoloring) {
    175:       ISColoring coloring;

    177:       PetscCall(DMCreateColoring(dm, dm->coloringtype, &coloring));
    178:       PetscCall(MatFDColoringCreate(B, coloring, &fdcoloring));
    179:       switch (dm->coloringtype) {
    180:       case IS_COLORING_GLOBAL:
    181:         PetscCall(MatFDColoringSetFunction(fdcoloring, (PetscErrorCode (*)(void))SNESComputeFunction_DMDA, dmdasnes));
    182:         break;
    183:       default:
    184:         SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
    185:       }
    186:       PetscCall(PetscObjectSetOptionsPrefix((PetscObject)fdcoloring, ((PetscObject)dm)->prefix));
    187:       PetscCall(MatFDColoringSetFromOptions(fdcoloring));
    188:       PetscCall(MatFDColoringSetUp(B, coloring, fdcoloring));
    189:       PetscCall(ISColoringDestroy(&coloring));
    190:       PetscCall(PetscObjectCompose((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject)fdcoloring));
    191:       PetscCall(PetscObjectDereference((PetscObject)fdcoloring));

    193:       /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
    194:        * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
    195:        * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
    196:        * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
    197:        * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
    198:        */
    199:       PetscCall(PetscObjectDereference((PetscObject)dm));
    200:     }
    201:     PetscCall(MatFDColoringApply(B, fdcoloring, X, snes));
    202:   }
    203:   /* This will be redundant if the user called both, but it's too common to forget. */
    204:   if (A != B) {
    205:     PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
    206:     PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
    207:   }
    208:   PetscFunctionReturn(PETSC_SUCCESS);
    209: }

    211: /*@C
    212:   DMDASNESSetFunctionLocal - set a local residual evaluation function for use with `DMDA`

    214:   Logically Collective

    216:   Input Parameters:
    217: + dm    - `DM` to associate callback with
    218: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
    219: . func  - local residual evaluation
    220: - ctx   - optional context for local residual evaluation

    222:   Calling sequence of `func`:
    223: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
    224: . x    - dimensional pointer to state at which to evaluate residual (e.g. PetscScalar *x or **x or ***x)
    225: . f    - dimensional pointer to residual, write the residual here (e.g. PetscScalar *f or **f or ***f)
    226: - ctx  - optional context passed above

    228:   Level: beginner

    230: .seealso: [](ch_snes), `DMDA`, `DMDASNESSetJacobianLocal()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
    231: @*/
    232: PetscErrorCode DMDASNESSetFunctionLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *info, void *x, void *f, void *ctx), void *ctx)
    233: {
    234:   DMSNES     sdm;
    235:   DMSNES_DA *dmdasnes;

    237:   PetscFunctionBegin;
    239:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    240:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    242:   dmdasnes->residuallocalimode = imode;
    243:   dmdasnes->residuallocal      = func;
    244:   dmdasnes->residuallocalctx   = ctx;

    246:   PetscCall(DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
    247:   if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
    248:     PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
    249:   }
    250:   PetscFunctionReturn(PETSC_SUCCESS);
    251: }

    253: /*@C
    254:   DMDASNESSetFunctionLocalVec - set a local residual evaluation function that operates on a local vector for `DMDA`

    256:   Logically Collective

    258:   Input Parameters:
    259: + dm    - `DM` to associate callback with
    260: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
    261: . func  - local residual evaluation
    262: - ctx   - optional context for local residual evaluation

    264:   Calling sequence of `func`:
    265: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
    266: . x    - state vector at which to evaluate residual
    267: . f    - residual vector
    268: - ctx  - optional context passed above

    270:   Level: beginner

    272: .seealso: [](ch_snes), `DMDA`, `DMDASNESSetFunctionLocal()`, `DMDASNESSetJacobianLocalVec()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
    273: @*/
    274: PetscErrorCode DMDASNESSetFunctionLocalVec(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *info, Vec x, Vec f, void *ctx), void *ctx)
    275: {
    276:   DMSNES     sdm;
    277:   DMSNES_DA *dmdasnes;

    279:   PetscFunctionBegin;
    281:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    282:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    284:   dmdasnes->residuallocalimode = imode;
    285:   dmdasnes->residuallocalvec   = func;
    286:   dmdasnes->residuallocalctx   = ctx;

    288:   PetscCall(DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
    289:   if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
    290:     PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
    291:   }
    292:   PetscFunctionReturn(PETSC_SUCCESS);
    293: }

    295: /*@C
    296:   DMDASNESSetJacobianLocal - set a local Jacobian evaluation function for use with `DMDA`

    298:   Logically Collective

    300:   Input Parameters:
    301: + dm   - `DM` to associate callback with
    302: . func - local Jacobian evaluation function
    303: - ctx  - optional context for local Jacobian evaluation

    305:   Calling sequence of `func`:
    306: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
    307: . x    - dimensional pointer to state at which to evaluate Jacobian (e.g. PetscScalar *x or **x or ***x)
    308: . J    - `Mat` object for the Jacobian
    309: . M    - `Mat` object for the Jacobian preconditioner matrix, often `J`
    310: - ctx  - optional context passed above

    312:   Level: beginner

    314:   Note:
    315:   The `J` and `M` matrices are created internally by `DMCreateMatrix()`

    317: .seealso: [](ch_snes), `DMDA`, `DMDASNESSetFunctionLocal()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
    318: @*/
    319: PetscErrorCode DMDASNESSetJacobianLocal(DM dm, PetscErrorCode (*func)(DMDALocalInfo *info, void *x, Mat J, Mat M, void *ctx), void *ctx)
    320: {
    321:   DMSNES     sdm;
    322:   DMSNES_DA *dmdasnes;

    324:   PetscFunctionBegin;
    326:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    327:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    329:   dmdasnes->jacobianlocal    = func;
    330:   dmdasnes->jacobianlocalctx = ctx;

    332:   PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
    333:   PetscFunctionReturn(PETSC_SUCCESS);
    334: }

    336: /*@C
    337:   DMDASNESSetJacobianLocalVec - set a local Jacobian evaluation function that operates on a local vector with `DMDA`

    339:   Logically Collective

    341:   Input Parameters:
    342: + dm   - `DM` to associate callback with
    343: . func - local Jacobian evaluation
    344: - ctx  - optional context for local Jacobian evaluation

    346:   Calling sequence of `func`:
    347: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
    348: . x    - state vector at which to evaluate Jacobian
    349: . J    - the Jacobian
    350: . M    - approximate Jacobian from which the preconditioner will be computed, often `J`
    351: - ctx  - optional context passed above

    353:   Level: beginner

    355: .seealso: [](ch_snes), `DMDA`, `DMDASNESSetJacobianLocal()`, `DMDASNESSetFunctionLocalVec()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
    356: @*/
    357: PetscErrorCode DMDASNESSetJacobianLocalVec(DM dm, PetscErrorCode (*func)(DMDALocalInfo *info, Vec x, Mat J, Mat M, void *), void *ctx)
    358: {
    359:   DMSNES     sdm;
    360:   DMSNES_DA *dmdasnes;

    362:   PetscFunctionBegin;
    364:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    365:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    367:   dmdasnes->jacobianlocalvec = func;
    368:   dmdasnes->jacobianlocalctx = ctx;

    370:   PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
    371:   PetscFunctionReturn(PETSC_SUCCESS);
    372: }

    374: /*@C
    375:   DMDASNESSetObjectiveLocal - set a local residual evaluation function to used with a `DMDA`

    377:   Logically Collective

    379:   Input Parameters:
    380: + dm   - `DM` to associate callback with
    381: . func - local objective evaluation, see `DMDASNESSetObjectiveLocal` for the calling sequence
    382: - ctx  - optional context for local residual evaluation

    384:   Calling sequence of `func`:
    385: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
    386: . x    - dimensional pointer to state at which to evaluate the objective (e.g. PetscScalar *x or **x or ***x)
    387: . obj  - returned objective value for the local subdomain
    388: - ctx  - optional context passed above

    390:   Level: beginner

    392: .seealso: [](ch_snes), `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocal()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`, `DMDASNESObjectiveFn`
    393: @*/
    394: PetscErrorCode DMDASNESSetObjectiveLocal(DM dm, PetscErrorCode (*func)(DMDALocalInfo *info, void *x, PetscReal *obj, void *), void *ctx)
    395: {
    396:   DMSNES     sdm;
    397:   DMSNES_DA *dmdasnes;

    399:   PetscFunctionBegin;
    401:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    402:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    404:   dmdasnes->objectivelocal    = func;
    405:   dmdasnes->objectivelocalctx = ctx;

    407:   PetscCall(DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes));
    408:   PetscFunctionReturn(PETSC_SUCCESS);
    409: }

    411: /*@C
    412:   DMDASNESSetObjectiveLocalVec - set a local residual evaluation function that operates on a local vector with `DMDA`

    414:   Logically Collective

    416:   Input Parameters:
    417: + dm   - `DM` to associate callback with
    418: . func - local objective evaluation, see `DMDASNESSetObjectiveLocalVec` for the calling sequence
    419: - ctx  - optional context for local residual evaluation

    421:   Calling sequence of `func`:
    422: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
    423: . x    - state vector at which to evaluate the objective
    424: . obj  - returned objective value for the local subdomain
    425: - ctx  - optional context passed above

    427:   Level: beginner

    429: .seealso: [](ch_snes), `DMDA`, `DMDASNESSetObjectiveLocal()`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocalVec()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`, `DMDASNESObjectiveVecFn`
    430: @*/
    431: PetscErrorCode DMDASNESSetObjectiveLocalVec(DM dm, PetscErrorCode (*func)(DMDALocalInfo *info, Vec x, PetscReal *obj, void *), void *ctx)
    432: {
    433:   DMSNES     sdm;
    434:   DMSNES_DA *dmdasnes;

    436:   PetscFunctionBegin;
    438:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    439:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    441:   dmdasnes->objectivelocalvec = func;
    442:   dmdasnes->objectivelocalctx = ctx;

    444:   PetscCall(DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes));
    445:   PetscFunctionReturn(PETSC_SUCCESS);
    446: }

    448: static PetscErrorCode SNESComputePicard_DMDA(SNES snes, Vec X, Vec F, void *ctx)
    449: {
    450:   DM            dm;
    451:   DMSNES_DA    *dmdasnes = (DMSNES_DA *)ctx;
    452:   DMDALocalInfo info;
    453:   Vec           Xloc;
    454:   void         *x, *f;

    456:   PetscFunctionBegin;
    460:   PetscCheck(dmdasnes->rhsplocal, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
    461:   PetscCall(SNESGetDM(snes, &dm));
    462:   PetscCall(DMGetLocalVector(dm, &Xloc));
    463:   PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
    464:   PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
    465:   PetscCall(DMDAGetLocalInfo(dm, &info));
    466:   PetscCall(DMDAVecGetArray(dm, Xloc, &x));
    467:   switch (dmdasnes->residuallocalimode) {
    468:   case INSERT_VALUES: {
    469:     PetscCall(DMDAVecGetArray(dm, F, &f));
    470:     PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
    471:     PetscCall(DMDAVecRestoreArray(dm, F, &f));
    472:   } break;
    473:   case ADD_VALUES: {
    474:     Vec Floc;
    475:     PetscCall(DMGetLocalVector(dm, &Floc));
    476:     PetscCall(VecZeroEntries(Floc));
    477:     PetscCall(DMDAVecGetArray(dm, Floc, &f));
    478:     PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
    479:     PetscCall(DMDAVecRestoreArray(dm, Floc, &f));
    480:     PetscCall(VecZeroEntries(F));
    481:     PetscCall(DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F));
    482:     PetscCall(DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F));
    483:     PetscCall(DMRestoreLocalVector(dm, &Floc));
    484:   } break;
    485:   default:
    486:     SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
    487:   }
    488:   PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
    489:   PetscCall(DMRestoreLocalVector(dm, &Xloc));
    490:   PetscFunctionReturn(PETSC_SUCCESS);
    491: }

    493: static PetscErrorCode SNESComputePicardJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
    494: {
    495:   DM            dm;
    496:   DMSNES_DA    *dmdasnes = (DMSNES_DA *)ctx;
    497:   DMDALocalInfo info;
    498:   Vec           Xloc;
    499:   void         *x;

    501:   PetscFunctionBegin;
    502:   PetscCheck(dmdasnes->jacobianplocal, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
    503:   PetscCall(SNESGetDM(snes, &dm));

    505:   PetscCall(DMGetLocalVector(dm, &Xloc));
    506:   PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
    507:   PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
    508:   PetscCall(DMDAGetLocalInfo(dm, &info));
    509:   PetscCall(DMDAVecGetArray(dm, Xloc, &x));
    510:   PetscCallBack("SNES Picard DMDA local callback Jacobian", (*dmdasnes->jacobianplocal)(&info, x, A, B, dmdasnes->picardlocalctx));
    511:   PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
    512:   PetscCall(DMRestoreLocalVector(dm, &Xloc));
    513:   if (A != B) {
    514:     PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
    515:     PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
    516:   }
    517:   PetscFunctionReturn(PETSC_SUCCESS);
    518: }

    520: /*@C
    521:   DMDASNESSetPicardLocal - set a local right-hand side and matrix evaluation function for Picard iteration with `DMDA`

    523:   Logically Collective

    525:   Input Parameters:
    526: + dm    - `DM` to associate callback with
    527: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
    528: . func  - local residual evaluation
    529: . jac   - function to compute Jacobian
    530: - ctx   - optional context for local residual evaluation

    532:   Calling sequence of `func`:
    533: + info - defines the subdomain to evaluate the residual on
    534: . x    - dimensional pointer to state at which to evaluate residual
    535: . f    - dimensional pointer to residual, write the residual here
    536: - ctx  - optional context passed above

    538:   Calling sequence of `jac`:
    539: + info - defines the subdomain to evaluate the residual on
    540: . x    - dimensional pointer to state at which to evaluate residual
    541: . jac  - the Jacobian
    542: . Jp   - approximation to the Jacobian used to compute the preconditioner, often `J`
    543: - ctx  - optional context passed above

    545:   Level: beginner

    547:   Note:
    548:   The user must use `SNESSetFunction`(`snes`,`NULL`,`SNESPicardComputeFunction`,&user));
    549:   in their code before calling this routine.

    551: .seealso: [](ch_snes), `SNES`, `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
    552: @*/
    553: PetscErrorCode DMDASNESSetPicardLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *info, void *x, void *f, void *ctx), PetscErrorCode (*jac)(DMDALocalInfo *info, void *x, Mat jac, Mat Jp, void *ctx), void *ctx)
    554: {
    555:   DMSNES     sdm;
    556:   DMSNES_DA *dmdasnes;

    558:   PetscFunctionBegin;
    560:   PetscCall(DMGetDMSNESWrite(dm, &sdm));
    561:   PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));

    563:   dmdasnes->residuallocalimode = imode;
    564:   dmdasnes->rhsplocal          = func;
    565:   dmdasnes->jacobianplocal     = jac;
    566:   dmdasnes->picardlocalctx     = ctx;

    568:   PetscCall(DMSNESSetPicard(dm, SNESComputePicard_DMDA, SNESComputePicardJacobian_DMDA, dmdasnes));
    569:   PetscCall(DMSNESSetMFFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
    570:   PetscFunctionReturn(PETSC_SUCCESS);
    571: }