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