Actual source code: nasm.c

  1: #include <petsc/private/snesimpl.h>
  2: #include <petscdm.h>

  4: typedef struct {
  5:   PetscInt    n;             /* local subdomains */
  6:   SNES       *subsnes;       /* nonlinear solvers for each subdomain */
  7:   Vec        *x;             /* solution vectors */
  8:   Vec        *xl;            /* solution local vectors */
  9:   Vec        *y;             /* step vectors */
 10:   Vec        *b;             /* rhs vectors */
 11:   Vec         weight;        /* weighting for adding updates on overlaps, in global space */
 12:   VecScatter *oscatter;      /* scatter from global space to the subdomain global space */
 13:   VecScatter *oscatter_copy; /* copy of the above */
 14:   VecScatter *iscatter;      /* scatter from global space to the nonoverlapping subdomain space */
 15:   VecScatter *gscatter;      /* scatter from global space to the subdomain local space */
 16:   PCASMType   type;          /* ASM type */
 17:   PetscBool   usesdm;        /* use the DM for setting up the subproblems */
 18:   PetscBool   finaljacobian; /* compute the jacobian of the converged solution */
 19:   PetscReal   damping;       /* damping parameter for updates from the blocks */
 20:   PetscBool   weight_set;    /* use a weight in the overlap updates */

 22:   /* logging events */
 23:   PetscLogEvent eventrestrictinterp;
 24:   PetscLogEvent eventsubsolve;

 26:   PetscInt fjtype; /* type of computed jacobian */
 27:   Vec      xinit;  /* initial solution in case the final jacobian type is computed as first */
 28: } SNES_NASM;

 30: const char *const SNESNASMTypes[]   = {"NONE", "RESTRICT", "INTERPOLATE", "BASIC", "PCASMType", "PC_ASM_", NULL};
 31: const char *const SNESNASMFJTypes[] = {"FINALOUTER", "FINALINNER", "INITIAL"};

 33: static PetscErrorCode SNESReset_NASM(SNES snes)
 34: {
 35:   SNES_NASM *nasm = (SNES_NASM *)snes->data;
 36:   PetscInt   i;

 38:   PetscFunctionBegin;
 39:   for (i = 0; i < nasm->n; i++) {
 40:     if (nasm->xl) PetscCall(VecDestroy(&nasm->xl[i]));
 41:     if (nasm->x) PetscCall(VecDestroy(&nasm->x[i]));
 42:     if (nasm->y) PetscCall(VecDestroy(&nasm->y[i]));
 43:     if (nasm->b) PetscCall(VecDestroy(&nasm->b[i]));

 45:     if (nasm->subsnes) PetscCall(SNESDestroy(&nasm->subsnes[i]));
 46:     if (nasm->oscatter) PetscCall(VecScatterDestroy(&nasm->oscatter[i]));
 47:     if (nasm->oscatter_copy) PetscCall(VecScatterDestroy(&nasm->oscatter_copy[i]));
 48:     if (nasm->iscatter) PetscCall(VecScatterDestroy(&nasm->iscatter[i]));
 49:     if (nasm->gscatter) PetscCall(VecScatterDestroy(&nasm->gscatter[i]));
 50:   }

 52:   PetscCall(PetscFree(nasm->x));
 53:   PetscCall(PetscFree(nasm->xl));
 54:   PetscCall(PetscFree(nasm->y));
 55:   PetscCall(PetscFree(nasm->b));

 57:   if (nasm->xinit) PetscCall(VecDestroy(&nasm->xinit));

 59:   PetscCall(PetscFree(nasm->subsnes));
 60:   PetscCall(PetscFree(nasm->oscatter));
 61:   PetscCall(PetscFree(nasm->oscatter_copy));
 62:   PetscCall(PetscFree(nasm->iscatter));
 63:   PetscCall(PetscFree(nasm->gscatter));

 65:   if (nasm->weight_set) PetscCall(VecDestroy(&nasm->weight));

 67:   nasm->eventrestrictinterp = 0;
 68:   nasm->eventsubsolve       = 0;
 69:   PetscFunctionReturn(PETSC_SUCCESS);
 70: }

 72: static PetscErrorCode SNESDestroy_NASM(SNES snes)
 73: {
 74:   PetscFunctionBegin;
 75:   PetscCall(SNESReset_NASM(snes));
 76:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetType_C", NULL));
 77:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetType_C", NULL));
 78:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetSubdomains_C", NULL));
 79:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetSubdomains_C", NULL));
 80:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetDamping_C", NULL));
 81:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetDamping_C", NULL));
 82:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetSubdomainVecs_C", NULL));
 83:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetComputeFinalJacobian_C", NULL));
 84:   PetscCall(PetscFree(snes->data));
 85:   PetscFunctionReturn(PETSC_SUCCESS);
 86: }

 88: static PetscErrorCode DMGlobalToLocalSubDomainDirichletHook_Private(DM dm, Vec g, InsertMode mode, Vec l, void *ctx)
 89: {
 90:   Vec bcs = (Vec)ctx;

 92:   PetscFunctionBegin;
 93:   PetscCall(VecCopy(bcs, l));
 94:   PetscFunctionReturn(PETSC_SUCCESS);
 95: }

 97: static PetscErrorCode SNESSetUp_NASM(SNES snes)
 98: {
 99:   SNES_NASM  *nasm = (SNES_NASM *)snes->data;
100:   DM          dm, subdm;
101:   DM         *subdms;
102:   PetscInt    i;
103:   const char *optionsprefix;
104:   Vec         F;
105:   PetscMPIInt size;
106:   KSP         ksp;
107:   PC          pc;

109:   PetscFunctionBegin;
110:   if (!nasm->subsnes) {
111:     PetscCall(SNESGetDM(snes, &dm));
112:     if (dm) {
113:       nasm->usesdm = PETSC_TRUE;
114:       PetscCall(DMCreateDomainDecomposition(dm, &nasm->n, NULL, NULL, NULL, &subdms));
115:       PetscCheck(subdms, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONGSTATE, "DM has no default decomposition defined.  Set subsolves manually with SNESNASMSetSubdomains().");
116:       PetscCall(DMCreateDomainDecompositionScatters(dm, nasm->n, subdms, &nasm->iscatter, &nasm->oscatter, &nasm->gscatter));
117:       PetscCall(PetscMalloc1(nasm->n, &nasm->oscatter_copy));
118:       for (i = 0; i < nasm->n; i++) PetscCall(VecScatterCopy(nasm->oscatter[i], &nasm->oscatter_copy[i]));

120:       PetscCall(SNESGetOptionsPrefix(snes, &optionsprefix));
121:       PetscCall(PetscMalloc1(nasm->n, &nasm->subsnes));
122:       for (i = 0; i < nasm->n; i++) {
123:         PetscCall(SNESCreate(PETSC_COMM_SELF, &nasm->subsnes[i]));
124:         PetscCall(PetscObjectIncrementTabLevel((PetscObject)nasm->subsnes[i], (PetscObject)snes, 1));
125:         PetscCall(SNESAppendOptionsPrefix(nasm->subsnes[i], optionsprefix));
126:         PetscCall(SNESAppendOptionsPrefix(nasm->subsnes[i], "sub_"));
127:         PetscCall(SNESSetDM(nasm->subsnes[i], subdms[i]));
128:         PetscCallMPI(MPI_Comm_size(PetscObjectComm((PetscObject)nasm->subsnes[i]), &size));
129:         if (size == 1) {
130:           PetscCall(SNESGetKSP(nasm->subsnes[i], &ksp));
131:           PetscCall(KSPGetPC(ksp, &pc));
132:           PetscCall(KSPSetType(ksp, KSPPREONLY));
133:           PetscCall(PCSetType(pc, PCLU));
134:         }
135:         PetscCall(SNESSetFromOptions(nasm->subsnes[i]));
136:         PetscCall(DMDestroy(&subdms[i]));
137:       }
138:       PetscCall(PetscFree(subdms));
139:     } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Cannot construct local problems automatically without a DM!");
140:   } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Must set subproblems manually if there is no DM!");
141:   /* allocate the global vectors */
142:   if (!nasm->x) PetscCall(PetscCalloc1(nasm->n, &nasm->x));
143:   if (!nasm->xl) PetscCall(PetscCalloc1(nasm->n, &nasm->xl));
144:   if (!nasm->y) PetscCall(PetscCalloc1(nasm->n, &nasm->y));
145:   if (!nasm->b) PetscCall(PetscCalloc1(nasm->n, &nasm->b));

147:   for (i = 0; i < nasm->n; i++) {
148:     PetscCall(SNESGetFunction(nasm->subsnes[i], &F, NULL, NULL));
149:     if (!nasm->x[i]) PetscCall(VecDuplicate(F, &nasm->x[i]));
150:     if (!nasm->y[i]) PetscCall(VecDuplicate(F, &nasm->y[i]));
151:     if (!nasm->b[i]) PetscCall(VecDuplicate(F, &nasm->b[i]));
152:     if (!nasm->xl[i]) {
153:       PetscCall(SNESGetDM(nasm->subsnes[i], &subdm));
154:       PetscCall(DMCreateLocalVector(subdm, &nasm->xl[i]));
155:       PetscCall(DMGlobalToLocalHookAdd(subdm, DMGlobalToLocalSubDomainDirichletHook_Private, NULL, nasm->xl[i]));
156:     }
157:   }
158:   if (nasm->finaljacobian) {
159:     PetscCall(SNESSetUpMatrices(snes));
160:     if (nasm->fjtype == 2) PetscCall(VecDuplicate(snes->vec_sol, &nasm->xinit));
161:     for (i = 0; i < nasm->n; i++) PetscCall(SNESSetUpMatrices(nasm->subsnes[i]));
162:   }
163:   PetscFunctionReturn(PETSC_SUCCESS);
164: }

