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