xref: /petsc/src/snes/impls/patch/snespatch.c (revision ab270f919ae7ceeb10c0d05c3e81e905961f49c1)
1 /*
2       Defines a SNES that can consist of a collection of SNESes on patches of the domain
3 */
4 #include <petsc/private/snesimpl.h> /*I "petscsnes.h" I*/
5 #include <petsc/private/pcpatchimpl.h> /* We need internal access to PCPatch right now, until that part is moved to Plex */
6 
7 typedef struct {
8   PC pc; /* The linear patch preconditioner */
9   SNESCompositeType type;
10 } SNES_Patch;
11 
12 static PetscErrorCode SNESPatchComputeResidual_Private(SNES snes, Vec x, Vec F, void *ctx)
13 {
14   PC             pc      = (PC) ctx;
15   PC_PATCH      *pcpatch = (PC_PATCH *) pc->data;
16   PetscErrorCode ierr;
17 
18   PetscFunctionBegin;
19   ierr = PCPatchComputeFunction_Internal(pc, x, F, pcpatch->currentPatch);CHKERRQ(ierr);
20   PetscFunctionReturn(0);
21 }
22 
23 static PetscErrorCode SNESPatchComputeJacobian_Private(SNES snes, Vec x, Mat J, Mat M, void *ctx)
24 {
25   PC             pc      = (PC) ctx;
26   PC_PATCH      *pcpatch = (PC_PATCH *) pc->data;
27   PetscErrorCode ierr;
28 
29   PetscFunctionBegin;
30   ierr = PCPatchComputeOperator_Internal(pc, x, M, pcpatch->currentPatch, PETSC_FALSE);CHKERRQ(ierr);
31   PetscFunctionReturn(0);
32 }
33 
34 static PetscErrorCode PCSetUp_PATCH_Nonlinear(PC pc)
35 {
36   PC_PATCH      *patch = (PC_PATCH *) pc->data;
37   const char    *prefix;
38   PetscInt       i;
39   PetscErrorCode ierr;
40 
41   PetscFunctionBegin;
42   if (!pc->setupcalled) {
43     ierr = PetscMalloc1(patch->npatch, &patch->solver);CHKERRQ(ierr);
44     ierr = PCGetOptionsPrefix(pc, &prefix);CHKERRQ(ierr);
45     for (i = 0; i < patch->npatch; ++i) {
46       SNES snes;
47       KSP  subksp;
48 
49       ierr = SNESCreate(PETSC_COMM_SELF, &snes);CHKERRQ(ierr);
50       ierr = SNESSetOptionsPrefix(snes, prefix);CHKERRQ(ierr);
51       ierr = SNESAppendOptionsPrefix(snes, "sub_");CHKERRQ(ierr);
52       ierr = PetscObjectIncrementTabLevel((PetscObject) snes, (PetscObject) pc, 1);CHKERRQ(ierr);
53       ierr = SNESGetKSP(snes, &subksp);CHKERRQ(ierr);
54       ierr = PetscObjectIncrementTabLevel((PetscObject) subksp, (PetscObject) pc, 1);CHKERRQ(ierr);
55       ierr = PetscLogObjectParent((PetscObject) pc, (PetscObject) snes);CHKERRQ(ierr);
56       patch->solver[i] = (PetscObject) snes;
57     }
58   }
59   for (i = 0; i < patch->npatch; ++i) {
60     SNES snes = (SNES) patch->solver[i];
61 
62     ierr = SNESSetFunction(snes, patch->patchX[i], SNESPatchComputeResidual_Private, pc);CHKERRQ(ierr);
63     ierr = SNESSetJacobian(snes, patch->mat[i], patch->mat[i], SNESPatchComputeJacobian_Private, pc);CHKERRQ(ierr);
64   }
65   if (!pc->setupcalled && patch->optionsSet) for (i = 0; i < patch->npatch; ++i) {ierr = SNESSetFromOptions((SNES) patch->solver[i]);CHKERRQ(ierr);}
66   PetscFunctionReturn(0);
67 }
68 
69 static PetscErrorCode PCApply_PATCH_Nonlinear(PC pc, PetscInt i, Vec x, Vec y)
70 {
71   PC_PATCH      *patch = (PC_PATCH *) pc->data;
72   PetscErrorCode ierr;
73 
74   PetscFunctionBegin;
75   patch->currentPatch = i;
76   ierr = PetscLogEventBegin(PC_Patch_Solve, pc, 0, 0, 0);CHKERRQ(ierr);
77   ierr = VecCopy(x, y);CHKERRQ(ierr);
78   ierr = SNESSolve((SNES) patch->solver[i], NULL, y);CHKERRQ(ierr); /* FIXME: RHS for FAS */
79   ierr = VecAXPY(y, -1.0, x);CHKERRQ(ierr);
80   ierr = PetscLogEventEnd(PC_Patch_Solve, pc, 0, 0, 0);CHKERRQ(ierr);
81   PetscFunctionReturn(0);
82 }
83 
84 static PetscErrorCode PCReset_PATCH_Nonlinear(PC pc)
85 {
86   PC_PATCH      *patch = (PC_PATCH *) pc->data;
87   PetscInt       i;
88   PetscErrorCode ierr;
89 
90   PetscFunctionBegin;
91   if (patch->solver) {
92     for (i = 0; i < patch->npatch; ++i) {ierr = SNESReset((SNES) patch->solver[i]);CHKERRQ(ierr);}
93   }
94   PetscFunctionReturn(0);
95 }
96 
97 static PetscErrorCode PCDestroy_PATCH_Nonlinear(PC pc)
98 {
99   PC_PATCH      *patch = (PC_PATCH *) pc->data;
100   PetscInt       i;
101   PetscErrorCode ierr;
102 
103   PetscFunctionBegin;
104   if (patch->solver) {
105     for (i = 0; i < patch->npatch; ++i) {ierr = SNESDestroy((SNES *) &patch->solver[i]);CHKERRQ(ierr);}
106     ierr = PetscFree(patch->solver);CHKERRQ(ierr);
107   }
108   PetscFunctionReturn(0);
109 }
110 
111 static PetscErrorCode SNESSetUp_Patch(SNES snes)
112 {
113   SNES_Patch    *patch = (SNES_Patch *) snes->data;
114   DM             dm;
115   Mat            dummy;
116   Vec            F;
117   PetscInt       n, N;
118   PetscErrorCode ierr;
119 
120   PetscFunctionBegin;
121   ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr);
122   ierr = PCSetDM(patch->pc, dm);CHKERRQ(ierr);
123   ierr = SNESGetFunction(snes, &F, NULL, NULL);CHKERRQ(ierr);
124   ierr = VecGetLocalSize(F, &n);CHKERRQ(ierr);
125   ierr = VecGetSize(F, &N);CHKERRQ(ierr);
126   ierr = MatCreateShell(PetscObjectComm((PetscObject) snes), n, n, N, N, (void *) snes, &dummy);CHKERRQ(ierr);
127   ierr = PCSetOperators(patch->pc, dummy, dummy);CHKERRQ(ierr);
128   ierr = MatDestroy(&dummy);CHKERRQ(ierr);
129   ierr = PCSetUp(patch->pc);CHKERRQ(ierr);
130   /* allocate workspace */
131   PetscFunctionReturn(0);
132 }
133 
134 static PetscErrorCode SNESReset_Patch(SNES snes)
135 {
136   SNES_Patch    *patch = (SNES_Patch *) snes->data;
137   PetscErrorCode ierr;
138 
139   PetscFunctionBegin;
140   ierr = PCReset(patch->pc);CHKERRQ(ierr);
141   PetscFunctionReturn(0);
142 }
143 
144 static PetscErrorCode SNESDestroy_Patch(SNES snes)
145 {
146   SNES_Patch    *patch = (SNES_Patch *) snes->data;
147   PetscErrorCode ierr;
148 
149   PetscFunctionBegin;
150   ierr = SNESReset_Patch(snes);CHKERRQ(ierr);
151   ierr = PCDestroy(&patch->pc);CHKERRQ(ierr);
152   ierr = PetscFree(snes->data);CHKERRQ(ierr);
153   PetscFunctionReturn(0);
154 }
155 
156 static PetscErrorCode SNESSetFromOptions_Patch(PetscOptionItems *PetscOptionsObject, SNES snes)
157 {
158   SNES_Patch    *patch = (SNES_Patch *) snes->data;
159   PetscBool      flg;
160   const char    *prefix;
161   PetscErrorCode ierr;
162 
163   PetscFunctionBegin;
164   ierr = PetscOptionsHead(PetscOptionsObject, "Patch nonlinear preconditioner options");CHKERRQ(ierr);
165   ierr = PetscOptionsEnum("-snes_patch_type", "Type of composition", "SNESPatchSetType", SNESCompositeTypes, (PetscEnum) patch->type, (PetscEnum *) &patch->type, &flg);CHKERRQ(ierr);
166   if (flg) {ierr = SNESPatchSetType(snes, patch->type);CHKERRQ(ierr);}
167   ierr = PetscOptionsTail();CHKERRQ(ierr);
168 
169   ierr = PetscObjectGetOptionsPrefix((PetscObject)snes, &prefix);CHKERRQ(ierr);
170   ierr = PetscObjectSetOptionsPrefix((PetscObject)patch->pc, prefix);CHKERRQ(ierr);
171   ierr = PCSetFromOptions(patch->pc);CHKERRQ(ierr);
172   PetscFunctionReturn(0);
173 }
174 
175 static PetscErrorCode SNESView_Patch(SNES snes,PetscViewer viewer)
176 {
177   SNES_Patch    *patch = (SNES_Patch *) snes->data;
178   PetscBool      iascii;
179   PetscErrorCode ierr;
180 
181   PetscFunctionBegin;
182   ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr);
183   if (iascii) {
184     ierr = PetscViewerASCIIPrintf(viewer,"  type - %s\n",SNESCompositeTypes[patch->type]);CHKERRQ(ierr);
185   }
186   ierr = PetscViewerASCIIPushTab(viewer);CHKERRQ(ierr);
187   ierr = PCView(patch->pc, viewer);CHKERRQ(ierr);
188   ierr = PetscViewerASCIIPopTab(viewer);CHKERRQ(ierr);
189   PetscFunctionReturn(0);
190 }
191 
192 static PetscErrorCode SNESSolve_Patch(SNES snes)
193 {
194   SNES_Patch *patch = (SNES_Patch *) snes->data;
195   Vec x, y;
196   PetscErrorCode ierr;
197 
198   PetscFunctionBegin;
199   ierr = SNESGetSolution(snes, &x);CHKERRQ(ierr);
200   ierr = SNESGetSolutionUpdate(snes, &y);CHKERRQ(ierr);
201   ierr = PCApply(patch->pc, x, y);
202   ierr = VecAXPY(x, 1.0, y);CHKERRQ(ierr); /* FIXME: line search */
203   ierr = SNESSetConvergedReason(snes, SNES_CONVERGED_ITS);CHKERRQ(ierr);
204   PetscFunctionReturn(0);
205 }
206 
207 static PetscErrorCode SNESPatchSetType_Patch(SNES snes, SNESCompositeType type)
208 {
209   SNES_Patch *patch = (SNES_Patch *) snes->data;
210 
211   PetscFunctionBegin;
212   patch->type = type;
213   PetscFunctionReturn(0);
214 }
215 
216 /*MC
217   SNESPATCH - Solve a nonlinear problem by composing together many nonlinear solvers on patches
218 
219   Options Database Keys:
220 . -snes_patch_type <type: one of additive, multiplicative, symmetric_multiplicative> - Sets composition type
221 
222   Level: intermediate
223 
224   Concepts: composing solvers
225 
226 .seealso:  SNESCreate(), SNESSetType(), SNESType (for list of available types), SNES,
227            PCPATCH
228 
229    References:
230 .  1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", SIAM Review, 57(4), 2015
231 
232 M*/
233 PETSC_EXTERN PetscErrorCode SNESCreate_Patch(SNES snes)
234 {
235   PetscErrorCode ierr;
236   SNES_Patch    *patch;
237 
238   PetscFunctionBegin;
239   ierr = PetscNewLog(snes, &patch);CHKERRQ(ierr);
240 
241   snes->ops->solve          = SNESSolve_Patch;
242   snes->ops->setup          = SNESSetUp_Patch;
243   snes->ops->reset          = SNESReset_Patch;
244   snes->ops->destroy        = SNESDestroy_Patch;
245   snes->ops->setfromoptions = SNESSetFromOptions_Patch;
246   snes->ops->view           = SNESView_Patch;
247 
248   snes->alwayscomputesfinalresidual = PETSC_FALSE;
249 
250   snes->data = (void *) patch;
251   patch->type = SNES_COMPOSITE_ADDITIVE;
252   ierr = PCCreate(PetscObjectComm((PetscObject) snes), &patch->pc);CHKERRQ(ierr);
253   ierr = PCSetType(patch->pc, PCPATCH);CHKERRQ(ierr);
254 
255   ((PC_PATCH *) patch->pc->data)->setupsolver   = PCSetUp_PATCH_Nonlinear;
256   ((PC_PATCH *) patch->pc->data)->applysolver   = PCApply_PATCH_Nonlinear;
257   ((PC_PATCH *) patch->pc->data)->resetsolver   = PCReset_PATCH_Nonlinear;
258   ((PC_PATCH *) patch->pc->data)->destroysolver = PCDestroy_PATCH_Nonlinear;
259 
260   ierr = PetscObjectComposeFunction((PetscObject) snes, "SNESPatchSetType_C", SNESPatchSetType_Patch);CHKERRQ(ierr);
261   PetscFunctionReturn(0);
262 }
263 
264 /*@
265   SNESPatchSetType - Sets the type of composition.
266 
267   Logically Collective on SNES
268 
269   Input Parameter:
270 + snes - the preconditioner context
271 - type - SNES_COMPOSITE_ADDITIVE (default), SNES_COMPOSITE_MULTIPLICATIVE, etc.
272 
273   Options Database Key:
274 . -snes_composite_type <type: one of additive, multiplicative, etc> - Sets composite preconditioner type
275 
276   Level: Developer
277 
278 .keywords: SNES, set, type, composite preconditioner, additive, multiplicative
279 @*/
280 PetscErrorCode SNESPatchSetType(SNES snes, SNESCompositeType type)
281 {
282   PetscErrorCode ierr;
283 
284   PetscFunctionBegin;
285   PetscValidHeaderSpecific(snes, SNES_CLASSID, 1);
286   PetscValidLogicalCollectiveEnum(snes,type, 2);
287   ierr = PetscTryMethod(snes, "SNESPatchSetType_C", (SNES,SNESCompositeType), (snes,type));CHKERRQ(ierr);
288   PetscFunctionReturn(0);
289 }
290 
291 PetscErrorCode SNESPatchSetDiscretisationInfo(SNES snes, PetscInt nsubspaces, DM *dms, PetscInt *bs, PetscInt *nodesPerCell, const PetscInt **cellNodeMap,
292                                             const PetscInt *subspaceOffsets, PetscInt numGhostBcs, const PetscInt *ghostBcNodes, PetscInt numGlobalBcs, const PetscInt *globalBcNodes)
293 {
294   SNES_Patch    *patch = (SNES_Patch *) snes->data;
295   PetscErrorCode ierr;
296   DM dm;
297 
298   PetscFunctionBegin;
299   ierr = SNESGetDM(snes, &dm);CHKERRQ(ierr);
300   if (!dm) SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_WRONGSTATE, "DM not yet set on patch SNES\n");
301   ierr = PCSetDM(patch->pc, dm);CHKERRQ(ierr);
302   ierr = PCPatchSetDiscretisationInfo(patch->pc, nsubspaces, dms, bs, nodesPerCell, cellNodeMap, subspaceOffsets, numGhostBcs, ghostBcNodes, numGlobalBcs, globalBcNodes);CHKERRQ(ierr);
303   PetscFunctionReturn(0);
304 }
305 
306 PetscErrorCode SNESPatchSetComputeOperator(SNES snes, PetscErrorCode (*func)(PC, PetscInt, Vec, Mat, IS, PetscInt, const PetscInt *, void *), void *ctx)
307 {
308   SNES_Patch    *patch = (SNES_Patch *) snes->data;
309   PetscErrorCode ierr;
310 
311   PetscFunctionBegin;
312   ierr = PCPatchSetComputeOperator(patch->pc, func, ctx);CHKERRQ(ierr);
313   PetscFunctionReturn(0);
314 }
315 
316 PetscErrorCode SNESPatchSetComputeFunction(SNES snes, PetscErrorCode (*func)(PC, PetscInt, Vec, Vec, IS, PetscInt, const PetscInt *, void *), void *ctx)
317 {
318   SNES_Patch    *patch = (SNES_Patch *) snes->data;
319   PetscErrorCode ierr;
320 
321   PetscFunctionBegin;
322   ierr = PCPatchSetComputeFunction(patch->pc, func, ctx);CHKERRQ(ierr);
323   PetscFunctionReturn(0);
324 }
325 
326 PetscErrorCode SNESPatchSetConstructType(SNES snes, PCPatchConstructType ctype, PetscErrorCode (*func)(PC, PetscInt *, IS **, IS *, void *), void *ctx)
327 {
328   SNES_Patch    *patch = (SNES_Patch *) snes->data;
329   PetscErrorCode ierr;
330 
331   PetscFunctionBegin;
332   ierr = PCPatchSetConstructType(patch->pc, ctype, func, ctx);CHKERRQ(ierr);
333   PetscFunctionReturn(0);
334 }
335 
336 PetscErrorCode SNESPatchSetCellNumbering(SNES snes, PetscSection cellNumbering)
337 {
338   SNES_Patch    *patch = (SNES_Patch *) snes->data;
339   PetscErrorCode ierr;
340 
341   PetscFunctionBegin;
342   ierr = PCPatchSetCellNumbering(patch->pc, cellNumbering);CHKERRQ(ierr);
343   PetscFunctionReturn(0);
344 }
345