166: static PetscErrorCode SNESSetFromOptions_NASM(SNES snes, PetscOptionItems *PetscOptionsObject)
167: {
168:   PCASMType  asmtype;
169:   PetscBool  flg, monflg;
170:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

172:   PetscFunctionBegin;
173:   PetscOptionsHeadBegin(PetscOptionsObject, "Nonlinear Additive Schwarz options");
174:   PetscCall(PetscOptionsEnum("-snes_nasm_type", "Type of restriction/extension", "", SNESNASMTypes, (PetscEnum)nasm->type, (PetscEnum *)&asmtype, &flg));
175:   if (flg) PetscCall(SNESNASMSetType(snes, asmtype));
176:   flg    = PETSC_FALSE;
177:   monflg = PETSC_TRUE;
178:   PetscCall(PetscOptionsReal("-snes_nasm_damping", "The new solution is obtained as old solution plus dmp times (sum of the solutions on the subdomains)", "SNESNASMSetDamping", nasm->damping, &nasm->damping, &flg));
179:   if (flg) PetscCall(SNESNASMSetDamping(snes, nasm->damping));
180:   PetscCall(PetscOptionsDeprecated("-snes_nasm_sub_view", NULL, "3.15", "Use -snes_view ::ascii_info_detail"));
181:   PetscCall(PetscOptionsBool("-snes_nasm_finaljacobian", "Compute the global jacobian of the final iterate (for ASPIN)", "", nasm->finaljacobian, &nasm->finaljacobian, NULL));
182:   PetscCall(PetscOptionsEList("-snes_nasm_finaljacobian_type", "The type of the final jacobian computed.", "", SNESNASMFJTypes, 3, SNESNASMFJTypes[0], &nasm->fjtype, NULL));
183:   PetscCall(PetscOptionsBool("-snes_nasm_log", "Log times for subSNES solves and restriction", "", monflg, &monflg, &flg));
184:   if (flg) {
185:     PetscCall(PetscLogEventRegister("SNESNASMSubSolve", ((PetscObject)snes)->classid, &nasm->eventsubsolve));
186:     PetscCall(PetscLogEventRegister("SNESNASMRestrict", ((PetscObject)snes)->classid, &nasm->eventrestrictinterp));
187:   }
188:   PetscOptionsHeadEnd();
189:   PetscFunctionReturn(PETSC_SUCCESS);
190: }

