Actual source code: dmsnes.c
petsc-3.3-p7 2013-05-11
1: #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/
2: #include <petscdm.h> /*I "petscdm.h" I*/
6: /*@C
7: SNESDMComputeFunction - This is a universal function evaluation routine that
8: may be used with SNESSetFunction() as long as the user context has a DM
9: as its first record and the user has called DMSetLocalFunction().
11: Collective on SNES
13: Input Parameters:
14: + snes - the SNES context
15: . X - input vector
16: . F - function vector
17: - ptr - pointer to a structure that must have a DM as its first entry.
18: This ptr must have been passed into SNESDMComputeFunction() as the context.
20: Level: intermediate
22: .seealso: DMSetLocalFunction(), DMSetLocalJacobian(), SNESSetFunction(), SNESSetJacobian()
23: @*/
24: PetscErrorCode SNESDMComputeFunction(SNES snes, Vec X, Vec F, void *ptr)
25: {
26: DM dm = *(DM*) ptr;
27: PetscErrorCode (*lf)(DM, Vec, Vec, void *);
28: Vec localX, localF;
29: PetscInt N, n;
30: PetscErrorCode ierr;
36: if (!dm) SETERRQ(((PetscObject)snes)->comm, PETSC_ERR_ARG_WRONGSTATE, "Looks like you called SNESSetFromFuntion(snes,SNESDMComputeFunction,) without the DM context");
39: /* determine whether X = localX */
40: DMGetLocalVector(dm, &localX);
41: DMGetLocalVector(dm, &localF);
42: VecGetSize(X, &N);
43: VecGetSize(localX, &n);
45: if (n != N){ /* X != localX */
46: /* Scatter ghost points to local vector, using the 2-step process
47: DMGlobalToLocalBegin(), DMGlobalToLocalEnd().
48: */
49: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, localX);
50: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, localX);
51: } else {
52: DMRestoreLocalVector(dm, &localX);
53: localX = X;
54: }
55: DMGetLocalFunction(dm, &lf);
56: (*lf)(dm, localX, localF, ptr);
57: if (n != N){
58: DMRestoreLocalVector(dm, &localX);
59: }
60: DMLocalToGlobalBegin(dm, localF, ADD_VALUES, F);
61: DMLocalToGlobalEnd(dm, localF, ADD_VALUES, F);
62: DMRestoreLocalVector(dm, &localF);
63: return(0);
64: }
68: /*
69: SNESDMComputeJacobian - This is a universal Jacobian evaluation routine that
70: may be used with SNESSetJacobian() as long as the user context has a DM
71: as its first record and the user has called DMSetLocalJacobian().
73: Collective on SNES
75: Input Parameters:
76: + snes - the SNES context
77: . X - input vector
78: . J - Jacobian
79: . B - Jacobian used in preconditioner (usally same as J)
80: . flag - indicates if the matrix changed its structure
81: - ptr - pointer to a structure that must have a DM as its first entry.
82: This ptr must have been passed into SNESDMComputeFunction() as the context.
84: Level: intermediate
86: .seealso: DMSetLocalFunction(), DMSetLocalJacobian(), SNESSetFunction(), SNESSetJacobian()
87: */
88: PetscErrorCode SNESDMComputeJacobian(SNES snes, Vec X, Mat *J, Mat *B, MatStructure *flag, void *ptr)
89: {
90: DM dm = *(DM*) ptr;
91: PetscErrorCode (*lj)(DM, Vec, Mat, Mat, void *);
92: Vec localX;
93: PetscErrorCode ierr;
96: DMGetLocalVector(dm, &localX);
97: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, localX);
98: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, localX);
99: DMGetLocalJacobian(dm, &lj);
100: (*lj)(dm, localX, *J, *B, ptr);
101: DMRestoreLocalVector(dm, &localX);
102: /* Assemble true Jacobian; if it is different */
103: if (*J != *B) {
104: MatAssemblyBegin(*J, MAT_FINAL_ASSEMBLY);
105: MatAssemblyEnd(*J, MAT_FINAL_ASSEMBLY);
106: }
107: MatSetOption(*B, MAT_NEW_NONZERO_LOCATION_ERR, PETSC_TRUE);
108: *flag = SAME_NONZERO_PATTERN;
109: return(0);
110: }
114: static PetscErrorCode PetscContainerDestroy_SNESDM(void *ctx)
115: {
117: SNESDM sdm = (SNESDM)ctx;
120: if (sdm->destroy) {(*sdm->destroy)(sdm);}
121: PetscFree(sdm);
122: return(0);
123: }
127: /* Attaches the SNESDM to the coarse level.
128: * Under what conditions should we copy versus duplicate?
129: */
130: static PetscErrorCode DMCoarsenHook_SNESDM(DM dm,DM dmc,void *ctx)
131: {
135: DMSNESCopyContext(dm,dmc);
136: return(0);
137: }
141: /* This could restrict auxiliary information to the coarse level.
142: */
143: static PetscErrorCode DMRestrictHook_SNESDM(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
144: {
147: return(0);
148: }
152: /*@C
153: DMSNESGetContext - get read-only private SNESDM context from a DM
155: Not Collective
157: Input Argument:
158: . dm - DM to be used with SNES
160: Output Argument:
161: . snesdm - private SNESDM context
163: Level: developer
165: Notes:
166: Use DMSNESGetContextWrite() if write access is needed. The DMSNESSetXXX API should be used wherever possible.
168: .seealso: DMSNESGetContextWrite()
169: @*/
170: PetscErrorCode DMSNESGetContext(DM dm,SNESDM *snesdm)
171: {
173: PetscContainer container;
174: SNESDM sdm;
178: PetscObjectQuery((PetscObject)dm,"SNESDM",(PetscObject*)&container);
179: if (container) {
180: PetscContainerGetPointer(container,(void**)snesdm);
181: } else {
182: PetscInfo(dm,"Creating new SNESDM\n");
183: PetscContainerCreate(((PetscObject)dm)->comm,&container);
184: PetscNewLog(dm,struct _n_SNESDM,&sdm);
185: PetscContainerSetPointer(container,sdm);
186: PetscContainerSetUserDestroy(container,PetscContainerDestroy_SNESDM);
187: PetscObjectCompose((PetscObject)dm,"SNESDM",(PetscObject)container);
188: DMCoarsenHookAdd(dm,DMCoarsenHook_SNESDM,DMRestrictHook_SNESDM,PETSC_NULL);
189: PetscContainerGetPointer(container,(void**)snesdm);
190: PetscContainerDestroy(&container);
191: }
192: return(0);
193: }
197: /*@C
198: DMSNESGetContextWrite - get write access to private SNESDM context from a DM
200: Not Collective
202: Input Argument:
203: . dm - DM to be used with SNES
205: Output Argument:
206: . snesdm - private SNESDM context
208: Level: developer
210: .seealso: DMSNESGetContext()
211: @*/
212: PetscErrorCode DMSNESGetContextWrite(DM dm,SNESDM *snesdm)
213: {
215: SNESDM sdm;
219: DMSNESGetContext(dm,&sdm);
220: if (!sdm->originaldm) sdm->originaldm = dm;
221: if (sdm->originaldm != dm) { /* Copy on write */
222: PetscContainer container;
223: SNESDM oldsdm = sdm;
224: PetscInfo(dm,"Copying SNESDM due to write\n");
225: PetscContainerCreate(((PetscObject)dm)->comm,&container);
226: PetscNewLog(dm,struct _n_SNESDM,&sdm);
227: PetscMemcpy(sdm,oldsdm,sizeof *sdm);
228: PetscContainerSetPointer(container,sdm);
229: PetscContainerSetUserDestroy(container,PetscContainerDestroy_SNESDM);
230: PetscObjectCompose((PetscObject)dm,"SNESDM",(PetscObject)container);
231: PetscContainerDestroy(&container);
232: }
233: *snesdm = sdm;
234: return(0);
235: }
239: /*@C
240: DMSNESCopyContext - copies a DM context to a new DM
242: Logically Collective
244: Input Arguments:
245: + dmsrc - DM to obtain context from
246: - dmdest - DM to add context to
248: Level: developer
250: Note:
251: The context is copied by reference. This function does not ensure that a context exists.
253: .seealso: DMSNESGetContext(), SNESSetDM()
254: @*/
255: PetscErrorCode DMSNESCopyContext(DM dmsrc,DM dmdest)
256: {
258: PetscContainer container;
263: PetscObjectQuery((PetscObject)dmsrc,"SNESDM",(PetscObject*)&container);
264: if (container) {
265: PetscObjectCompose((PetscObject)dmdest,"SNESDM",(PetscObject)container);
266: DMCoarsenHookAdd(dmdest,DMCoarsenHook_SNESDM,DMRestrictHook_SNESDM,PETSC_NULL);
267: }
268: return(0);
269: }
273: /*@C
274: DMSNESSetFunction - set SNES residual evaluation function
276: Not Collective
278: Input Arguments:
279: + dm - DM to be used with SNES
280: . func - residual evaluation function, see SNESSetFunction() for calling sequence
281: - ctx - context for residual evaluation
283: Level: advanced
285: Note:
286: SNESSetFunction() is normally used, but it calls this function internally because the user context is actually
287: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
288: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
290: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
291: @*/
292: PetscErrorCode DMSNESSetFunction(DM dm,PetscErrorCode (*func)(SNES,Vec,Vec,void*),void *ctx)
293: {
295: SNESDM sdm;
299: DMSNESGetContextWrite(dm,&sdm);
300: if (func) sdm->computefunction = func;
301: if (ctx) sdm->functionctx = ctx;
302: return(0);
303: }
307: /*@C
308: DMSNESGetFunction - get SNES residual evaluation function
310: Not Collective
312: Input Argument:
313: . dm - DM to be used with SNES
315: Output Arguments:
316: + func - residual evaluation function, see SNESSetFunction() for calling sequence
317: - ctx - context for residual evaluation
319: Level: advanced
321: Note:
322: SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
323: associated with the DM.
325: .seealso: DMSNESSetContext(), DMSNESSetFunction(), SNESSetFunction()
326: @*/
327: PetscErrorCode DMSNESGetFunction(DM dm,PetscErrorCode (**func)(SNES,Vec,Vec,void*),void **ctx)
328: {
330: SNESDM sdm;
334: DMSNESGetContext(dm,&sdm);
335: if (func) *func = sdm->computefunction;
336: if (ctx) *ctx = sdm->functionctx;
337: return(0);
338: }
342: /*@C
343: DMSNESSetGS - set SNES Gauss-Seidel relaxation function
345: Not Collective
347: Input Argument:
348: + dm - DM to be used with SNES
349: . func - relaxation function, see SNESSetGS() for calling sequence
350: - ctx - context for residual evaluation
352: Level: advanced
354: Note:
355: SNESSetGS() is normally used, but it calls this function internally because the user context is actually
356: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
357: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
359: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), DMSNESSetFunction()
360: @*/
361: PetscErrorCode DMSNESSetGS(DM dm,PetscErrorCode (*func)(SNES,Vec,Vec,void*),void *ctx)
362: {
364: SNESDM sdm;
368: DMSNESGetContextWrite(dm,&sdm);
369: if (func) sdm->computegs = func;
370: if (ctx) sdm->gsctx = ctx;
371: return(0);
372: }
376: /*@C
377: DMSNESGetGS - get SNES Gauss-Seidel relaxation function
379: Not Collective
381: Input Argument:
382: . dm - DM to be used with SNES
384: Output Arguments:
385: + func - relaxation function, see SNESSetGS() for calling sequence
386: - ctx - context for residual evaluation
388: Level: advanced
390: Note:
391: SNESGetGS() is normally used, but it calls this function internally because the user context is actually
392: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
393: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
395: .seealso: DMSNESSetContext(), SNESGetGS(), DMSNESGetJacobian(), DMSNESGetFunction()
396: @*/
397: PetscErrorCode DMSNESGetGS(DM dm,PetscErrorCode (**func)(SNES,Vec,Vec,void*),void **ctx)
398: {
400: SNESDM sdm;
404: DMSNESGetContext(dm,&sdm);
405: if (func) *func = sdm->computegs;
406: if (ctx) *ctx = sdm->gsctx;
407: return(0);
408: }
412: /*@C
413: DMSNESSetJacobian - set SNES Jacobian evaluation function
415: Not Collective
417: Input Argument:
418: + dm - DM to be used with SNES
419: . func - Jacobian evaluation function, see SNESSetJacobian() for calling sequence
420: - ctx - context for residual evaluation
422: Level: advanced
424: Note:
425: SNESSetJacobian() is normally used, but it calls this function internally because the user context is actually
426: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
427: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
429: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESGetJacobian(), SNESSetJacobian()
430: @*/
431: PetscErrorCode DMSNESSetJacobian(DM dm,PetscErrorCode (*func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
432: {
434: SNESDM sdm;
438: DMSNESGetContextWrite(dm,&sdm);
439: if (func) sdm->computejacobian = func;
440: if (ctx) sdm->jacobianctx = ctx;
441: return(0);
442: }
446: /*@C
447: DMSNESGetJacobian - get SNES Jacobian evaluation function
449: Not Collective
451: Input Argument:
452: . dm - DM to be used with SNES
454: Output Arguments:
455: + func - Jacobian evaluation function, see SNESSetJacobian() for calling sequence
456: - ctx - context for residual evaluation
458: Level: advanced
460: Note:
461: SNESGetJacobian() is normally used, but it calls this function internally because the user context is actually
462: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
463: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
465: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
466: @*/
467: PetscErrorCode DMSNESGetJacobian(DM dm,PetscErrorCode (**func)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void **ctx)
468: {
470: SNESDM sdm;
474: DMSNESGetContext(dm,&sdm);
475: if (func) *func = sdm->computejacobian;
476: if (ctx) *ctx = sdm->jacobianctx;
477: return(0);
478: }
482: /*@C
483: DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions.
485: Not Collective
487: Input Argument:
488: + dm - DM to be used with SNES
489: . func - RHS evaluation function, see SNESSetFunction() for calling sequence
490: . pjac - Picard matrix evaluation function, see SNESSetJacobian() for calling sequence
491: - ctx - context for residual evaluation
493: Level: advanced
495: .seealso: SNESSetPicard(), DMSNESSetFunction(), DMSNESSetJacobian()
496: @*/
497: PetscErrorCode DMSNESSetPicard(DM dm,PetscErrorCode (*pfunc)(SNES,Vec,Vec,void*),PetscErrorCode (*pjac)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx)
498: {
500: SNESDM sdm;
504: DMSNESGetContext(dm,&sdm);
505: if (pfunc) sdm->computepfunction = pfunc;
506: if (pjac) sdm->computepjacobian = pjac;
507: if (ctx) sdm->pctx = ctx;
508: return(0);
509: }
514: /*@C
515: DMSNESGetPicard - get SNES Picard iteration evaluation functions
517: Not Collective
519: Input Argument:
520: . dm - DM to be used with SNES
522: Output Arguments:
523: + pfunc - Jacobian evaluation function, see SNESSetJacobian() for calling sequence
524: . pjac - RHS evaluation function, see SNESSetFunction() for calling sequence
525: - ctx - context for residual evaluation
527: Level: advanced
529: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
530: @*/
531: PetscErrorCode DMSNESGetPicard(DM dm,PetscErrorCode (**pfunc)(SNES,Vec,Vec,void*),PetscErrorCode (**pjac)(SNES,Vec,Mat*,Mat*,MatStructure*,void*),void **ctx)
532: {
534: SNESDM sdm;
538: DMSNESGetContext(dm,&sdm);
539: if (pfunc) *pfunc = sdm->computepfunction;
540: if (pjac) *pjac = sdm->computepjacobian;
541: if (ctx) *ctx = sdm->pctx;
542: return(0);
543: }
548: static PetscErrorCode SNESDefaultComputeFunction_DMLegacy(SNES snes,Vec X,Vec F,void *ctx)
549: {
551: DM dm;
554: SNESGetDM(snes,&dm);
555: DMComputeFunction(dm,X,F);
556: return(0);
557: }
561: static PetscErrorCode SNESDefaultComputeJacobian_DMLegacy(SNES snes,Vec X,Mat *A,Mat *B,MatStructure *mstr,void *ctx)
562: {
564: DM dm;
567: SNESGetDM(snes,&dm);
568: DMComputeJacobian(dm,X,*A,*B,mstr);
569: return(0);
570: }
574: /* Sets up calling of legacy DM routines */
575: PetscErrorCode DMSNESSetUpLegacy(DM dm)
576: {
578: SNESDM sdm;
581: DMSNESGetContext(dm,&sdm);
582: if (!sdm->computefunction) {DMSNESSetFunction(dm,SNESDefaultComputeFunction_DMLegacy,PETSC_NULL);}
583: DMSNESGetContext(dm,&sdm);
584: if (!sdm->computejacobian) {DMSNESSetJacobian(dm,SNESDefaultComputeJacobian_DMLegacy,PETSC_NULL);}
585: return(0);
586: }