Actual source code: dmsnes.c
petsc-3.5.4 2015-05-23
1: #include <petsc-private/snesimpl.h> /*I "petscsnes.h" I*/
2: #include <petsc-private/dmimpl.h> /*I "petscdm.h" I*/
6: static PetscErrorCode DMSNESDestroy(DMSNES *kdm)
7: {
11: if (!*kdm) return(0);
13: if (--((PetscObject)(*kdm))->refct > 0) {*kdm = 0; return(0);}
14: if ((*kdm)->ops->destroy) {((*kdm)->ops->destroy)(*kdm);}
15: PetscHeaderDestroy(kdm);
16: return(0);
17: }
21: PetscErrorCode DMSNESLoad(DMSNES kdm,PetscViewer viewer)
22: {
26: PetscViewerBinaryRead(viewer,&kdm->ops->computefunction,1,PETSC_FUNCTION);
27: PetscViewerBinaryRead(viewer,&kdm->ops->computejacobian,1,PETSC_FUNCTION);
28: return(0);
29: }
33: PetscErrorCode DMSNESView(DMSNES kdm,PetscViewer viewer)
34: {
36: PetscBool isascii,isbinary;
39: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
40: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
41: if (isascii) {
42: #if defined(PETSC_SERIALIZE_FUNCTIONS)
43: const char *fname;
45: PetscFPTFind(kdm->ops->computefunction,&fname);
46: if (fname) {
47: PetscViewerASCIIPrintf(viewer,"Function used by SNES: %s\n",fname);
48: }
49: PetscFPTFind(kdm->ops->computejacobian,&fname);
50: if (fname) {
51: PetscViewerASCIIPrintf(viewer,"Jacobian function used by SNES: %s\n",fname);
52: }
53: #endif
54: } else if (isbinary) {
55: struct {
56: PetscErrorCode (*func)(SNES,Vec,Vec,void*);
57: } funcstruct;
58: struct {
59: PetscErrorCode (*jac)(SNES,Vec,Mat,Mat,void*);
60: } jacstruct;
61: funcstruct.func = kdm->ops->computefunction;
62: jacstruct.jac = kdm->ops->computejacobian;
63: PetscViewerBinaryWrite(viewer,&funcstruct,1,PETSC_FUNCTION,PETSC_FALSE);
64: PetscViewerBinaryWrite(viewer,&jacstruct,1,PETSC_FUNCTION,PETSC_FALSE);
65: }
66: return(0);
67: }
71: static PetscErrorCode DMSNESCreate(MPI_Comm comm,DMSNES *kdm)
72: {
76: SNESInitializePackage();
77: PetscHeaderCreate(*kdm, _p_DMSNES, struct _DMSNESOps, DMSNES_CLASSID, "DMSNES", "DMSNES", "DMSNES", comm, DMSNESDestroy, DMSNESView);
78: PetscMemzero((*kdm)->ops, sizeof(struct _DMSNESOps));
79: return(0);
80: }
84: /* Attaches the DMSNES to the coarse level.
85: * Under what conditions should we copy versus duplicate?
86: */
87: static PetscErrorCode DMCoarsenHook_DMSNES(DM dm,DM dmc,void *ctx)
88: {
92: DMCopyDMSNES(dm,dmc);
93: return(0);
94: }
98: /* This could restrict auxiliary information to the coarse level.
99: */
100: static PetscErrorCode DMRestrictHook_DMSNES(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
101: {
104: return(0);
105: }
109: /* Attaches the DMSNES to the subdomain. */
110: static PetscErrorCode DMSubDomainHook_DMSNES(DM dm,DM subdm,void *ctx)
111: {
115: DMCopyDMSNES(dm,subdm);
116: return(0);
117: }
121: /* This could restrict auxiliary information to the coarse level.
122: */
123: static PetscErrorCode DMSubDomainRestrictHook_DMSNES(DM dm,VecScatter gscat,VecScatter lscat,DM subdm,void *ctx)
124: {
127: return(0);
128: }
132: static PetscErrorCode DMRefineHook_DMSNES(DM dm,DM dmf,void *ctx)
133: {
137: DMCopyDMSNES(dm,dmf);
138: return(0);
139: }
143: /* This could restrict auxiliary information to the coarse level.
144: */
145: static PetscErrorCode DMInterpolateHook_DMSNES(DM dm,Mat Interp,DM dmf,void *ctx)
146: {
149: return(0);
150: }
154: /*@C
155: DMSNESCopy - copies the information in a DMSNES to another DMSNES
157: Not Collective
159: Input Argument:
160: + kdm - Original DMSNES
161: - nkdm - DMSNES to receive the data, should have been created with DMSNESCreate()
163: Level: developer
165: .seealso: DMSNESCreate(), DMSNESDestroy()
166: @*/
167: PetscErrorCode DMSNESCopy(DMSNES kdm,DMSNES nkdm)
168: {
174: nkdm->ops->computefunction = kdm->ops->computefunction;
175: nkdm->ops->computejacobian = kdm->ops->computejacobian;
176: nkdm->ops->computegs = kdm->ops->computegs;
177: nkdm->ops->computeobjective = kdm->ops->computeobjective;
178: nkdm->ops->computepjacobian = kdm->ops->computepjacobian;
179: nkdm->ops->computepfunction = kdm->ops->computepfunction;
180: nkdm->ops->destroy = kdm->ops->destroy;
181: nkdm->ops->duplicate = kdm->ops->duplicate;
183: nkdm->functionctx = kdm->functionctx;
184: nkdm->gsctx = kdm->gsctx;
185: nkdm->pctx = kdm->pctx;
186: nkdm->jacobianctx = kdm->jacobianctx;
187: nkdm->objectivectx = kdm->objectivectx;
188: nkdm->data = kdm->data;
190: /*
191: nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
192: nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
193: nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
194: */
196: /* implementation specific copy hooks */
197: if (kdm->ops->duplicate) {(*kdm->ops->duplicate)(kdm,nkdm);}
198: return(0);
199: }
203: /*@C
204: DMGetDMSNES - get read-only private DMSNES context from a DM
206: Not Collective
208: Input Argument:
209: . dm - DM to be used with SNES
211: Output Argument:
212: . snesdm - private DMSNES context
214: Level: developer
216: Notes:
217: Use DMGetDMSNESWrite() if write access is needed. The DMSNESSetXXX API should be used wherever possible.
219: .seealso: DMGetDMSNESWrite()
220: @*/
221: PetscErrorCode DMGetDMSNES(DM dm,DMSNES *snesdm)
222: {
227: *snesdm = (DMSNES) dm->dmsnes;
228: if (!*snesdm) {
229: PetscInfo(dm,"Creating new DMSNES\n");
230: DMSNESCreate(PetscObjectComm((PetscObject)dm),snesdm);
232: dm->dmsnes = (PetscObject) *snesdm;
234: DMCoarsenHookAdd(dm,DMCoarsenHook_DMSNES,DMRestrictHook_DMSNES,NULL);
235: DMRefineHookAdd(dm,DMRefineHook_DMSNES,DMInterpolateHook_DMSNES,NULL);
236: DMSubDomainHookAdd(dm,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
237: }
238: return(0);
239: }
243: /*@C
244: DMGetDMSNESWrite - get write access to private DMSNES context from a DM
246: Not Collective
248: Input Argument:
249: . dm - DM to be used with SNES
251: Output Argument:
252: . snesdm - private DMSNES context
254: Level: developer
256: .seealso: DMGetDMSNES()
257: @*/
258: PetscErrorCode DMGetDMSNESWrite(DM dm,DMSNES *snesdm)
259: {
261: DMSNES sdm;
265: DMGetDMSNES(dm,&sdm);
266: if (!sdm->originaldm) sdm->originaldm = dm;
267: if (sdm->originaldm != dm) { /* Copy on write */
268: DMSNES oldsdm = sdm;
269: PetscInfo(dm,"Copying DMSNES due to write\n");
270: DMSNESCreate(PetscObjectComm((PetscObject)dm),&sdm);
271: DMSNESCopy(oldsdm,sdm);
272: DMSNESDestroy((DMSNES*)&dm->dmsnes);
273: dm->dmsnes = (PetscObject)sdm;
274: }
275: *snesdm = sdm;
276: return(0);
277: }
281: /*@C
282: DMCopyDMSNES - copies a DM context to a new DM
284: Logically Collective
286: Input Arguments:
287: + dmsrc - DM to obtain context from
288: - dmdest - DM to add context to
290: Level: developer
292: Note:
293: The context is copied by reference. This function does not ensure that a context exists.
295: .seealso: DMGetDMSNES(), SNESSetDM()
296: @*/
297: PetscErrorCode DMCopyDMSNES(DM dmsrc,DM dmdest)
298: {
304: DMSNESDestroy((DMSNES*)&dmdest->dmsnes);
306: dmdest->dmsnes = dmsrc->dmsnes;
308: PetscObjectReference(dmdest->dmsnes);
309: DMCoarsenHookAdd(dmdest,DMCoarsenHook_DMSNES,NULL,NULL);
310: DMRefineHookAdd(dmdest,DMRefineHook_DMSNES,NULL,NULL);
311: DMSubDomainHookAdd(dmdest,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
312: return(0);
313: }
317: /*@C
318: DMSNESSetFunction - set SNES residual evaluation function
320: Not Collective
322: Input Arguments:
323: + dm - DM to be used with SNES
324: . f - residual evaluation function; see SNESFunction for details
325: - ctx - context for residual evaluation
327: Level: advanced
329: Note:
330: SNESSetFunction() is normally used, but it calls this function internally because the user context is actually
331: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
332: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
334: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESFunction
335: @*/
336: PetscErrorCode DMSNESSetFunction(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
337: {
339: DMSNES sdm;
343: if (f || ctx) {
344: DMGetDMSNESWrite(dm,&sdm);
345: }
346: if (f) sdm->ops->computefunction = f;
347: if (ctx) sdm->functionctx = ctx;
348: return(0);
349: }
353: /*@C
354: DMSNESGetFunction - get SNES residual evaluation function
356: Not Collective
358: Input Argument:
359: . dm - DM to be used with SNES
361: Output Arguments:
362: + f - residual evaluation function; see SNESFunction for details
363: - ctx - context for residual evaluation
365: Level: advanced
367: Note:
368: SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
369: associated with the DM.
371: .seealso: DMSNESSetContext(), DMSNESSetFunction(), SNESSetFunction(), SNESFunction
372: @*/
373: PetscErrorCode DMSNESGetFunction(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
374: {
376: DMSNES sdm;
380: DMGetDMSNES(dm,&sdm);
381: if (f) *f = sdm->ops->computefunction;
382: if (ctx) *ctx = sdm->functionctx;
383: return(0);
384: }
388: /*@C
389: DMSNESSetObjective - set SNES objective evaluation function
391: Not Collective
393: Input Arguments:
394: + dm - DM to be used with SNES
395: . obj - objective evaluation function; see SNESObjectiveFunction for details
396: - ctx - context for residual evaluation
398: Level: advanced
400: .seealso: DMSNESSetContext(), SNESGetObjective(), DMSNESSetFunction()
401: @*/
402: PetscErrorCode DMSNESSetObjective(DM dm,PetscErrorCode (*obj)(SNES,Vec,PetscReal*,void*),void *ctx)
403: {
405: DMSNES sdm;
409: if (obj || ctx) {
410: DMGetDMSNESWrite(dm,&sdm);
411: }
412: if (obj) sdm->ops->computeobjective = obj;
413: if (ctx) sdm->objectivectx = ctx;
414: return(0);
415: }
419: /*@C
420: DMSNESGetObjective - get SNES objective evaluation function
422: Not Collective
424: Input Argument:
425: . dm - DM to be used with SNES
427: Output Arguments:
428: + obj- residual evaluation function; see SNESObjectiveFunction for details
429: - ctx - context for residual evaluation
431: Level: advanced
433: Note:
434: SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
435: associated with the DM.
437: .seealso: DMSNESSetContext(), DMSNESSetObjective(), SNESSetFunction()
438: @*/
439: PetscErrorCode DMSNESGetObjective(DM dm,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
440: {
442: DMSNES sdm;
446: DMGetDMSNES(dm,&sdm);
447: if (obj) *obj = sdm->ops->computeobjective;
448: if (ctx) *ctx = sdm->objectivectx;
449: return(0);
450: }
454: /*@C
455: DMSNESSetNGS - set SNES Gauss-Seidel relaxation function
457: Not Collective
459: Input Argument:
460: + dm - DM to be used with SNES
461: . f - relaxation function, see SNESGSFunction
462: - ctx - context for residual evaluation
464: Level: advanced
466: Note:
467: SNESSetNGS() is normally used, but it calls this function internally because the user context is actually
468: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
469: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
471: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), DMSNESSetFunction(), SNESGSFunction
472: @*/
473: PetscErrorCode DMSNESSetNGS(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
474: {
476: DMSNES sdm;
480: if (f || ctx) {
481: DMGetDMSNESWrite(dm,&sdm);
482: }
483: if (f) sdm->ops->computegs = f;
484: if (ctx) sdm->gsctx = ctx;
485: return(0);
486: }
490: /*@C
491: DMSNESGetNGS - get SNES Gauss-Seidel relaxation function
493: Not Collective
495: Input Argument:
496: . dm - DM to be used with SNES
498: Output Arguments:
499: + f - relaxation function which performs Gauss-Seidel sweeps, see SNESGSFunction
500: - ctx - context for residual evaluation
502: Level: advanced
504: Note:
505: SNESGetNGS() is normally used, but it calls this function internally because the user context is actually
506: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
507: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
509: .seealso: DMSNESSetContext(), SNESGetNGS(), DMSNESGetJacobian(), DMSNESGetFunction(), SNESNGSFunction
510: @*/
511: PetscErrorCode DMSNESGetNGS(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
512: {
514: DMSNES sdm;
518: DMGetDMSNES(dm,&sdm);
519: if (f) *f = sdm->ops->computegs;
520: if (ctx) *ctx = sdm->gsctx;
521: return(0);
522: }
526: /*@C
527: DMSNESSetJacobian - set SNES Jacobian evaluation function
529: Not Collective
531: Input Argument:
532: + dm - DM to be used with SNES
533: . J - Jacobian evaluation function
534: - ctx - context for residual evaluation
536: Level: advanced
538: Note:
539: SNESSetJacobian() is normally used, but it calls this function internally because the user context is actually
540: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
541: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
543: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESGetJacobian(), SNESSetJacobian(), SNESJacobianFunction
544: @*/
545: PetscErrorCode DMSNESSetJacobian(DM dm,PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
546: {
548: DMSNES sdm;
552: if (J || ctx) {
553: DMGetDMSNESWrite(dm,&sdm);
554: }
555: if (J) sdm->ops->computejacobian = J;
556: if (ctx) sdm->jacobianctx = ctx;
557: return(0);
558: }
562: /*@C
563: DMSNESGetJacobian - get SNES Jacobian evaluation function
565: Not Collective
567: Input Argument:
568: . dm - DM to be used with SNES
570: Output Arguments:
571: + J - Jacobian evaluation function; see SNESJacobianFunction for all calling sequence
572: - ctx - context for residual evaluation
574: Level: advanced
576: Note:
577: SNESGetJacobian() is normally used, but it calls this function internally because the user context is actually
578: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
579: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
581: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESJacobianFunction
582: @*/
583: PetscErrorCode DMSNESGetJacobian(DM dm,PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
584: {
586: DMSNES sdm;
590: DMGetDMSNES(dm,&sdm);
591: if (J) *J = sdm->ops->computejacobian;
592: if (ctx) *ctx = sdm->jacobianctx;
593: return(0);
594: }
598: /*@C
599: DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions.
601: Not Collective
603: Input Argument:
604: + dm - DM to be used with SNES
605: . b - RHS evaluation function
606: . J - Picard matrix evaluation function
607: - ctx - context for residual evaluation
609: Level: advanced
611: .seealso: SNESSetPicard(), DMSNESSetFunction(), DMSNESSetJacobian()
612: @*/
613: PetscErrorCode DMSNESSetPicard(DM dm,PetscErrorCode (*b)(SNES,Vec,Vec,void*),PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
614: {
616: DMSNES sdm;
620: DMGetDMSNES(dm,&sdm);
621: if (b) sdm->ops->computepfunction = b;
622: if (J) sdm->ops->computepjacobian = J;
623: if (ctx) sdm->pctx = ctx;
624: return(0);
625: }
629: /*@C
630: DMSNESGetPicard - get SNES Picard iteration evaluation functions
632: Not Collective
634: Input Argument:
635: . dm - DM to be used with SNES
637: Output Arguments:
638: + b - RHS evaluation function; see SNESFunction for details
639: . J - RHS evaluation function; see SNESJacobianFunction for detailsa
640: - ctx - context for residual evaluation
642: Level: advanced
644: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
645: @*/
646: PetscErrorCode DMSNESGetPicard(DM dm,PetscErrorCode (**b)(SNES,Vec,Vec,void*),PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
647: {
649: DMSNES sdm;
653: DMGetDMSNES(dm,&sdm);
654: if (b) *b = sdm->ops->computepfunction;
655: if (J) *J = sdm->ops->computepjacobian;
656: if (ctx) *ctx = sdm->pctx;
657: return(0);
658: }