192: static PetscErrorCode SNESView_NASM(SNES snes, PetscViewer viewer)
193: {
194:   SNES_NASM        *nasm = (SNES_NASM *)snes->data;
195:   PetscMPIInt       rank, size;
196:   PetscInt          i, N, bsz;
197:   PetscBool         iascii, isstring;
198:   PetscViewer       sviewer;
199:   MPI_Comm          comm;
200:   PetscViewerFormat format;
201:   const char       *prefix;

203:   PetscFunctionBegin;
204:   PetscCall(PetscObjectGetComm((PetscObject)snes, &comm));
205:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
206:   PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERSTRING, &isstring));
207:   PetscCallMPI(MPI_Comm_rank(comm, &rank));
208:   PetscCallMPI(MPI_Comm_size(comm, &size));
209:   PetscCall(MPIU_Allreduce(&nasm->n, &N, 1, MPIU_INT, MPI_SUM, comm));
210:   if (iascii) {
211:     PetscCall(PetscViewerASCIIPrintf(viewer, "  total subdomain blocks = %" PetscInt_FMT "\n", N));
212:     PetscCall(PetscViewerGetFormat(viewer, &format));
213:     if (format != PETSC_VIEWER_ASCII_INFO_DETAIL) {
214:       if (nasm->subsnes) {
215:         PetscCall(PetscViewerASCIIPrintf(viewer, "  Local solver information for first block on rank 0:\n"));
216:         PetscCall(SNESGetOptionsPrefix(snes, &prefix));
217:         PetscCall(PetscViewerASCIIPrintf(viewer, "  Use -%ssnes_view ::ascii_info_detail to display information for all blocks\n", prefix ? prefix : ""));
218:         PetscCall(PetscViewerASCIIPushTab(viewer));
219:         PetscCall(PetscViewerGetSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
220:         if (rank == 0) {
221:           PetscCall(PetscViewerASCIIPushTab(viewer));
222:           PetscCall(SNESView(nasm->subsnes[0], sviewer));
223:           PetscCall(PetscViewerASCIIPopTab(viewer));
224:         }
225:         PetscCall(PetscViewerRestoreSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
226:         PetscCall(PetscViewerASCIIPopTab(viewer));
227:       }
228:     } else {
229:       /* print the solver on each block */
230:       PetscCall(PetscViewerASCIIPushSynchronized(viewer));
231:       PetscCall(PetscViewerASCIISynchronizedPrintf(viewer, "  [%d] number of local blocks = %" PetscInt_FMT "\n", (int)rank, nasm->n));
232:       PetscCall(PetscViewerFlush(viewer));
233:       PetscCall(PetscViewerASCIIPopSynchronized(viewer));
234:       PetscCall(PetscViewerASCIIPrintf(viewer, "  Local solver information for each block is in the following SNES objects:\n"));
235:       PetscCall(PetscViewerASCIIPushTab(viewer));
236:       PetscCall(PetscViewerASCIIPrintf(viewer, "- - - - - - - - - - - - - - - - - -\n"));
237:       PetscCall(PetscViewerGetSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
238:       for (i = 0; i < nasm->n; i++) {
239:         PetscCall(VecGetLocalSize(nasm->x[i], &bsz));
240:         PetscCall(PetscViewerASCIIPrintf(sviewer, "[%d] local block number %" PetscInt_FMT ", size = %" PetscInt_FMT "\n", (int)rank, i, bsz));
241:         PetscCall(SNESView(nasm->subsnes[i], sviewer));
242:         PetscCall(PetscViewerASCIIPrintf(sviewer, "- - - - - - - - - - - - - - - - - -\n"));
243:       }
244:       PetscCall(PetscViewerRestoreSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
245:       PetscCall(PetscViewerFlush(viewer));
246:       PetscCall(PetscViewerASCIIPopTab(viewer));
247:     }
248:   } else if (isstring) {
249:     PetscCall(PetscViewerStringSPrintf(viewer, " blocks=%" PetscInt_FMT ",type=%s", N, SNESNASMTypes[nasm->type]));
250:     PetscCall(PetscViewerGetSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
251:     if (nasm->subsnes && rank == 0) PetscCall(SNESView(nasm->subsnes[0], sviewer));
252:     PetscCall(PetscViewerRestoreSubViewer(viewer, PETSC_COMM_SELF, &sviewer));
253:   }
254:   PetscFunctionReturn(PETSC_SUCCESS);
255: }

257: /*@
258:   SNESNASMSetType - Set the type of subdomain update used for the nonlinear additive Schwarz solver `SNESNASM`

260:   Logically Collective

262:   Input Parameters:
263: + snes - the `SNES` context
264: - type - the type of update, `PC_ASM_BASIC` or `PC_ASM_RESTRICT`

266:   Options Database Key:
267: . -snes_nasm_type <basic,restrict> - type of subdomain update used

269:   Level: intermediate

271: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMGetType()`, `PCASMSetType()`, `PC_ASM_BASIC`, `PC_ASM_RESTRICT`, `PCASMType`
272: @*/
273: PetscErrorCode SNESNASMSetType(SNES snes, PCASMType type)
274: {
275:   PetscErrorCode (*f)(SNES, PCASMType);

277:   PetscFunctionBegin;
278:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMSetType_C", &f));
279:   if (f) PetscCall((f)(snes, type));
280:   PetscFunctionReturn(PETSC_SUCCESS);
281: }

283: static PetscErrorCode SNESNASMSetType_NASM(SNES snes, PCASMType type)
284: {
285:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

287:   PetscFunctionBegin;
288:   PetscCheck(type == PC_ASM_BASIC || type == PC_ASM_RESTRICT, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "SNESNASM only supports basic and restrict types");
289:   nasm->type = type;
290:   PetscFunctionReturn(PETSC_SUCCESS);
291: }

293: /*@
294:   SNESNASMGetType - Get the type of subdomain update used for the nonlinear additive Schwarz solver `SNESNASM`

296:   Logically Collective

298:   Input Parameter:
299: . snes - the `SNES` context

301:   Output Parameter:
302: . type - the type of update

304:   Level: intermediate

306: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMSetType()`, `PCASMGetType()`, `PC_ASM_BASIC`, `PC_ASM_RESTRICT`, `PCASMType`
307: @*/
308: PetscErrorCode SNESNASMGetType(SNES snes, PCASMType *type)
309: {
310:   PetscFunctionBegin;
311:   PetscUseMethod(snes, "SNESNASMGetType_C", (SNES, PCASMType *), (snes, type));
312:   PetscFunctionReturn(PETSC_SUCCESS);
313: }

315: static PetscErrorCode SNESNASMGetType_NASM(SNES snes, PCASMType *type)
316: {
317:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

319:   PetscFunctionBegin;
320:   *type = nasm->type;
321:   PetscFunctionReturn(PETSC_SUCCESS);
322: }

324: /*@
325:   SNESNASMSetSubdomains - Manually Set the context required to restrict and solve subdomain problems in the nonlinear additive Schwarz solver

327:   Logically Collective

329:   Input Parameters:
330: + snes     - the `SNES` context
331: . n        - the number of local subdomains
332: . subsnes  - solvers defined on the local subdomains
333: . iscatter - scatters into the nonoverlapping portions of the local subdomains
334: . oscatter - scatters into the overlapping portions of the local subdomains
335: - gscatter - scatters into the (ghosted) local vector of the local subdomain

337:   Level: intermediate

339: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMGetSubdomains()`
340: @*/
341: PetscErrorCode SNESNASMSetSubdomains(SNES snes, PetscInt n, SNES subsnes[], VecScatter iscatter[], VecScatter oscatter[], VecScatter gscatter[])
342: {
343:   PetscErrorCode (*f)(SNES, PetscInt, SNES *, VecScatter *, VecScatter *, VecScatter *);

345:   PetscFunctionBegin;
346:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMSetSubdomains_C", &f));
347:   if (f) PetscCall((f)(snes, n, subsnes, iscatter, oscatter, gscatter));
348:   PetscFunctionReturn(PETSC_SUCCESS);
349: }

351: static PetscErrorCode SNESNASMSetSubdomains_NASM(SNES snes, PetscInt n, SNES subsnes[], VecScatter iscatter[], VecScatter oscatter[], VecScatter gscatter[])
352: {
353:   PetscInt   i;
354:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

356:   PetscFunctionBegin;
357:   PetscCheck(!snes->setupcalled, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNESNASMSetSubdomains() should be called before calling SNESSetUp().");

359:   /* tear down the previously set things */
360:   PetscCall(SNESReset(snes));

362:   nasm->n = n;
363:   if (oscatter) {
364:     for (i = 0; i < n; i++) PetscCall(PetscObjectReference((PetscObject)oscatter[i]));
365:   }
366:   if (iscatter) {
367:     for (i = 0; i < n; i++) PetscCall(PetscObjectReference((PetscObject)iscatter[i]));
368:   }
369:   if (gscatter) {
370:     for (i = 0; i < n; i++) PetscCall(PetscObjectReference((PetscObject)gscatter[i]));
371:   }
372:   if (oscatter) {
373:     PetscCall(PetscMalloc1(n, &nasm->oscatter));
374:     PetscCall(PetscMalloc1(n, &nasm->oscatter_copy));
375:     for (i = 0; i < n; i++) {
376:       nasm->oscatter[i] = oscatter[i];
377:       PetscCall(VecScatterCopy(oscatter[i], &nasm->oscatter_copy[i]));
378:     }
379:   }
380:   if (iscatter) {
381:     PetscCall(PetscMalloc1(n, &nasm->iscatter));
382:     for (i = 0; i < n; i++) nasm->iscatter[i] = iscatter[i];
383:   }
384:   if (gscatter) {
385:     PetscCall(PetscMalloc1(n, &nasm->gscatter));
386:     for (i = 0; i < n; i++) nasm->gscatter[i] = gscatter[i];
387:   }

389:   if (subsnes) {
390:     PetscCall(PetscMalloc1(n, &nasm->subsnes));
391:     for (i = 0; i < n; i++) nasm->subsnes[i] = subsnes[i];
392:   }
393:   PetscFunctionReturn(PETSC_SUCCESS);
394: }

396: /*@
397:   SNESNASMGetSubdomains - Get the local subdomain contexts for the nonlinear additive Schwarz solver

399:   Not Collective but some of the objects returned will be parallel

401:   Input Parameter:
402: . snes - the `SNES` context

404:   Output Parameters:
405: + n        - the number of local subdomains
406: . subsnes  - solvers defined on the local subdomains
407: . iscatter - scatters into the nonoverlapping portions of the local subdomains
408: . oscatter - scatters into the overlapping portions of the local subdomains
409: - gscatter - scatters into the (ghosted) local vector of the local subdomain

411:   Level: intermediate

413: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMSetSubdomains()`
414: @*/
415: PetscErrorCode SNESNASMGetSubdomains(SNES snes, PetscInt *n, SNES *subsnes[], VecScatter *iscatter[], VecScatter *oscatter[], VecScatter *gscatter[])
416: {
417:   PetscErrorCode (*f)(SNES, PetscInt *, SNES **, VecScatter **, VecScatter **, VecScatter **);

419:   PetscFunctionBegin;
420:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMGetSubdomains_C", &f));
421:   if (f) PetscCall((f)(snes, n, subsnes, iscatter, oscatter, gscatter));
422:   PetscFunctionReturn(PETSC_SUCCESS);
423: }

425: static PetscErrorCode SNESNASMGetSubdomains_NASM(SNES snes, PetscInt *n, SNES *subsnes[], VecScatter *iscatter[], VecScatter *oscatter[], VecScatter *gscatter[])
426: {
427:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

429:   PetscFunctionBegin;
430:   if (n) *n = nasm->n;
431:   if (oscatter) *oscatter = nasm->oscatter;
432:   if (iscatter) *iscatter = nasm->iscatter;
433:   if (gscatter) *gscatter = nasm->gscatter;
434:   if (subsnes) *subsnes = nasm->subsnes;
435:   PetscFunctionReturn(PETSC_SUCCESS);
436: }

438: /*@
439:   SNESNASMGetSubdomainVecs - Get the processor-local subdomain vectors for the nonlinear additive Schwarz solver

441:   Not Collective

443:   Input Parameter:
444: . snes - the `SNES` context

446:   Output Parameters:
447: + n  - the number of local subdomains
448: . x  - The subdomain solution vector
449: . y  - The subdomain step vector
450: . b  - The subdomain RHS vector
451: - xl - The subdomain local vectors (ghosted)

453:   Level: developer

455: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMGetSubdomains()`
456: @*/
457: PetscErrorCode SNESNASMGetSubdomainVecs(SNES snes, PetscInt *n, Vec **x, Vec **y, Vec **b, Vec **xl)
458: {
459:   PetscErrorCode (*f)(SNES, PetscInt *, Vec **, Vec **, Vec **, Vec **);

461:   PetscFunctionBegin;
462:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMGetSubdomainVecs_C", &f));
463:   if (f) PetscCall((f)(snes, n, x, y, b, xl));
464:   PetscFunctionReturn(PETSC_SUCCESS);
465: }

467: static PetscErrorCode SNESNASMGetSubdomainVecs_NASM(SNES snes, PetscInt *n, Vec **x, Vec **y, Vec **b, Vec **xl)
468: {
469:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

471:   PetscFunctionBegin;
472:   if (n) *n = nasm->n;
473:   if (x) *x = nasm->x;
474:   if (y) *y = nasm->y;
475:   if (b) *b = nasm->b;
476:   if (xl) *xl = nasm->xl;
477:   PetscFunctionReturn(PETSC_SUCCESS);
478: }

480: /*@
481:   SNESNASMSetComputeFinalJacobian - Schedules the computation of the global and subdomain Jacobians upon convergence for the
482:   nonlinear additive Schwarz solver

484:   Collective

486:   Input Parameters:
487: + snes - the SNES context
488: - flg  - `PETSC_TRUE` to compute the Jacobians

490:   Level: developer

492:   Notes:
493:   This is used almost exclusively in the implementation of `SNESASPIN`, where the converged subdomain and global Jacobian
494:   is needed at each linear iteration.

496: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMGetSubdomains()`
497: @*/
498: PetscErrorCode SNESNASMSetComputeFinalJacobian(SNES snes, PetscBool flg)
499: {
500:   PetscErrorCode (*f)(SNES, PetscBool);

502:   PetscFunctionBegin;
503:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMSetComputeFinalJacobian_C", &f));
504:   if (f) PetscCall((f)(snes, flg));
505:   PetscFunctionReturn(PETSC_SUCCESS);
506: }

508: static PetscErrorCode SNESNASMSetComputeFinalJacobian_NASM(SNES snes, PetscBool flg)
509: {
510:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

512:   PetscFunctionBegin;
513:   nasm->finaljacobian = flg;
514:   PetscFunctionReturn(PETSC_SUCCESS);
515: }

517: /*@
518:   SNESNASMSetDamping - Sets the update damping for `SNESNASM` the nonlinear additive Schwarz solver

520:   Logically Collective

522:   Input Parameters:
523: + snes - the `SNES` context
524: - dmp  - damping

526:   Options Database Key:
527: . -snes_nasm_damping <dmp> - the new solution is obtained as old solution plus `dmp` times (sum of the solutions on the subdomains)

529:   Level: intermediate

531:   Note:
532:   The new solution is obtained as old solution plus dmp times (sum of the solutions on the subdomains)

534: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMGetDamping()`
535: @*/
536: PetscErrorCode SNESNASMSetDamping(SNES snes, PetscReal dmp)
537: {
538:   PetscErrorCode (*f)(SNES, PetscReal);

540:   PetscFunctionBegin;
541:   PetscCall(PetscObjectQueryFunction((PetscObject)snes, "SNESNASMSetDamping_C", (void (**)(void)) & f));
542:   if (f) PetscCall((f)(snes, dmp));
543:   PetscFunctionReturn(PETSC_SUCCESS);
544: }

546: static PetscErrorCode SNESNASMSetDamping_NASM(SNES snes, PetscReal dmp)
547: {
548:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

550:   PetscFunctionBegin;
551:   nasm->damping = dmp;
552:   PetscFunctionReturn(PETSC_SUCCESS);
553: }

555: /*@
556:   SNESNASMGetDamping - Gets the update damping for `SNESNASM` the nonlinear additive Schwarz solver

558:   Not Collective

560:   Input Parameter:
561: . snes - the `SNES` context

563:   Output Parameter:
564: . dmp - damping

566:   Level: intermediate

568: .seealso: [](ch_snes), `SNES`, `SNESNASM`, `SNESNASMSetDamping()`
569: @*/
570: PetscErrorCode SNESNASMGetDamping(SNES snes, PetscReal *dmp)
571: {
572:   PetscFunctionBegin;
573:   PetscUseMethod(snes, "SNESNASMGetDamping_C", (SNES, PetscReal *), (snes, dmp));
574:   PetscFunctionReturn(PETSC_SUCCESS);
575: }

577: static PetscErrorCode SNESNASMGetDamping_NASM(SNES snes, PetscReal *dmp)
578: {
579:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

581:   PetscFunctionBegin;
582:   *dmp = nasm->damping;
583:   PetscFunctionReturn(PETSC_SUCCESS);
584: }

586: /*
587:   Input Parameters:
588: + snes - The solver
589: . B - The RHS vector
590: - X - The initial guess

592:   Output Parameter:
593: . Y - The solution update

595:   TODO: All scatters should be packed into one
596: */
597: static PetscErrorCode SNESNASMSolveLocal_Private(SNES snes, Vec B, Vec Y, Vec X)
598: {
599:   SNES_NASM *nasm = (SNES_NASM *)snes->data;
600:   SNES       subsnes;
601:   PetscInt   i;
602:   PetscReal  dmp;
603:   Vec        Xl, Bl, Yl, Xlloc;
604:   VecScatter iscat, oscat, gscat, oscat_copy;
605:   DM         dm, subdm;
606:   PCASMType  type;

608:   PetscFunctionBegin;
609:   PetscCall(SNESNASMGetType(snes, &type));
610:   PetscCall(SNESGetDM(snes, &dm));
611:   PetscCall(VecSet(Y, 0));
612:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventBegin(nasm->eventrestrictinterp, snes, 0, 0, 0));
613:   for (i = 0; i < nasm->n; i++) {
614:     /* scatter the solution to the global solution and the local solution */
615:     Xl         = nasm->x[i];
616:     Xlloc      = nasm->xl[i];
617:     oscat      = nasm->oscatter[i];
618:     oscat_copy = nasm->oscatter_copy[i];
619:     gscat      = nasm->gscatter[i];
620:     PetscCall(VecScatterBegin(oscat, X, Xl, INSERT_VALUES, SCATTER_FORWARD));
621:     PetscCall(VecScatterBegin(gscat, X, Xlloc, INSERT_VALUES, SCATTER_FORWARD));
622:     if (B) {
623:       /* scatter the RHS to the local RHS */
624:       Bl = nasm->b[i];
625:       PetscCall(VecScatterBegin(oscat_copy, B, Bl, INSERT_VALUES, SCATTER_FORWARD));
626:     }
627:   }
628:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventEnd(nasm->eventrestrictinterp, snes, 0, 0, 0));

630:   if (nasm->eventsubsolve) PetscCall(PetscLogEventBegin(nasm->eventsubsolve, snes, 0, 0, 0));
631:   for (i = 0; i < nasm->n; i++) {
632:     Xl      = nasm->x[i];
633:     Xlloc   = nasm->xl[i];
634:     Yl      = nasm->y[i];
635:     subsnes = nasm->subsnes[i];
636:     PetscCall(SNESGetDM(subsnes, &subdm));
637:     iscat      = nasm->iscatter[i];
638:     oscat      = nasm->oscatter[i];
639:     oscat_copy = nasm->oscatter_copy[i];
640:     gscat      = nasm->gscatter[i];
641:     PetscCall(VecScatterEnd(oscat, X, Xl, INSERT_VALUES, SCATTER_FORWARD));
642:     PetscCall(VecScatterEnd(gscat, X, Xlloc, INSERT_VALUES, SCATTER_FORWARD));
643:     if (B) {
644:       Bl = nasm->b[i];
645:       PetscCall(VecScatterEnd(oscat_copy, B, Bl, INSERT_VALUES, SCATTER_FORWARD));
646:     } else Bl = NULL;

648:     PetscCall(DMSubDomainRestrict(dm, oscat, gscat, subdm));
649:     PetscCall(VecCopy(Xl, Yl));
650:     PetscCall(SNESSolve(subsnes, Bl, Xl));
651:     PetscCall(VecAYPX(Yl, -1.0, Xl));
652:     PetscCall(VecScale(Yl, nasm->damping));
653:     if (type == PC_ASM_BASIC) {
654:       PetscCall(VecScatterBegin(oscat, Yl, Y, ADD_VALUES, SCATTER_REVERSE));
655:       PetscCall(VecScatterEnd(oscat, Yl, Y, ADD_VALUES, SCATTER_REVERSE));
656:     } else if (type == PC_ASM_RESTRICT) {
657:       PetscCall(VecScatterBegin(iscat, Yl, Y, ADD_VALUES, SCATTER_REVERSE));
658:       PetscCall(VecScatterEnd(iscat, Yl, Y, ADD_VALUES, SCATTER_REVERSE));
659:     } else SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "Only basic and restrict types are supported for SNESNASM");
660:   }
661:   if (nasm->eventsubsolve) PetscCall(PetscLogEventEnd(nasm->eventsubsolve, snes, 0, 0, 0));
662:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventBegin(nasm->eventrestrictinterp, snes, 0, 0, 0));
663:   if (nasm->weight_set) PetscCall(VecPointwiseMult(Y, Y, nasm->weight));
664:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventEnd(nasm->eventrestrictinterp, snes, 0, 0, 0));
665:   PetscCall(SNESNASMGetDamping(snes, &dmp));
666:   PetscCall(VecAXPY(X, dmp, Y));
667:   PetscFunctionReturn(PETSC_SUCCESS);
668: }

