Actual source code: dmsnes.c
petsc-3.11.0 2019-03-29
1: #include <petsc/private/snesimpl.h>
2: #include <petsc/private/dmimpl.h>
4: static PetscErrorCode DMSNESDestroy(DMSNES *kdm)
5: {
9: if (!*kdm) return(0);
11: if (--((PetscObject)(*kdm))->refct > 0) {*kdm = 0; return(0);}
12: if ((*kdm)->ops->destroy) {((*kdm)->ops->destroy)(*kdm);}
13: PetscHeaderDestroy(kdm);
14: return(0);
15: }
17: PetscErrorCode DMSNESLoad(DMSNES kdm,PetscViewer viewer)
18: {
22: PetscViewerBinaryRead(viewer,&kdm->ops->computefunction,1,NULL,PETSC_FUNCTION);
23: PetscViewerBinaryRead(viewer,&kdm->ops->computejacobian,1,NULL,PETSC_FUNCTION);
24: return(0);
25: }
27: PetscErrorCode DMSNESView(DMSNES kdm,PetscViewer viewer)
28: {
30: PetscBool isascii,isbinary;
33: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&isascii);
34: PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERBINARY,&isbinary);
35: if (isascii) {
36: #if defined(PETSC_SERIALIZE_FUNCTIONS) && defined(PETSC_SERIALIZE_FUNCTIONS_VIEW)
37: const char *fname;
39: PetscFPTFind(kdm->ops->computefunction,&fname);
40: if (fname) {
41: PetscViewerASCIIPrintf(viewer,"Function used by SNES: %s\n",fname);
42: }
43: PetscFPTFind(kdm->ops->computejacobian,&fname);
44: if (fname) {
45: PetscViewerASCIIPrintf(viewer,"Jacobian function used by SNES: %s\n",fname);
46: }
47: #endif
48: } else if (isbinary) {
49: struct {
50: PetscErrorCode (*func)(SNES,Vec,Vec,void*);
51: } funcstruct;
52: struct {
53: PetscErrorCode (*jac)(SNES,Vec,Mat,Mat,void*);
54: } jacstruct;
55: funcstruct.func = kdm->ops->computefunction;
56: jacstruct.jac = kdm->ops->computejacobian;
57: PetscViewerBinaryWrite(viewer,&funcstruct,1,PETSC_FUNCTION,PETSC_FALSE);
58: PetscViewerBinaryWrite(viewer,&jacstruct,1,PETSC_FUNCTION,PETSC_FALSE);
59: }
60: return(0);
61: }
63: static PetscErrorCode DMSNESCreate(MPI_Comm comm,DMSNES *kdm)
64: {
68: SNESInitializePackage();
69: PetscHeaderCreate(*kdm, DMSNES_CLASSID, "DMSNES", "DMSNES", "DMSNES", comm, DMSNESDestroy, DMSNESView);
70: return(0);
71: }
73: /* Attaches the DMSNES to the coarse level.
74: * Under what conditions should we copy versus duplicate?
75: */
76: static PetscErrorCode DMCoarsenHook_DMSNES(DM dm,DM dmc,void *ctx)
77: {
81: DMCopyDMSNES(dm,dmc);
82: return(0);
83: }
85: /* This could restrict auxiliary information to the coarse level.
86: */
87: static PetscErrorCode DMRestrictHook_DMSNES(DM dm,Mat Restrict,Vec rscale,Mat Inject,DM dmc,void *ctx)
88: {
91: return(0);
92: }
94: /* Attaches the DMSNES to the subdomain. */
95: static PetscErrorCode DMSubDomainHook_DMSNES(DM dm,DM subdm,void *ctx)
96: {
100: DMCopyDMSNES(dm,subdm);
101: return(0);
102: }
104: /* This could restrict auxiliary information to the coarse level.
105: */
106: static PetscErrorCode DMSubDomainRestrictHook_DMSNES(DM dm,VecScatter gscat,VecScatter lscat,DM subdm,void *ctx)
107: {
110: return(0);
111: }
113: static PetscErrorCode DMRefineHook_DMSNES(DM dm,DM dmf,void *ctx)
114: {
118: DMCopyDMSNES(dm,dmf);
119: return(0);
120: }
122: /* This could restrict auxiliary information to the coarse level.
123: */
124: static PetscErrorCode DMInterpolateHook_DMSNES(DM dm,Mat Interp,DM dmf,void *ctx)
125: {
128: return(0);
129: }
131: /*@C
132: DMSNESCopy - copies the information in a DMSNES to another DMSNES
134: Not Collective
136: Input Argument:
137: + kdm - Original DMSNES
138: - nkdm - DMSNES to receive the data, should have been created with DMSNESCreate()
140: Level: developer
142: .seealso: DMSNESCreate(), DMSNESDestroy()
143: @*/
144: PetscErrorCode DMSNESCopy(DMSNES kdm,DMSNES nkdm)
145: {
151: nkdm->ops->computefunction = kdm->ops->computefunction;
152: nkdm->ops->computejacobian = kdm->ops->computejacobian;
153: nkdm->ops->computegs = kdm->ops->computegs;
154: nkdm->ops->computeobjective = kdm->ops->computeobjective;
155: nkdm->ops->computepjacobian = kdm->ops->computepjacobian;
156: nkdm->ops->computepfunction = kdm->ops->computepfunction;
157: nkdm->ops->destroy = kdm->ops->destroy;
158: nkdm->ops->duplicate = kdm->ops->duplicate;
160: nkdm->functionctx = kdm->functionctx;
161: nkdm->gsctx = kdm->gsctx;
162: nkdm->pctx = kdm->pctx;
163: nkdm->jacobianctx = kdm->jacobianctx;
164: nkdm->objectivectx = kdm->objectivectx;
166: /*
167: nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
168: nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
169: nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
170: */
172: /* implementation specific copy hooks */
173: if (kdm->ops->duplicate) {(*kdm->ops->duplicate)(kdm,nkdm);}
174: return(0);
175: }
177: /*@C
178: DMGetDMSNES - get read-only private DMSNES context from a DM
180: Not Collective
182: Input Argument:
183: . dm - DM to be used with SNES
185: Output Argument:
186: . snesdm - private DMSNES context
188: Level: developer
190: Notes:
191: Use DMGetDMSNESWrite() if write access is needed. The DMSNESSetXXX API should be used wherever possible.
193: .seealso: DMGetDMSNESWrite()
194: @*/
195: PetscErrorCode DMGetDMSNES(DM dm,DMSNES *snesdm)
196: {
201: *snesdm = (DMSNES) dm->dmsnes;
202: if (!*snesdm) {
203: PetscInfo(dm,"Creating new DMSNES\n");
204: DMSNESCreate(PetscObjectComm((PetscObject)dm),snesdm);
206: dm->dmsnes = (PetscObject) *snesdm;
208: DMCoarsenHookAdd(dm,DMCoarsenHook_DMSNES,DMRestrictHook_DMSNES,NULL);
209: DMRefineHookAdd(dm,DMRefineHook_DMSNES,DMInterpolateHook_DMSNES,NULL);
210: DMSubDomainHookAdd(dm,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
211: }
212: return(0);
213: }
215: /*@C
216: DMGetDMSNESWrite - get write access to private DMSNES context from a DM
218: Not Collective
220: Input Argument:
221: . dm - DM to be used with SNES
223: Output Argument:
224: . snesdm - private DMSNES context
226: Level: developer
228: .seealso: DMGetDMSNES()
229: @*/
230: PetscErrorCode DMGetDMSNESWrite(DM dm,DMSNES *snesdm)
231: {
233: DMSNES sdm;
237: DMGetDMSNES(dm,&sdm);
238: if (!sdm->originaldm) sdm->originaldm = dm;
239: if (sdm->originaldm != dm) { /* Copy on write */
240: DMSNES oldsdm = sdm;
241: PetscInfo(dm,"Copying DMSNES due to write\n");
242: DMSNESCreate(PetscObjectComm((PetscObject)dm),&sdm);
243: DMSNESCopy(oldsdm,sdm);
244: DMSNESDestroy((DMSNES*)&dm->dmsnes);
245: dm->dmsnes = (PetscObject)sdm;
246: }
247: *snesdm = sdm;
248: return(0);
249: }
251: /*@C
252: DMCopyDMSNES - copies a DM context to a new DM
254: Logically Collective
256: Input Arguments:
257: + dmsrc - DM to obtain context from
258: - dmdest - DM to add context to
260: Level: developer
262: Note:
263: The context is copied by reference. This function does not ensure that a context exists.
265: .seealso: DMGetDMSNES(), SNESSetDM()
266: @*/
267: PetscErrorCode DMCopyDMSNES(DM dmsrc,DM dmdest)
268: {
274: if (!dmdest->dmsnes) {DMSNESCreate(PetscObjectComm((PetscObject) dmdest), (DMSNES *) &dmdest->dmsnes);}
275: DMSNESCopy((DMSNES) dmsrc->dmsnes, (DMSNES) dmdest->dmsnes);
276: DMCoarsenHookAdd(dmdest,DMCoarsenHook_DMSNES,NULL,NULL);
277: DMRefineHookAdd(dmdest,DMRefineHook_DMSNES,NULL,NULL);
278: DMSubDomainHookAdd(dmdest,DMSubDomainHook_DMSNES,DMSubDomainRestrictHook_DMSNES,NULL);
279: return(0);
280: }
282: /*@C
283: DMSNESSetFunction - set SNES residual evaluation function
285: Not Collective
287: Input Arguments:
288: + dm - DM to be used with SNES
289: . f - residual evaluation function; see SNESFunction for details
290: - ctx - context for residual evaluation
292: Level: advanced
294: Note:
295: SNESSetFunction() is normally used, but it calls this function internally because the user context is actually
296: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
297: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
299: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESFunction
300: @*/
301: PetscErrorCode DMSNESSetFunction(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
302: {
304: DMSNES sdm;
308: if (f || ctx) {
309: DMGetDMSNESWrite(dm,&sdm);
310: }
311: if (f) sdm->ops->computefunction = f;
312: if (ctx) sdm->functionctx = ctx;
313: return(0);
314: }
316: /*@C
317: DMSNESGetFunction - get SNES residual evaluation function
319: Not Collective
321: Input Argument:
322: . dm - DM to be used with SNES
324: Output Arguments:
325: + f - residual evaluation function; see SNESFunction for details
326: - ctx - context for residual evaluation
328: Level: advanced
330: Note:
331: SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
332: associated with the DM.
334: .seealso: DMSNESSetContext(), DMSNESSetFunction(), SNESSetFunction(), SNESFunction
335: @*/
336: PetscErrorCode DMSNESGetFunction(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
337: {
339: DMSNES sdm;
343: DMGetDMSNES(dm,&sdm);
344: if (f) *f = sdm->ops->computefunction;
345: if (ctx) *ctx = sdm->functionctx;
346: return(0);
347: }
349: /*@C
350: DMSNESSetObjective - set SNES objective evaluation function
352: Not Collective
354: Input Arguments:
355: + dm - DM to be used with SNES
356: . obj - objective evaluation function; see SNESObjectiveFunction for details
357: - ctx - context for residual evaluation
359: Level: advanced
361: .seealso: DMSNESSetContext(), SNESGetObjective(), DMSNESSetFunction()
362: @*/
363: PetscErrorCode DMSNESSetObjective(DM dm,PetscErrorCode (*obj)(SNES,Vec,PetscReal*,void*),void *ctx)
364: {
366: DMSNES sdm;
370: if (obj || ctx) {
371: DMGetDMSNESWrite(dm,&sdm);
372: }
373: if (obj) sdm->ops->computeobjective = obj;
374: if (ctx) sdm->objectivectx = ctx;
375: return(0);
376: }
378: /*@C
379: DMSNESGetObjective - get SNES objective evaluation function
381: Not Collective
383: Input Argument:
384: . dm - DM to be used with SNES
386: Output Arguments:
387: + obj- residual evaluation function; see SNESObjectiveFunction for details
388: - ctx - context for residual evaluation
390: Level: advanced
392: Note:
393: SNESGetFunction() is normally used, but it calls this function internally because the user context is actually
394: associated with the DM.
396: .seealso: DMSNESSetContext(), DMSNESSetObjective(), SNESSetFunction()
397: @*/
398: PetscErrorCode DMSNESGetObjective(DM dm,PetscErrorCode (**obj)(SNES,Vec,PetscReal*,void*),void **ctx)
399: {
401: DMSNES sdm;
405: DMGetDMSNES(dm,&sdm);
406: if (obj) *obj = sdm->ops->computeobjective;
407: if (ctx) *ctx = sdm->objectivectx;
408: return(0);
409: }
411: /*@C
412: DMSNESSetNGS - set SNES Gauss-Seidel relaxation function
414: Not Collective
416: Input Argument:
417: + dm - DM to be used with SNES
418: . f - relaxation function, see SNESGSFunction
419: - ctx - context for residual evaluation
421: Level: advanced
423: Note:
424: SNESSetNGS() is normally used, but it calls this function internally because the user context is actually
425: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
426: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
428: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), DMSNESSetFunction(), SNESGSFunction
429: @*/
430: PetscErrorCode DMSNESSetNGS(DM dm,PetscErrorCode (*f)(SNES,Vec,Vec,void*),void *ctx)
431: {
433: DMSNES sdm;
437: if (f || ctx) {
438: DMGetDMSNESWrite(dm,&sdm);
439: }
440: if (f) sdm->ops->computegs = f;
441: if (ctx) sdm->gsctx = ctx;
442: return(0);
443: }
445: /*@C
446: DMSNESGetNGS - get SNES Gauss-Seidel relaxation function
448: Not Collective
450: Input Argument:
451: . dm - DM to be used with SNES
453: Output Arguments:
454: + f - relaxation function which performs Gauss-Seidel sweeps, see SNESGSFunction
455: - ctx - context for residual evaluation
457: Level: advanced
459: Note:
460: SNESGetNGS() is normally used, but it calls this function internally because the user context is actually
461: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
462: not. If DM took a more central role at some later date, this could become the primary method of setting the residual.
464: .seealso: DMSNESSetContext(), SNESGetNGS(), DMSNESGetJacobian(), DMSNESGetFunction(), SNESNGSFunction
465: @*/
466: PetscErrorCode DMSNESGetNGS(DM dm,PetscErrorCode (**f)(SNES,Vec,Vec,void*),void **ctx)
467: {
469: DMSNES sdm;
473: DMGetDMSNES(dm,&sdm);
474: if (f) *f = sdm->ops->computegs;
475: if (ctx) *ctx = sdm->gsctx;
476: return(0);
477: }
479: /*@C
480: DMSNESSetJacobian - set SNES Jacobian evaluation function
482: Not Collective
484: Input Argument:
485: + dm - DM to be used with SNES
486: . J - Jacobian evaluation function
487: - ctx - context for residual evaluation
489: Level: advanced
491: Note:
492: SNESSetJacobian() is normally used, but it calls this function internally because the user context is actually
493: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
494: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
496: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESGetJacobian(), SNESSetJacobian(), SNESJacobianFunction
497: @*/
498: PetscErrorCode DMSNESSetJacobian(DM dm,PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
499: {
501: DMSNES sdm;
505: if (J || ctx) {
506: DMGetDMSNESWrite(dm,&sdm);
507: }
508: if (J) sdm->ops->computejacobian = J;
509: if (ctx) sdm->jacobianctx = ctx;
510: return(0);
511: }
513: /*@C
514: DMSNESGetJacobian - get SNES Jacobian evaluation function
516: Not Collective
518: Input Argument:
519: . dm - DM to be used with SNES
521: Output Arguments:
522: + J - Jacobian evaluation function; see SNESJacobianFunction for all calling sequence
523: - ctx - context for residual evaluation
525: Level: advanced
527: Note:
528: SNESGetJacobian() is normally used, but it calls this function internally because the user context is actually
529: associated with the DM. This makes the interface consistent regardless of whether the user interacts with a DM or
530: not. If DM took a more central role at some later date, this could become the primary method of setting the Jacobian.
532: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian(), SNESJacobianFunction
533: @*/
534: PetscErrorCode DMSNESGetJacobian(DM dm,PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
535: {
537: DMSNES sdm;
541: DMGetDMSNES(dm,&sdm);
542: if (J) *J = sdm->ops->computejacobian;
543: if (ctx) *ctx = sdm->jacobianctx;
544: return(0);
545: }
547: /*@C
548: DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions.
550: Not Collective
552: Input Argument:
553: + dm - DM to be used with SNES
554: . b - RHS evaluation function
555: . J - Picard matrix evaluation function
556: - ctx - context for residual evaluation
558: Level: advanced
560: .seealso: SNESSetPicard(), DMSNESSetFunction(), DMSNESSetJacobian()
561: @*/
562: PetscErrorCode DMSNESSetPicard(DM dm,PetscErrorCode (*b)(SNES,Vec,Vec,void*),PetscErrorCode (*J)(SNES,Vec,Mat,Mat,void*),void *ctx)
563: {
565: DMSNES sdm;
569: DMGetDMSNES(dm,&sdm);
570: if (b) sdm->ops->computepfunction = b;
571: if (J) sdm->ops->computepjacobian = J;
572: if (ctx) sdm->pctx = ctx;
573: return(0);
574: }
576: /*@C
577: DMSNESGetPicard - get SNES Picard iteration evaluation functions
579: Not Collective
581: Input Argument:
582: . dm - DM to be used with SNES
584: Output Arguments:
585: + b - RHS evaluation function; see SNESFunction for details
586: . J - RHS evaluation function; see SNESJacobianFunction for detailsa
587: - ctx - context for residual evaluation
589: Level: advanced
591: .seealso: DMSNESSetContext(), SNESSetFunction(), DMSNESSetJacobian()
592: @*/
593: PetscErrorCode DMSNESGetPicard(DM dm,PetscErrorCode (**b)(SNES,Vec,Vec,void*),PetscErrorCode (**J)(SNES,Vec,Mat,Mat,void*),void **ctx)
594: {
596: DMSNES sdm;
600: DMGetDMSNES(dm,&sdm);
601: if (b) *b = sdm->ops->computepfunction;
602: if (J) *J = sdm->ops->computepjacobian;
603: if (ctx) *ctx = sdm->pctx;
604: return(0);
605: }