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