670: static PetscErrorCode SNESNASMComputeFinalJacobian_Private(SNES snes, Vec Xfinal)
671: {
672:   Vec        X    = Xfinal;
673:   SNES_NASM *nasm = (SNES_NASM *)snes->data;
674:   SNES       subsnes;
675:   PetscInt   i, lag = 1;
676:   Vec        Xlloc, Xl, Fl, F;
677:   VecScatter oscat, gscat;
678:   DM         dm, subdm;

680:   PetscFunctionBegin;
681:   if (nasm->fjtype == 2) X = nasm->xinit;
682:   F = snes->vec_func;
683:   if (snes->normschedule == SNES_NORM_NONE) PetscCall(SNESComputeFunction(snes, X, F));
684:   PetscCall(SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre));
685:   PetscCall(SNESGetDM(snes, &dm));
686:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventBegin(nasm->eventrestrictinterp, snes, 0, 0, 0));
687:   if (nasm->fjtype != 1) {
688:     for (i = 0; i < nasm->n; i++) {
689:       Xlloc = nasm->xl[i];
690:       gscat = nasm->gscatter[i];
691:       PetscCall(VecScatterBegin(gscat, X, Xlloc, INSERT_VALUES, SCATTER_FORWARD));
692:     }
693:   }
694:   if (nasm->eventrestrictinterp) PetscCall(PetscLogEventEnd(nasm->eventrestrictinterp, snes, 0, 0, 0));
695:   for (i = 0; i < nasm->n; i++) {
696:     Fl      = nasm->subsnes[i]->vec_func;
697:     Xl      = nasm->x[i];
698:     Xlloc   = nasm->xl[i];
699:     subsnes = nasm->subsnes[i];
700:     oscat   = nasm->oscatter[i];
701:     gscat   = nasm->gscatter[i];
702:     if (nasm->fjtype != 1) PetscCall(VecScatterEnd(gscat, X, Xlloc, INSERT_VALUES, SCATTER_FORWARD));
703:     PetscCall(SNESGetDM(subsnes, &subdm));
704:     PetscCall(DMSubDomainRestrict(dm, oscat, gscat, subdm));
705:     if (nasm->fjtype != 1) {
706:       PetscCall(DMLocalToGlobalBegin(subdm, Xlloc, INSERT_VALUES, Xl));
707:       PetscCall(DMLocalToGlobalEnd(subdm, Xlloc, INSERT_VALUES, Xl));
708:     }
709:     if (subsnes->lagjacobian == -1) subsnes->lagjacobian = -2;
710:     else if (subsnes->lagjacobian > 1) lag = subsnes->lagjacobian;
711:     PetscCall(SNESComputeFunction(subsnes, Xl, Fl));
712:     PetscCall(SNESComputeJacobian(subsnes, Xl, subsnes->jacobian, subsnes->jacobian_pre));
713:     if (lag > 1) subsnes->lagjacobian = lag;
714:   }
715:   PetscFunctionReturn(PETSC_SUCCESS);
716: }

