Actual source code: ms.c
petsc-3.14.6 2021-03-30
1: #include <petsc/private/snesimpl.h>
3: static SNESMSType SNESMSDefault = SNESMSM62;
4: static PetscBool SNESMSRegisterAllCalled;
5: static PetscBool SNESMSPackageInitialized;
7: typedef struct _SNESMSTableau *SNESMSTableau;
8: struct _SNESMSTableau {
9: char *name;
10: PetscInt nstages; /* Number of stages */
11: PetscInt nregisters; /* Number of registers */
12: PetscReal stability; /* Scaled stability region */
13: PetscReal *gamma; /* Coefficients of 3S* method */
14: PetscReal *delta; /* Coefficients of 3S* method */
15: PetscReal *betasub; /* Subdiagonal of beta in Shu-Osher form */
16: };
18: typedef struct _SNESMSTableauLink *SNESMSTableauLink;
19: struct _SNESMSTableauLink {
20: struct _SNESMSTableau tab;
21: SNESMSTableauLink next;
22: };
23: static SNESMSTableauLink SNESMSTableauList;
25: typedef struct {
26: SNESMSTableau tableau; /* Tableau in low-storage form */
27: PetscReal damping; /* Damping parameter, like length of (pseudo) time step */
28: PetscBool norms; /* Compute norms, usually only for monitoring purposes */
29: } SNES_MS;
31: /*@C
32: SNESMSRegisterAll - Registers all of the multi-stage methods in SNESMS
34: Not Collective, but should be called by all processes which will need the schemes to be registered
36: Level: advanced
38: .seealso: SNESMSRegisterDestroy()
39: @*/
40: PetscErrorCode SNESMSRegisterAll(void)
41: {
45: if (SNESMSRegisterAllCalled) return(0);
46: SNESMSRegisterAllCalled = PETSC_TRUE;
48: {
49: PetscReal alpha[1] = {1.0};
50: SNESMSRegister(SNESMSEULER,1,1,1.0,NULL,NULL,alpha);
51: }
53: {
54: PetscReal gamma[3][6] = {
55: {0.0000000000000000E+00, -7.0304722367110606E-01, -1.9836719667506464E-01, -1.6023843981863788E+00, 9.4483822882855284E-02, -1.4204296130641869E-01},
56: {1.0000000000000000E+00, 1.1111025767083920E+00, 5.6150921583923230E-01, 7.4151723494934041E-01, 3.1714538168600587E-01, 4.6479276238548706E-01},
57: {0.0000000000000000E+00, 0.0000000000000000E+00, 0.0000000000000000E+00, 6.7968174970583317E-01, -4.1755042846051737E-03, -1.9115668129923846E-01}
58: };
59: PetscReal delta[6] = {1.0000000000000000E+00, 5.3275427433201750E-01, 6.0143544663985238E-01, 4.5874077053842177E-01, 2.7544386906104651E-01, 0.0000000000000000E+00};
60: PetscReal betasub[6] = {8.4753115429481929E-01, 7.4018896368655618E-01, 6.5963574086583309E-03, 4.6747795645517759E-01, 1.3314545813643919E-01, 5.3260800028018784E-01};
61: SNESMSRegister(SNESMSM62,6,3,1.0,&gamma[0][0],delta,betasub);
62: }
64: { /* Jameson (1983) */
65: PetscReal alpha[4] = {0.25, 0.5, 0.55, 1.0};
66: SNESMSRegister(SNESMSJAMESON83,4,1,1.0,NULL,NULL,alpha);
67: }
69: { /* Van Leer, Tai, and Powell (1989) 1 stage, order 1 */
70: PetscReal alpha[1] = {1.0};
71: SNESMSRegister(SNESMSVLTP11,1,1,0.5,NULL,NULL,alpha);
72: }
73: { /* Van Leer, Tai, and Powell (1989) 2 stage, order 1 */
74: PetscReal alpha[2] = {0.3333, 1.0};
75: SNESMSRegister(SNESMSVLTP21,2,1,1.0,NULL,NULL,alpha);
76: }
77: { /* Van Leer, Tai, and Powell (1989) 3 stage, order 1 */
78: PetscReal alpha[3] = {0.1481, 0.4000, 1.0};
79: SNESMSRegister(SNESMSVLTP31,3,1,1.5,NULL,NULL,alpha);
80: }
81: { /* Van Leer, Tai, and Powell (1989) 4 stage, order 1 */
82: PetscReal alpha[4] = {0.0833, 0.2069, 0.4265, 1.0};
83: SNESMSRegister(SNESMSVLTP41,4,1,2.0,NULL,NULL,alpha);
84: }
85: { /* Van Leer, Tai, and Powell (1989) 5 stage, order 1 */
86: PetscReal alpha[5] = {0.0533, 0.1263, 0.2375, 0.4414,1.0};
87: SNESMSRegister(SNESMSVLTP51,5,1,2.5,NULL,NULL,alpha);
88: }
89: { /* Van Leer, Tai, and Powell (1989) 6 stage, order 1 */
90: PetscReal alpha[6] = {0.0370, 0.0851, 0.1521, 0.2562, 0.4512, 1.0};
91: SNESMSRegister(SNESMSVLTP61,6,1,3.0,NULL,NULL,alpha);
92: }
93: return(0);
94: }
96: /*@C
97: SNESMSRegisterDestroy - Frees the list of schemes that were registered by SNESMSRegister().
99: Not Collective
101: Level: advanced
103: .seealso: SNESMSRegister(), SNESMSRegisterAll()
104: @*/
105: PetscErrorCode SNESMSRegisterDestroy(void)
106: {
107: PetscErrorCode ierr;
108: SNESMSTableauLink link;
111: while ((link = SNESMSTableauList)) {
112: SNESMSTableau t = &link->tab;
113: SNESMSTableauList = link->next;
115: PetscFree(t->name);
116: PetscFree(t->gamma);
117: PetscFree(t->delta);
118: PetscFree(t->betasub);
119: PetscFree(link);
120: }
121: SNESMSRegisterAllCalled = PETSC_FALSE;
122: return(0);
123: }
125: /*@C
126: SNESMSInitializePackage - This function initializes everything in the SNESMS package. It is called
127: from SNESInitializePackage().
129: Level: developer
131: .seealso: PetscInitialize()
132: @*/
133: PetscErrorCode SNESMSInitializePackage(void)
134: {
138: if (SNESMSPackageInitialized) return(0);
139: SNESMSPackageInitialized = PETSC_TRUE;
141: SNESMSRegisterAll();
142: PetscRegisterFinalize(SNESMSFinalizePackage);
143: return(0);
144: }
146: /*@C
147: SNESMSFinalizePackage - This function destroys everything in the SNESMS package. It is
148: called from PetscFinalize().
150: Level: developer
152: .seealso: PetscFinalize()
153: @*/
154: PetscErrorCode SNESMSFinalizePackage(void)
155: {
159: SNESMSPackageInitialized = PETSC_FALSE;
161: SNESMSRegisterDestroy();
162: return(0);
163: }
165: /*@C
166: SNESMSRegister - register a multistage scheme
168: Not Collective, but the same schemes should be registered on all processes on which they will be used
170: Input Parameters:
171: + name - identifier for method
172: . nstages - number of stages
173: . nregisters - number of registers used by low-storage implementation
174: . stability - scaled stability region
175: . gamma - coefficients, see Ketcheson's paper
176: . delta - coefficients, see Ketcheson's paper
177: - betasub - subdiagonal of Shu-Osher form
179: Notes:
180: The notation is described in Ketcheson (2010) Runge-Kutta methods with minimum storage implementations.
182: Many multistage schemes are of the form
183: $ X_0 = X^{(n)}
184: $ X_k = X_0 + \alpha_k * F(X_{k-1}), k = 1,\ldots,s
185: $ X^{(n+1)} = X_s
186: These methods can be registered with
187: .vb
188: SNESMSRegister("name",s,1,stability,NULL,NULL,alpha);
189: .ve
191: Level: advanced
193: .seealso: SNESMS
194: @*/
195: PetscErrorCode SNESMSRegister(SNESMSType name,PetscInt nstages,PetscInt nregisters,PetscReal stability,const PetscReal gamma[],const PetscReal delta[],const PetscReal betasub[])
196: {
197: PetscErrorCode ierr;
198: SNESMSTableauLink link;
199: SNESMSTableau t;
203: if (nstages < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must have at least one stage");
204: if (gamma || delta) {
205: if (nregisters != 3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Only support for methods written in 3-register form");
208: } else {
209: if (nregisters != 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Only support for methods written in 1-register form");
210: }
213: SNESMSInitializePackage();
214: PetscNew(&link);
215: t = &link->tab;
216: PetscStrallocpy(name,&t->name);
217: t->nstages = nstages;
218: t->nregisters = nregisters;
219: t->stability = stability;
221: if (gamma && delta) {
222: PetscMalloc1(nstages*nregisters,&t->gamma);
223: PetscMalloc1(nstages,&t->delta);
224: PetscArraycpy(t->gamma,gamma,nstages*nregisters);
225: PetscArraycpy(t->delta,delta,nstages);
226: }
227: PetscMalloc1(nstages,&t->betasub);
228: PetscArraycpy(t->betasub,betasub,nstages);
230: link->next = SNESMSTableauList;
231: SNESMSTableauList = link;
232: return(0);
233: }
235: /*
236: X - initial state, updated in-place.
237: F - residual, computed at the initial X on input
238: */
239: static PetscErrorCode SNESMSStep_3Sstar(SNES snes,Vec X,Vec F)
240: {
241: PetscErrorCode ierr;
242: SNES_MS *ms = (SNES_MS*)snes->data;
243: SNESMSTableau t = ms->tableau;
244: const PetscReal *gamma = t->gamma,*delta = t->delta,*betasub = t->betasub;
245: Vec S1,S2,S3,Y;
246: PetscInt i,nstages = t->nstages;
249: Y = snes->work[0];
250: S1 = X;
251: S2 = snes->work[1];
252: S3 = snes->work[2];
253: VecZeroEntries(S2);
254: VecCopy(X,S3);
255: for (i = 0; i < nstages; i++) {
256: Vec Ss[4];
257: PetscScalar scoeff[4];
259: Ss[0] = S1; Ss[1] = S2; Ss[2] = S3; Ss[3] = Y;
261: scoeff[0] = gamma[0*nstages+i] - 1;
262: scoeff[1] = gamma[1*nstages+i];
263: scoeff[2] = gamma[2*nstages+i];
264: scoeff[3] = -betasub[i]*ms->damping;
266: VecAXPY(S2,delta[i],S1);
267: if (i > 0) {
268: SNESComputeFunction(snes,S1,F);
269: }
270: KSPSolve(snes->ksp,F,Y);
271: VecMAXPY(S1,4,scoeff,Ss);
272: }
273: return(0);
274: }
276: /*
277: X - initial state, updated in-place.
278: F - residual, computed at the initial X on input
279: */
280: static PetscErrorCode SNESMSStep_Basic(SNES snes,Vec X,Vec F)
281: {
282: PetscErrorCode ierr;
283: SNES_MS *ms = (SNES_MS*)snes->data;
284: SNESMSTableau tab = ms->tableau;
285: const PetscReal *alpha = tab->betasub, h = ms->damping;
286: PetscInt i,nstages = tab->nstages;
287: Vec X0 = snes->work[0];
290: VecCopy(X,X0);
291: for (i = 0; i < nstages; i++) {
292: if (i > 0) {
293: SNESComputeFunction(snes,X,F);
294: }
295: KSPSolve(snes->ksp,F,X);
296: VecAYPX(X,-alpha[i]*h,X0);
297: }
298: return(0);
299: }
301: static PetscErrorCode SNESMSStep_Step(SNES snes,Vec X,Vec F)
302: {
303: SNES_MS *ms = (SNES_MS*)snes->data;
304: SNESMSTableau tab = ms->tableau;
308: if (tab->gamma && tab->delta) {
309: SNESMSStep_3Sstar(snes,X,F);
310: } else {
311: SNESMSStep_Basic(snes,X,F);
312: }
313: return(0);
314: }
316: static PetscErrorCode SNESMSStep_Norms(SNES snes,PetscInt iter,Vec F)
317: {
318: SNES_MS *ms = (SNES_MS*)snes->data;
319: PetscReal fnorm;
323: if (ms->norms) {
324: VecNorm(F,NORM_2,&fnorm); /* fnorm <- ||F|| */
325: SNESCheckFunctionNorm(snes,fnorm);
326: /* Monitor convergence */
327: PetscObjectSAWsTakeAccess((PetscObject)snes);
328: snes->iter = iter;
329: snes->norm = fnorm;
330: PetscObjectSAWsGrantAccess((PetscObject)snes);
331: SNESLogConvergenceHistory(snes,snes->norm,0);
332: SNESMonitor(snes,snes->iter,snes->norm);
333: /* Test for convergence */
334: (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);
335: } else if (iter > 0) {
336: PetscObjectSAWsTakeAccess((PetscObject)snes);
337: snes->iter = iter;
338: PetscObjectSAWsGrantAccess((PetscObject)snes);
339: }
340: return(0);
341: }
343: static PetscErrorCode SNESSolve_MS(SNES snes)
344: {
345: SNES_MS *ms = (SNES_MS*)snes->data;
346: Vec X = snes->vec_sol,F = snes->vec_func;
347: PetscInt i;
351: if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name);
352: PetscCitationsRegister(SNESCitation,&SNEScite);
354: snes->reason = SNES_CONVERGED_ITERATING;
355: PetscObjectSAWsTakeAccess((PetscObject)snes);
356: snes->iter = 0;
357: snes->norm = 0;
358: PetscObjectSAWsGrantAccess((PetscObject)snes);
360: if (!snes->vec_func_init_set) {
361: SNESComputeFunction(snes,X,F);
362: } else snes->vec_func_init_set = PETSC_FALSE;
364: SNESMSStep_Norms(snes,0,F);
365: if (snes->reason) return(0);
367: for (i = 0; i < snes->max_its; i++) {
369: /* Call general purpose update function */
370: if (snes->ops->update) {
371: (*snes->ops->update)(snes,snes->iter);
372: }
374: if (i == 0 && snes->jacobian) {
375: /* This method does not require a Jacobian, but it is usually preconditioned by PBJacobi */
376: SNESComputeJacobian(snes,snes->vec_sol,snes->jacobian,snes->jacobian_pre);
377: SNESCheckJacobianDomainerror(snes);
378: KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);
379: }
381: SNESMSStep_Step(snes,X,F);
383: if (i < snes->max_its-1 || ms->norms) {
384: SNESComputeFunction(snes,X,F);
385: }
387: SNESMSStep_Norms(snes,i+1,F);
388: if (snes->reason) return(0);
389: }
390: if (!snes->reason) snes->reason = SNES_CONVERGED_ITS;
391: return(0);
392: }
394: static PetscErrorCode SNESSetUp_MS(SNES snes)
395: {
396: SNES_MS *ms = (SNES_MS*)snes->data;
397: SNESMSTableau tab = ms->tableau;
398: PetscInt nwork = tab->nregisters;
402: SNESSetWorkVecs(snes,nwork);
403: SNESSetUpMatrices(snes);
404: return(0);
405: }
407: static PetscErrorCode SNESReset_MS(SNES snes)
408: {
410: return(0);
411: }
413: static PetscErrorCode SNESDestroy_MS(SNES snes)
414: {
418: SNESReset_MS(snes);
419: PetscFree(snes->data);
420: PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetType_C",NULL);
421: PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",NULL);
422: PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetDamping_C",NULL);
423: PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetDamping_C",NULL);
424: return(0);
425: }
427: static PetscErrorCode SNESView_MS(SNES snes,PetscViewer viewer)
428: {
429: PetscBool iascii;
431: SNES_MS *ms = (SNES_MS*)snes->data;
432: SNESMSTableau tab = ms->tableau;
435: PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);
436: if (iascii) {
437: PetscViewerASCIIPrintf(viewer," multi-stage method type: %s\n",tab->name);
438: }
439: return(0);
440: }
442: static PetscErrorCode SNESSetFromOptions_MS(PetscOptionItems *PetscOptionsObject,SNES snes)
443: {
444: SNES_MS *ms = (SNES_MS*)snes->data;
448: PetscOptionsHead(PetscOptionsObject,"SNES MS options");
449: {
450: SNESMSTableauLink link;
451: PetscInt count,choice;
452: PetscBool flg;
453: const char **namelist;
454: SNESMSType mstype;
455: PetscReal damping;
457: SNESMSGetType(snes,&mstype);
458: for (link=SNESMSTableauList,count=0; link; link=link->next,count++) ;
459: PetscMalloc1(count,(char***)&namelist);
460: for (link=SNESMSTableauList,count=0; link; link=link->next,count++) namelist[count] = link->tab.name;
461: PetscOptionsEList("-snes_ms_type","Multistage smoother type","SNESMSSetType",(const char*const*)namelist,count,mstype,&choice,&flg);
462: if (flg) {SNESMSSetType(snes,namelist[choice]);}
463: PetscFree(namelist);
464: SNESMSGetDamping(snes,&damping);
465: PetscOptionsReal("-snes_ms_damping","Damping for multistage method","SNESMSSetDamping",damping,&damping,&flg);
466: if (flg) {SNESMSSetDamping(snes,damping);}
467: PetscOptionsBool("-snes_ms_norms","Compute norms for monitoring","none",ms->norms,&ms->norms,NULL);
468: }
469: PetscOptionsTail();
470: return(0);
471: }
473: static PetscErrorCode SNESMSGetType_MS(SNES snes,SNESMSType *mstype)
474: {
475: SNES_MS *ms = (SNES_MS*)snes->data;
476: SNESMSTableau tab = ms->tableau;
479: *mstype = tab->name;
480: return(0);
481: }
483: static PetscErrorCode SNESMSSetType_MS(SNES snes,SNESMSType mstype)
484: {
485: SNES_MS *ms = (SNES_MS*)snes->data;
486: SNESMSTableauLink link;
487: PetscBool match;
488: PetscErrorCode ierr;
491: if (ms->tableau) {
492: PetscStrcmp(ms->tableau->name,mstype,&match);
493: if (match) return(0);
494: }
495: for (link = SNESMSTableauList; link; link=link->next) {
496: PetscStrcmp(link->tab.name,mstype,&match);
497: if (match) {
498: if (snes->setupcalled) {SNESReset_MS(snes);}
499: ms->tableau = &link->tab;
500: if (snes->setupcalled) {SNESSetUp_MS(snes);}
501: return(0);
502: }
503: }
504: SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_UNKNOWN_TYPE,"Could not find '%s'",mstype);
505: }
507: /*@C
508: SNESMSGetType - Get the type of multistage smoother
510: Not collective
512: Input Parameter:
513: . snes - nonlinear solver context
515: Output Parameter:
516: . mstype - type of multistage method
518: Level: beginner
520: .seealso: SNESMSSetType(), SNESMSType, SNESMS
521: @*/
522: PetscErrorCode SNESMSGetType(SNES snes,SNESMSType *mstype)
523: {
529: PetscUseMethod(snes,"SNESMSGetType_C",(SNES,SNESMSType*),(snes,mstype));
530: return(0);
531: }
533: /*@C
534: SNESMSSetType - Set the type of multistage smoother
536: Logically collective
538: Input Parameter:
539: + snes - nonlinear solver context
540: - mstype - type of multistage method
542: Level: beginner
544: .seealso: SNESMSGetType(), SNESMSType, SNESMS
545: @*/
546: PetscErrorCode SNESMSSetType(SNES snes,SNESMSType mstype)
547: {
553: PetscTryMethod(snes,"SNESMSSetType_C",(SNES,SNESMSType),(snes,mstype));
554: return(0);
555: }
557: static PetscErrorCode SNESMSGetDamping_MS(SNES snes,PetscReal *damping)
558: {
559: SNES_MS *ms = (SNES_MS*)snes->data;
562: *damping = ms->damping;
563: return(0);
564: }
566: static PetscErrorCode SNESMSSetDamping_MS(SNES snes,PetscReal damping)
567: {
568: SNES_MS *ms = (SNES_MS*)snes->data;
571: ms->damping = damping;
572: return(0);
573: }
575: /*@C
576: SNESMSGetDamping - Get the damping parameter
578: Not collective
580: Input Parameter:
581: . snes - nonlinear solver context
583: Output Parameter:
584: . damping - damping parameter
586: Level: advanced
588: .seealso: SNESMSSetDamping(), SNESMS
589: @*/
590: PetscErrorCode SNESMSGetDamping(SNES snes,PetscReal *damping)
591: {
597: PetscUseMethod(snes,"SNESMSGetDamping_C",(SNES,PetscReal*),(snes,damping));
598: return(0);
599: }
601: /*@C
602: SNESMSSetDamping - Set the damping parameter
604: Logically collective
606: Input Parameter:
607: + snes - nonlinear solver context
608: - damping - damping parameter
610: Level: advanced
612: .seealso: SNESMSGetDamping(), SNESMS
613: @*/
614: PetscErrorCode SNESMSSetDamping(SNES snes,PetscReal damping)
615: {
621: PetscTryMethod(snes,"SNESMSSetDamping_C",(SNES,PetscReal),(snes,damping));
622: return(0);
623: }
625: /* -------------------------------------------------------------------------- */
626: /*MC
627: SNESMS - multi-stage smoothers
629: Options Database:
631: + -snes_ms_type - type of multi-stage smoother
632: - -snes_ms_damping - damping for multi-stage method
634: Notes:
635: These multistage methods are explicit Runge-Kutta methods that are often used as smoothers for
636: FAS multigrid for transport problems. In the linear case, these are equivalent to polynomial smoothers (such as Chebyshev).
638: Multi-stage smoothers should usually be preconditioned by point-block Jacobi to ensure proper scaling and to normalize the wave speeds.
640: The methods are specified in low storage form (Ketcheson 2010). New methods can be registered with SNESMSRegister().
642: References:
643: + 1. - Ketcheson (2010) Runge Kutta methods with minimum storage implementations (https://doi.org/10.1016/j.jcp.2009.11.006).
644: . 2. - Jameson (1983) Solution of the Euler equations for two dimensional transonic flow by a multigrid method (https://doi.org/10.1016/0096-3003(83)90019-X).
645: . 3. - Pierce and Giles (1997) Preconditioned multigrid methods for compressible flow calculations on stretched meshes (https://doi.org/10.1006/jcph.1997.5772).
646: - 4. - Van Leer, Tai, and Powell (1989) Design of optimally smoothing multi-stage schemes for the Euler equations (https://doi.org/10.2514/6.1989-1933).
648: Level: beginner
650: .seealso: SNESCreate(), SNES, SNESSetType(), SNESMS, SNESFAS, KSPCHEBYSHEV
652: M*/
653: PETSC_EXTERN PetscErrorCode SNESCreate_MS(SNES snes)
654: {
656: SNES_MS *ms;
659: SNESMSInitializePackage();
661: snes->ops->setup = SNESSetUp_MS;
662: snes->ops->solve = SNESSolve_MS;
663: snes->ops->destroy = SNESDestroy_MS;
664: snes->ops->setfromoptions = SNESSetFromOptions_MS;
665: snes->ops->view = SNESView_MS;
666: snes->ops->reset = SNESReset_MS;
668: snes->usesnpc = PETSC_FALSE;
669: snes->usesksp = PETSC_TRUE;
671: snes->alwayscomputesfinalresidual = PETSC_FALSE;
673: PetscNewLog(snes,&ms);
674: snes->data = (void*)ms;
675: ms->damping = 0.9;
676: ms->norms = PETSC_FALSE;
678: PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetType_C",SNESMSGetType_MS);
679: PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",SNESMSSetType_MS);
680: PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetDamping_C",SNESMSGetDamping_MS);
681: PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetDamping_C",SNESMSSetDamping_MS);
683: SNESMSSetType(snes,SNESMSDefault);
684: return(0);
685: }