xref: /petsc/src/snes/impls/ms/ms.c (revision 55e7fe800d976e85ed2b5cd8bfdef564daa37bd9)
1 #include <petsc/private/snesimpl.h>   /*I "petscsnes.h" I*/
2 
3 static SNESMSType SNESMSDefault = SNESMSM62;
4 static PetscBool  SNESMSRegisterAllCalled;
5 static PetscBool  SNESMSPackageInitialized;
6 
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 };
17 
18 typedef struct _SNESMSTableauLink *SNESMSTableauLink;
19 struct _SNESMSTableauLink {
20   struct _SNESMSTableau tab;
21   SNESMSTableauLink     next;
22 };
23 static SNESMSTableauLink SNESMSTableauList;
24 
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;
30 
31 /*@C
32   SNESMSRegisterAll - Registers all of the multi-stage methods in SNESMS
33 
34   Not Collective, but should be called by all processes which will need the schemes to be registered
35 
36   Level: advanced
37 
38 .keywords: SNES, SNESMS, register, all
39 
40 .seealso:  SNESMSRegisterDestroy()
41 @*/
42 PetscErrorCode SNESMSRegisterAll(void)
43 {
44   PetscErrorCode ierr;
45 
46   PetscFunctionBegin;
47   if (SNESMSRegisterAllCalled) PetscFunctionReturn(0);
48   SNESMSRegisterAllCalled = PETSC_TRUE;
49 
50   {
51     PetscReal gamma[3][1] = {{1.0},{0.0},{0.0}};
52     PetscReal delta[1]    = {0.0};
53     PetscReal betasub[1]  = {1.0};
54     ierr = SNESMSRegister(SNESMSEULER,1,3,1.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
55   }
56 
57   {
58     PetscReal gamma[3][6] = {
59       {0.0000000000000000E+00,  -7.0304722367110606E-01, -1.9836719667506464E-01, -1.6023843981863788E+00, 9.4483822882855284E-02,  -1.4204296130641869E-01},
60       {1.0000000000000000E+00,  1.1111025767083920E+00,  5.6150921583923230E-01,  7.4151723494934041E-01,  3.1714538168600587E-01,  4.6479276238548706E-01},
61       {0.0000000000000000E+00,  0.0000000000000000E+00,  0.0000000000000000E+00,  6.7968174970583317E-01,  -4.1755042846051737E-03, -1.9115668129923846E-01}
62     };
63     PetscReal delta[6]   = {1.0000000000000000E+00, 5.3275427433201750E-01, 6.0143544663985238E-01, 4.5874077053842177E-01, 2.7544386906104651E-01, 0.0000000000000000E+00};
64     PetscReal betasub[6] = {8.4753115429481929E-01, 7.4018896368655618E-01, 6.5963574086583309E-03, 4.6747795645517759E-01, 1.3314545813643919E-01, 5.3260800028018784E-01};
65     ierr = SNESMSRegister(SNESMSM62,6,3,1.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
66   }
67   {
68     PetscReal gamma[3][4] = {{0,0,0,0},{0,0,0,0},{1,1,1,1}};
69     PetscReal delta[4]    = {0,0,0,0};
70     PetscReal betasub[4]  = {0.25, 0.5, 0.55, 1.0};
71     ierr = SNESMSRegister(SNESMSJAMESON83,4,3,1.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
72   }
73   {                             /* Van Leer, Tai, and Powell (1989) 2 stage, order 1 */
74     PetscReal gamma[3][2] = {{0,0},{0,0},{1,1}};
75     PetscReal delta[2]    = {0,0};
76     PetscReal betasub[2]  = {0.3333,1.0};
77     ierr = SNESMSRegister(SNESMSVLTP21,2,3,1.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
78   }
79   {                             /* Van Leer, Tai, and Powell (1989) 3 stage, order 1 */
80     PetscReal gamma[3][3] = {{0,0,0},{0,0,0},{1,1,1}};
81     PetscReal delta[3]    = {0,0,0};
82     PetscReal betasub[3]  = {0.1481,0.4000,1.0};
83     ierr = SNESMSRegister(SNESMSVLTP31,3,3,1.5,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
84   }
85   {                             /* Van Leer, Tai, and Powell (1989) 4 stage, order 1 */
86     PetscReal gamma[3][4] = {{0,0,0,0},{0,0,0,0},{1,1,1,1}};
87     PetscReal delta[4]    = {0,0,0,0};
88     PetscReal betasub[4]  = {0.0833,0.2069,0.4265,1.0};
89     ierr = SNESMSRegister(SNESMSVLTP41,4,3,2.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
90   }
91   {                             /* Van Leer, Tai, and Powell (1989) 5 stage, order 1 */
92     PetscReal gamma[3][5] = {{0,0,0,0,0},{0,0,0,0,0},{1,1,1,1,1}};
93     PetscReal delta[5]    = {0,0,0,0,0};
94     PetscReal betasub[5]  = {0.0533,0.1263,0.2375,0.4414,1.0};
95     ierr = SNESMSRegister(SNESMSVLTP51,5,3,2.5,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
96   }
97   {                             /* Van Leer, Tai, and Powell (1989) 6 stage, order 1 */
98     PetscReal gamma[3][6] = {{0,0,0,0,0,0},{0,0,0,0,0,0},{1,1,1,1,1,1}};
99     PetscReal delta[6]    = {0,0,0,0,0,0};
100     PetscReal betasub[6]  = {0.0370,0.0851,0.1521,0.2562,0.4512,1.0};
101     ierr = SNESMSRegister(SNESMSVLTP61,6,3,3.0,&gamma[0][0],delta,betasub);CHKERRQ(ierr);
102   }
103   PetscFunctionReturn(0);
104 }
105 
106 /*@C
107    SNESMSRegisterDestroy - Frees the list of schemes that were registered by TSRosWRegister().
108 
109    Not Collective
110 
111    Level: advanced
112 
113 .keywords: TSRosW, register, destroy
114 .seealso: TSRosWRegister(), TSRosWRegisterAll(), TSRosWRegister()
115 @*/
116 PetscErrorCode SNESMSRegisterDestroy(void)
117 {
118   PetscErrorCode    ierr;
119   SNESMSTableauLink link;
120 
121   PetscFunctionBegin;
122   while ((link = SNESMSTableauList)) {
123     SNESMSTableau t = &link->tab;
124     SNESMSTableauList = link->next;
125 
126     ierr = PetscFree3(t->gamma,t->delta,t->betasub);CHKERRQ(ierr);
127     ierr = PetscFree(t->name);CHKERRQ(ierr);
128     ierr = PetscFree(link);CHKERRQ(ierr);
129   }
130   SNESMSRegisterAllCalled = PETSC_FALSE;
131   PetscFunctionReturn(0);
132 }
133 
134 /*@C
135   SNESMSInitializePackage - This function initializes everything in the SNESMS package. It is called
136   from PetscDLLibraryRegister() when using dynamic libraries, and on the first call to SNESCreate_MS()
137   when using static libraries.
138 
139   Level: developer
140 
141 .keywords: SNES, SNESMS, initialize, package
142 .seealso: PetscInitialize()
143 @*/
144 PetscErrorCode SNESMSInitializePackage(void)
145 {
146   PetscErrorCode ierr;
147 
148   PetscFunctionBegin;
149   if (SNESMSPackageInitialized) PetscFunctionReturn(0);
150   SNESMSPackageInitialized = PETSC_TRUE;
151 
152   ierr = SNESMSRegisterAll();CHKERRQ(ierr);
153   ierr = PetscRegisterFinalize(SNESMSFinalizePackage);CHKERRQ(ierr);
154   PetscFunctionReturn(0);
155 }
156 
157 /*@C
158   SNESMSFinalizePackage - This function destroys everything in the SNESMS package. It is
159   called from PetscFinalize().
160 
161   Level: developer
162 
163 .keywords: Petsc, destroy, package
164 .seealso: PetscFinalize()
165 @*/
166 PetscErrorCode SNESMSFinalizePackage(void)
167 {
168   PetscErrorCode ierr;
169 
170   PetscFunctionBegin;
171   SNESMSPackageInitialized = PETSC_FALSE;
172 
173   ierr = SNESMSRegisterDestroy();CHKERRQ(ierr);
174   PetscFunctionReturn(0);
175 }
176 
177 /*@C
178    SNESMSRegister - register a multistage scheme
179 
180    Not Collective, but the same schemes should be registered on all processes on which they will be used
181 
182    Input Parameters:
183 +  name - identifier for method
184 .  nstages - number of stages
185 .  nregisters - number of registers used by low-storage implementation
186 .  gamma - coefficients, see Ketcheson's paper
187 .  delta - coefficients, see Ketcheson's paper
188 -  betasub - subdiagonal of Shu-Osher form
189 
190    Notes:
191    The notation is described in Ketcheson (2010) Runge-Kutta methods with minimum storage implementations.
192 
193    Level: advanced
194 
195 .keywords: SNES, register
196 
197 .seealso: SNESMS
198 @*/
199 PetscErrorCode SNESMSRegister(SNESMSType name,PetscInt nstages,PetscInt nregisters,PetscReal stability,const PetscReal gamma[],const PetscReal delta[],const PetscReal betasub[])
200 {
201   PetscErrorCode    ierr;
202   SNESMSTableauLink link;
203   SNESMSTableau     t;
204 
205   PetscFunctionBegin;
206   PetscValidCharPointer(name,1);
207   if (nstages < 1) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_OUTOFRANGE,"Must have at least one stage");
208   if (nregisters != 3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_SUP,"Only support for methods written in 3-register form");
209   PetscValidPointer(gamma,4);
210   PetscValidPointer(delta,5);
211   PetscValidPointer(betasub,6);
212 
213   ierr          = SNESMSInitializePackage();CHKERRQ(ierr);
214   ierr          = PetscNew(&link);CHKERRQ(ierr);
215   t             = &link->tab;
216   ierr          = PetscStrallocpy(name,&t->name);CHKERRQ(ierr);
217   t->nstages    = nstages;
218   t->nregisters = nregisters;
219   t->stability  = stability;
220 
221   ierr = PetscMalloc3(nstages*nregisters,&t->gamma,nstages,&t->delta,nstages,&t->betasub);CHKERRQ(ierr);
222   ierr = PetscMemcpy(t->gamma,gamma,nstages*nregisters*sizeof(PetscReal));CHKERRQ(ierr);
223   ierr = PetscMemcpy(t->delta,delta,nstages*sizeof(PetscReal));CHKERRQ(ierr);
224   ierr = PetscMemcpy(t->betasub,betasub,nstages*sizeof(PetscReal));CHKERRQ(ierr);
225 
226   link->next        = SNESMSTableauList;
227   SNESMSTableauList = link;
228   PetscFunctionReturn(0);
229 }
230 
231 /*
232   X - initial state, updated in-place.
233   F - residual, computed at the initial X on input
234 */
235 static PetscErrorCode SNESMSStep_3Sstar(SNES snes,Vec X,Vec F)
236 {
237   PetscErrorCode  ierr;
238   SNES_MS         *ms    = (SNES_MS*)snes->data;
239   SNESMSTableau   t      = ms->tableau;
240   const PetscReal *gamma = t->gamma,*delta = t->delta,*betasub = t->betasub;
241   Vec             S1,S2,S3,Y;
242   PetscInt        i,nstages = t->nstages;;
243 
244 
245   PetscFunctionBegin;
246   Y    = snes->work[0];
247   S1   = X;
248   S2   = snes->work[1];
249   S3   = snes->work[2];
250   ierr = VecZeroEntries(S2);CHKERRQ(ierr);
251   ierr = VecCopy(X,S3);CHKERRQ(ierr);
252   for (i=0; i<nstages; i++) {
253     Vec         Ss[4];
254     PetscScalar scoeff[4];
255 
256     Ss[0] = S1; Ss[1] = S2; Ss[2] = S3; Ss[3] = Y;
257 
258     scoeff[0] = gamma[0*nstages+i]-(PetscReal)1.0;
259     scoeff[1] = gamma[1*nstages+i];
260     scoeff[2] = gamma[2*nstages+i];
261     scoeff[3] = -betasub[i]*ms->damping;
262 
263     ierr = VecAXPY(S2,delta[i],S1);CHKERRQ(ierr);
264     if (i>0) {
265       ierr = SNESComputeFunction(snes,S1,F);CHKERRQ(ierr);
266     }
267     ierr = KSPSolve(snes->ksp,F,Y);CHKERRQ(ierr);
268     ierr = VecMAXPY(S1,4,scoeff,Ss);CHKERRQ(ierr);
269   }
270   PetscFunctionReturn(0);
271 }
272 
273 static PetscErrorCode SNESSolve_MS(SNES snes)
274 {
275   SNES_MS        *ms = (SNES_MS*)snes->data;
276   Vec            X   = snes->vec_sol,F = snes->vec_func;
277   PetscReal      fnorm;
278   PetscInt       i;
279   PetscErrorCode ierr;
280 
281   PetscFunctionBegin;
282   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);
283 
284   ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr);
285   snes->reason = SNES_CONVERGED_ITERATING;
286   ierr         = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
287   snes->iter   = 0;
288   snes->norm   = 0.;
289   ierr         = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
290   if (!snes->vec_func_init_set) {
291     ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
292   } else snes->vec_func_init_set = PETSC_FALSE;
293 
294   if (snes->jacobian) {         /* This method does not require a Jacobian, but it is usually preconditioned by PBJacobi */
295     ierr = SNESComputeJacobian(snes,snes->vec_sol,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr);
296   }
297   if (ms->norms) {
298     ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); /* fnorm <- ||F||  */
299     SNESCheckFunctionNorm(snes,fnorm);
300     /* Monitor convergence */
301     ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
302     snes->iter = 0;
303     snes->norm = fnorm;
304     ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
305     ierr       = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr);
306     ierr       = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr);
307 
308     /* Test for convergence */
309     ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
310     if (snes->reason) PetscFunctionReturn(0);
311   }
312 
313   /* Call general purpose update function */
314   if (snes->ops->update) {
315     ierr = (*snes->ops->update)(snes,snes->iter);CHKERRQ(ierr);
316   }
317   for (i = 0; i < snes->max_its; i++) {
318     ierr = SNESMSStep_3Sstar(snes,X,F);CHKERRQ(ierr);
319 
320     if (i+1 < snes->max_its || ms->norms) {
321       ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr);
322     }
323 
324     if (ms->norms) {
325       ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); /* fnorm <- ||F||  */
326       SNESCheckFunctionNorm(snes,fnorm);
327 
328       /* Monitor convergence */
329       ierr       = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr);
330       snes->iter = i+1;
331       snes->norm = fnorm;
332       ierr       = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr);
333       ierr       = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr);
334       ierr       = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr);
335 
336       /* Test for convergence */
337       ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr);
338       if (snes->reason) PetscFunctionReturn(0);
339     }
340 
341     /* Call general purpose update function */
342     if (snes->ops->update) {
343       ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr);
344     }
345   }
346   if (!snes->reason) snes->reason = SNES_CONVERGED_ITS;
347   PetscFunctionReturn(0);
348 }
349 
350 static PetscErrorCode SNESSetUp_MS(SNES snes)
351 {
352   SNES_MS        *ms = (SNES_MS*)snes->data;
353   PetscErrorCode ierr;
354 
355   PetscFunctionBegin;
356   if (!ms->tableau) {ierr = SNESMSSetType(snes,SNESMSDefault);CHKERRQ(ierr);}
357   ierr = SNESSetWorkVecs(snes,3);CHKERRQ(ierr);
358   ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr);
359   PetscFunctionReturn(0);
360 }
361 
362 static PetscErrorCode SNESReset_MS(SNES snes)
363 {
364 
365   PetscFunctionBegin;
366   PetscFunctionReturn(0);
367 }
368 
369 static PetscErrorCode SNESDestroy_MS(SNES snes)
370 {
371   PetscErrorCode ierr;
372 
373   PetscFunctionBegin;
374   ierr = PetscFree(snes->data);CHKERRQ(ierr);
375   ierr = PetscObjectComposeFunction((PetscObject)snes,"",NULL);CHKERRQ(ierr);
376   PetscFunctionReturn(0);
377 }
378 
379 static PetscErrorCode SNESView_MS(SNES snes,PetscViewer viewer)
380 {
381   PetscBool      iascii;
382   PetscErrorCode ierr;
383   SNES_MS        *ms = (SNES_MS*)snes->data;
384 
385   PetscFunctionBegin;
386   ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr);
387   if (iascii) {
388     SNESMSTableau tab = ms->tableau;
389     ierr = PetscViewerASCIIPrintf(viewer,"  multi-stage method type: %s\n",tab ? tab->name : "not yet set");CHKERRQ(ierr);
390   }
391   PetscFunctionReturn(0);
392 }
393 
394 static PetscErrorCode SNESSetFromOptions_MS(PetscOptionItems *PetscOptionsObject,SNES snes)
395 {
396   SNES_MS        *ms = (SNES_MS*)snes->data;
397   PetscErrorCode ierr;
398 
399   PetscFunctionBegin;
400   ierr = PetscOptionsHead(PetscOptionsObject,"SNES MS options");CHKERRQ(ierr);
401   {
402     SNESMSTableauLink link;
403     PetscInt          count,choice;
404     PetscBool         flg;
405     const char        **namelist;
406     char              mstype[256];
407 
408     ierr = PetscStrncpy(mstype,SNESMSDefault,sizeof(mstype));CHKERRQ(ierr);
409     for (link=SNESMSTableauList,count=0; link; link=link->next,count++) ;
410     ierr = PetscMalloc1(count,(char***)&namelist);CHKERRQ(ierr);
411     for (link=SNESMSTableauList,count=0; link; link=link->next,count++) namelist[count] = link->tab.name;
412     ierr = PetscOptionsEList("-snes_ms_type","Multistage smoother type","SNESMSSetType",(const char*const*)namelist,count,mstype,&choice,&flg);CHKERRQ(ierr);
413     ierr = SNESMSSetType(snes,flg ? namelist[choice] : mstype);CHKERRQ(ierr);
414     ierr = PetscFree(namelist);CHKERRQ(ierr);
415     ierr = PetscOptionsReal("-snes_ms_damping","Damping for multistage method","SNESMSSetDamping",ms->damping,&ms->damping,NULL);CHKERRQ(ierr);
416     ierr = PetscOptionsBool("-snes_ms_norms","Compute norms for monitoring","none",ms->norms,&ms->norms,NULL);CHKERRQ(ierr);
417   }
418   ierr = PetscOptionsTail();CHKERRQ(ierr);
419   PetscFunctionReturn(0);
420 }
421 
422 PetscErrorCode  SNESMSSetType_MS(SNES snes,SNESMSType mstype)
423 {
424   PetscErrorCode    ierr;
425   SNES_MS           *ms = (SNES_MS*)snes->data;
426   SNESMSTableauLink link;
427   PetscBool         match;
428 
429   PetscFunctionBegin;
430   if (ms->tableau) {
431     ierr = PetscStrcmp(ms->tableau->name,mstype,&match);CHKERRQ(ierr);
432     if (match) PetscFunctionReturn(0);
433   }
434   for (link = SNESMSTableauList; link; link=link->next) {
435     ierr = PetscStrcmp(link->tab.name,mstype,&match);CHKERRQ(ierr);
436     if (match) {
437       ierr        = SNESReset_MS(snes);CHKERRQ(ierr);
438       ms->tableau = &link->tab;
439       PetscFunctionReturn(0);
440     }
441   }
442   SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_UNKNOWN_TYPE,"Could not find '%s'",mstype);
443   PetscFunctionReturn(0);
444 }
445 
446 /*@C
447   SNESMSSetType - Set the type of multistage smoother
448 
449   Logically collective
450 
451   Input Parameter:
452 +  snes - nonlinear solver context
453 -  mstype - type of multistage method
454 
455   Level: beginner
456 
457 .seealso: SNESMSGetType(), SNESMS
458 @*/
459 PetscErrorCode SNESMSSetType(SNES snes,SNESMSType rostype)
460 {
461   PetscErrorCode ierr;
462 
463   PetscFunctionBegin;
464   PetscValidHeaderSpecific(snes,SNES_CLASSID,1);
465   ierr = PetscTryMethod(snes,"SNESMSSetType_C",(SNES,SNESMSType),(snes,rostype));CHKERRQ(ierr);
466   PetscFunctionReturn(0);
467 }
468 
469 /* -------------------------------------------------------------------------- */
470 /*MC
471       SNESMS - multi-stage smoothers
472 
473       Options Database:
474 
475 +     -snes_ms_type - type of multi-stage smoother
476 -     -snes_ms_damping - damping for multi-stage method
477 
478       Notes:
479       These multistage methods are explicit Runge-Kutta methods that are often used as smoothers for
480       FAS multigrid for transport problems. In the linear case, these are equivalent to polynomial smoothers (such as Chebyshev).
481 
482       Multi-stage smoothers should usually be preconditioned by point-block Jacobi to ensure proper scaling and to normalize the wave speeds.
483 
484       The methods are specified in low storage form (Ketcheson 2010). New methods can be registered with SNESMSRegister().
485 
486       References:
487 +   1. -   Ketcheson (2010) Runge Kutta methods with minimum storage implementations.
488 .   2. -   Jameson (1983) Solution of the Euler equations for two dimensional transonic flow by a multigrid method.
489 -   3. -   Pierce and Giles (1997) Preconditioned multigrid methods for compressible flow calculations on stretched meshes.
490 
491       Level: beginner
492 
493 .seealso:  SNESCreate(), SNES, SNESSetType(), SNESMS, SNESFAS, KSPCHEBYSHEV
494 
495 M*/
496 PETSC_EXTERN PetscErrorCode SNESCreate_MS(SNES snes)
497 {
498   PetscErrorCode ierr;
499   SNES_MS        *ms;
500 
501   PetscFunctionBegin;
502   ierr = SNESMSInitializePackage();CHKERRQ(ierr);
503 
504   snes->ops->setup          = SNESSetUp_MS;
505   snes->ops->solve          = SNESSolve_MS;
506   snes->ops->destroy        = SNESDestroy_MS;
507   snes->ops->setfromoptions = SNESSetFromOptions_MS;
508   snes->ops->view           = SNESView_MS;
509   snes->ops->reset          = SNESReset_MS;
510 
511   snes->usesnpc = PETSC_FALSE;
512   snes->usesksp = PETSC_TRUE;
513 
514   snes->alwayscomputesfinalresidual = PETSC_FALSE;
515 
516   ierr        = PetscNewLog(snes,&ms);CHKERRQ(ierr);
517   snes->data  = (void*)ms;
518   ms->damping = 0.9;
519   ms->norms   = PETSC_FALSE;
520 
521   ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",SNESMSSetType_MS);CHKERRQ(ierr);
522   PetscFunctionReturn(0);
523 }
524 
525