718: static PetscErrorCode SNESSolve_NASM(SNES snes)
719: {
720:   Vec              F;
721:   Vec              X;
722:   Vec              B;
723:   Vec              Y;
724:   PetscInt         i;
725:   PetscReal        fnorm = 0.0;
726:   SNESNormSchedule normschedule;
727:   SNES_NASM       *nasm = (SNES_NASM *)snes->data;

729:   PetscFunctionBegin;

731:   PetscCheck(!snes->xl & !snes->xu && !snes->ops->computevariablebounds, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);

733:   PetscCall(PetscCitationsRegister(SNESCitation, &SNEScite));
734:   X = snes->vec_sol;
735:   Y = snes->vec_sol_update;
736:   F = snes->vec_func;
737:   B = snes->vec_rhs;

739:   PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
740:   snes->iter = 0;
741:   snes->norm = 0.;
742:   PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
743:   snes->reason = SNES_CONVERGED_ITERATING;
744:   PetscCall(SNESGetNormSchedule(snes, &normschedule));
745:   if (normschedule == SNES_NORM_ALWAYS || normschedule == SNES_NORM_INITIAL_ONLY || normschedule == SNES_NORM_INITIAL_FINAL_ONLY) {
746:     /* compute the initial function and preconditioned update delX */
747:     if (!snes->vec_func_init_set) {
748:       PetscCall(SNESComputeFunction(snes, X, F));
749:     } else snes->vec_func_init_set = PETSC_FALSE;

751:     PetscCall(VecNorm(F, NORM_2, &fnorm)); /* fnorm <- ||F||  */
752:     SNESCheckFunctionNorm(snes, fnorm);
753:     PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
754:     snes->iter = 0;
755:     snes->norm = fnorm;
756:     PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
757:     PetscCall(SNESLogConvergenceHistory(snes, snes->norm, 0));

759:     /* test convergence */
760:     PetscCall(SNESConverged(snes, 0, 0.0, 0.0, fnorm));
761:     PetscCall(SNESMonitor(snes, 0, snes->norm));
762:     if (snes->reason) PetscFunctionReturn(PETSC_SUCCESS);
763:   } else {
764:     PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
765:     PetscCall(SNESLogConvergenceHistory(snes, snes->norm, 0));
766:     PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
767:   }

769:   /* Call general purpose update function */
770:   PetscTryTypeMethod(snes, update, snes->iter);
771:   /* copy the initial solution over for later */
772:   if (nasm->fjtype == 2) PetscCall(VecCopy(X, nasm->xinit));

774:   for (i = 0; i < snes->max_its; i++) {
775:     PetscCall(SNESNASMSolveLocal_Private(snes, B, Y, X));
776:     if (normschedule == SNES_NORM_ALWAYS || ((i == snes->max_its - 1) && (normschedule == SNES_NORM_INITIAL_FINAL_ONLY || normschedule == SNES_NORM_FINAL_ONLY))) {
777:       PetscCall(SNESComputeFunction(snes, X, F));
778:       PetscCall(VecNorm(F, NORM_2, &fnorm)); /* fnorm <- ||F||  */
779:       SNESCheckFunctionNorm(snes, fnorm);
780:     }
781:     /* Monitor convergence */
782:     PetscCall(PetscObjectSAWsTakeAccess((PetscObject)snes));
783:     snes->iter = i + 1;
784:     snes->norm = fnorm;
785:     PetscCall(PetscObjectSAWsGrantAccess((PetscObject)snes));
786:     PetscCall(SNESLogConvergenceHistory(snes, snes->norm, 0));
787:     /* Test for convergence */
788:     PetscCall(SNESConverged(snes, snes->iter, 0.0, 0.0, fnorm));
789:     PetscCall(SNESMonitor(snes, snes->iter, snes->norm));
790:     if (snes->reason) break;
791:     /* Call general purpose update function */
792:     PetscTryTypeMethod(snes, update, snes->iter);
793:   }
794:   if (nasm->finaljacobian) {
795:     PetscCall(SNESNASMComputeFinalJacobian_Private(snes, X));
796:     SNESCheckJacobianDomainerror(snes);
797:   }
798:   PetscFunctionReturn(PETSC_SUCCESS);
799: }

801: /*MC
802:   SNESNASM - Nonlinear Additive Schwarz solver {cite}`ck02`, {cite}`bruneknepleysmithtu15`

804:    Options Database Keys:
805: +  -snes_nasm_log                                                - enable logging events for the communication and solve stages
806: .  -snes_nasm_type <basic,restrict>                              - type of subdomain update used
807: .  -snes_nasm_damping <dmp>                                      - the new solution is obtained as old solution plus dmp times (sum of the solutions on the subdomains)
808: .  -snes_nasm_finaljacobian                                      - compute the local and global Jacobians of the final iterate
809: .  -snes_nasm_finaljacobian_type <finalinner,finalouter,initial> - pick state the Jacobian is calculated at
810: .  -sub_snes_                                                    - options prefix of the subdomain nonlinear solves
811: .  -sub_ksp_                                                     - options prefix of the subdomain Krylov solver
812: -  -sub_pc_                                                      - options prefix of the subdomain preconditioner

814:    Level: advanced

816:    Note:
817:    This is not often used directly as a solver, it converges too slowly. However it works well as a nonlinear preconditioner for
818:    the `SNESASPIN` solver

820:    Developer Note:
821:    This is a non-Newton based nonlinear solver that does not directly require a Jacobian; hence the flag snes->usesksp is set to
822:    false and `SNESView()` and -snes_view do not display a `KSP` object. However, if the flag nasm->finaljacobian is set (for example, if
823:    `SNESNASM` is used as a nonlinear preconditioner for  `SNESASPIN`) then `SNESSetUpMatrices()` is called to generate the
824:    Jacobian (needed by `SNESASPIN`)
825:    and this utilizes the inner `KSP` object for storing the matrices, but the `KSP` is never used for solving a linear system. When `SNESNASM` is
826:    used by `SNESASPIN` they share the same Jacobian matrices because `SNESSetUp()` (called on the outer `SNESASPIN`) causes the inner `SNES`
827:    object (in this case `SNESNASM`) to inherit the outer Jacobian matrices.

829: .seealso: [](ch_snes), `SNESCreate()`, `SNES`, `SNESSetType()`, `SNESType`, `SNESNASMSetType()`, `SNESNASMGetType()`, `SNESNASMSetSubdomains()`, `SNESNASMGetSubdomains()`,
830:           `SNESNASMGetSubdomainVecs()`, `SNESNASMSetComputeFinalJacobian()`, `SNESNASMSetDamping()`, `SNESNASMGetDamping()`, `SNESNASMSetWeight()`,
831:           `SNESNASMGetSNES()`, `SNESNASMGetNumber()`
832: M*/

834: PETSC_EXTERN PetscErrorCode SNESCreate_NASM(SNES snes)
835: {
836:   SNES_NASM *nasm;

838:   PetscFunctionBegin;
839:   PetscCall(PetscNew(&nasm));
840:   snes->data = (void *)nasm;

842:   nasm->n             = PETSC_DECIDE;
843:   nasm->subsnes       = NULL;
844:   nasm->x             = NULL;
845:   nasm->xl            = NULL;
846:   nasm->y             = NULL;
847:   nasm->b             = NULL;
848:   nasm->oscatter      = NULL;
849:   nasm->oscatter_copy = NULL;
850:   nasm->iscatter      = NULL;
851:   nasm->gscatter      = NULL;
852:   nasm->damping       = 1.;

854:   nasm->type          = PC_ASM_BASIC;
855:   nasm->finaljacobian = PETSC_FALSE;
856:   nasm->weight_set    = PETSC_FALSE;

858:   snes->ops->destroy        = SNESDestroy_NASM;
859:   snes->ops->setup          = SNESSetUp_NASM;
860:   snes->ops->setfromoptions = SNESSetFromOptions_NASM;
861:   snes->ops->view           = SNESView_NASM;
862:   snes->ops->solve          = SNESSolve_NASM;
863:   snes->ops->reset          = SNESReset_NASM;

865:   snes->usesksp = PETSC_FALSE;
866:   snes->usesnpc = PETSC_FALSE;

868:   snes->alwayscomputesfinalresidual = PETSC_FALSE;

870:   nasm->fjtype              = 0;
871:   nasm->xinit               = NULL;
872:   nasm->eventrestrictinterp = 0;
873:   nasm->eventsubsolve       = 0;

875:   if (!snes->tolerancesset) {
876:     snes->max_its   = 10000;
877:     snes->max_funcs = 10000;
878:   }

880:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetType_C", SNESNASMSetType_NASM));
881:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetType_C", SNESNASMGetType_NASM));
882:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetSubdomains_C", SNESNASMSetSubdomains_NASM));
883:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetSubdomains_C", SNESNASMGetSubdomains_NASM));
884:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetDamping_C", SNESNASMSetDamping_NASM));
885:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetDamping_C", SNESNASMGetDamping_NASM));
886:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMGetSubdomainVecs_C", SNESNASMGetSubdomainVecs_NASM));
887:   PetscCall(PetscObjectComposeFunction((PetscObject)snes, "SNESNASMSetComputeFinalJacobian_C", SNESNASMSetComputeFinalJacobian_NASM));
888:   PetscFunctionReturn(PETSC_SUCCESS);
889: }

891: /*@
892:   SNESNASMGetSNES - Gets a subsolver

894:   Not Collective

896:   Input Parameters:
897: + snes - the `SNES` context
898: - i    - the number of the subsnes to get

900:   Output Parameter:
901: . subsnes - the subsolver context

903:   Level: intermediate

905: .seealso: [](ch_snes), `SNESNASM`, `SNESNASMGetNumber()`
906: @*/
907: PetscErrorCode SNESNASMGetSNES(SNES snes, PetscInt i, SNES *subsnes)
908: {
909:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

911:   PetscFunctionBegin;
912:   PetscCheck(i >= 0 && i < nasm->n, PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_OUTOFRANGE, "No such subsolver");
913:   *subsnes = nasm->subsnes[i];
914:   PetscFunctionReturn(PETSC_SUCCESS);
915: }

917: /*@
918:   SNESNASMGetNumber - Gets number of subsolvers

920:   Not Collective

922:   Input Parameter:
923: . snes - the `SNES` context

925:   Output Parameter:
926: . n - the number of subsolvers

928:   Level: intermediate

930: .seealso: [](ch_snes), `SNESNASM`, `SNESNASMGetSNES()`
931: @*/
932: PetscErrorCode SNESNASMGetNumber(SNES snes, PetscInt *n)
933: {
934:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

936:   PetscFunctionBegin;
937:   *n = nasm->n;
938:   PetscFunctionReturn(PETSC_SUCCESS);
939: }

941: /*@
942:   SNESNASMSetWeight - Sets weight to use when adding overlapping updates

944:   Collective

946:   Input Parameters:
947: + snes   - the `SNES` context
948: - weight - the weights to use (typically 1/N for each dof, where N is the number of patches it appears in)

950:   Level: intermediate

952: .seealso: [](ch_snes), `SNESNASM`
953: @*/
954: PetscErrorCode SNESNASMSetWeight(SNES snes, Vec weight)
955: {
956:   SNES_NASM *nasm = (SNES_NASM *)snes->data;

958:   PetscFunctionBegin;

960:   PetscCall(VecDestroy(&nasm->weight));
961:   nasm->weight_set = PETSC_TRUE;
962:   nasm->weight     = weight;
963:   PetscCall(PetscObjectReference((PetscObject)nasm->weight));

965:   PetscFunctionReturn(PETSC_SUCCESS);
966: }