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 PetscValidRealPointer(gamma,4); 202 PetscValidRealPointer(delta,5); 203 PetscValidRealPointer(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] - 1; 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 SNESMSStep_Norms(SNES snes,PetscInt iter,Vec F) 266 { 267 SNES_MS *ms = (SNES_MS*)snes->data; 268 PetscReal fnorm; 269 PetscErrorCode ierr; 270 271 PetscFunctionBegin; 272 if (ms->norms) { 273 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 274 SNESCheckFunctionNorm(snes,fnorm); 275 /* Monitor convergence */ 276 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 277 snes->iter = iter; 278 snes->norm = fnorm; 279 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 280 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 281 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 282 /* Test for convergence */ 283 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 284 } else if (iter > 0) { 285 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 286 snes->iter = iter; 287 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 288 } 289 PetscFunctionReturn(0); 290 } 291 292 293 static PetscErrorCode SNESSolve_MS(SNES snes) 294 { 295 SNES_MS *ms = (SNES_MS*)snes->data; 296 Vec X = snes->vec_sol,F = snes->vec_func; 297 PetscInt i; 298 PetscErrorCode ierr; 299 300 PetscFunctionBegin; 301 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); 302 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 303 304 snes->reason = SNES_CONVERGED_ITERATING; 305 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 306 snes->iter = 0; 307 snes->norm = 0; 308 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 309 310 if (!snes->vec_func_init_set) { 311 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 312 } else snes->vec_func_init_set = PETSC_FALSE; 313 314 ierr = SNESMSStep_Norms(snes,0,F);CHKERRQ(ierr); 315 if (snes->reason) PetscFunctionReturn(0); 316 317 for (i = 0; i < snes->max_its; i++) { 318 319 /* Call general purpose update function */ 320 if (snes->ops->update) { 321 ierr = (*snes->ops->update)(snes,snes->iter);CHKERRQ(ierr); 322 } 323 324 if (i == 0 && snes->jacobian) { 325 /* This method does not require a Jacobian, but it is usually preconditioned by PBJacobi */ 326 ierr = SNESComputeJacobian(snes,snes->vec_sol,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 327 SNESCheckJacobianDomainerror(snes); 328 ierr = KSPSetOperators(snes->ksp,snes->jacobian,snes->jacobian_pre);CHKERRQ(ierr); 329 } 330 331 ierr = SNESMSStep_3Sstar(snes,X,F);CHKERRQ(ierr); 332 333 if (i < snes->max_its-1 || ms->norms) { 334 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 335 } 336 337 ierr = SNESMSStep_Norms(snes,i+1,F);CHKERRQ(ierr); 338 if (snes->reason) PetscFunctionReturn(0); 339 } 340 if (!snes->reason) snes->reason = SNES_CONVERGED_ITS; 341 PetscFunctionReturn(0); 342 } 343 344 static PetscErrorCode SNESSetUp_MS(SNES snes) 345 { 346 SNES_MS *ms = (SNES_MS*)snes->data; 347 PetscErrorCode ierr; 348 349 PetscFunctionBegin; 350 if (!ms->tableau) {ierr = SNESMSSetType(snes,SNESMSDefault);CHKERRQ(ierr);} 351 ierr = SNESSetWorkVecs(snes,3);CHKERRQ(ierr); 352 ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 353 PetscFunctionReturn(0); 354 } 355 356 static PetscErrorCode SNESReset_MS(SNES snes) 357 { 358 359 PetscFunctionBegin; 360 PetscFunctionReturn(0); 361 } 362 363 static PetscErrorCode SNESDestroy_MS(SNES snes) 364 { 365 PetscErrorCode ierr; 366 367 PetscFunctionBegin; 368 ierr = SNESReset_MS(snes);CHKERRQ(ierr); 369 ierr = PetscFree(snes->data);CHKERRQ(ierr); 370 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",NULL);CHKERRQ(ierr); 371 PetscFunctionReturn(0); 372 } 373 374 static PetscErrorCode SNESView_MS(SNES snes,PetscViewer viewer) 375 { 376 PetscBool iascii; 377 PetscErrorCode ierr; 378 SNES_MS *ms = (SNES_MS*)snes->data; 379 380 PetscFunctionBegin; 381 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 382 if (iascii) { 383 SNESMSTableau tab = ms->tableau; 384 ierr = PetscViewerASCIIPrintf(viewer," multi-stage method type: %s\n",tab ? tab->name : "not yet set");CHKERRQ(ierr); 385 } 386 PetscFunctionReturn(0); 387 } 388 389 static PetscErrorCode SNESSetFromOptions_MS(PetscOptionItems *PetscOptionsObject,SNES snes) 390 { 391 SNES_MS *ms = (SNES_MS*)snes->data; 392 PetscErrorCode ierr; 393 394 PetscFunctionBegin; 395 ierr = PetscOptionsHead(PetscOptionsObject,"SNES MS options");CHKERRQ(ierr); 396 { 397 SNESMSTableauLink link; 398 PetscInt count,choice; 399 PetscBool flg; 400 const char **namelist; 401 char mstype[256]; 402 403 ierr = PetscStrncpy(mstype,SNESMSDefault,sizeof(mstype));CHKERRQ(ierr); 404 for (link=SNESMSTableauList,count=0; link; link=link->next,count++) ; 405 ierr = PetscMalloc1(count,(char***)&namelist);CHKERRQ(ierr); 406 for (link=SNESMSTableauList,count=0; link; link=link->next,count++) namelist[count] = link->tab.name; 407 ierr = PetscOptionsEList("-snes_ms_type","Multistage smoother type","SNESMSSetType",(const char*const*)namelist,count,mstype,&choice,&flg);CHKERRQ(ierr); 408 ierr = SNESMSSetType(snes,flg ? namelist[choice] : mstype);CHKERRQ(ierr); 409 ierr = PetscFree(namelist);CHKERRQ(ierr); 410 ierr = PetscOptionsReal("-snes_ms_damping","Damping for multistage method","SNESMSSetDamping",ms->damping,&ms->damping,NULL);CHKERRQ(ierr); 411 ierr = PetscOptionsBool("-snes_ms_norms","Compute norms for monitoring","none",ms->norms,&ms->norms,NULL);CHKERRQ(ierr); 412 } 413 ierr = PetscOptionsTail();CHKERRQ(ierr); 414 PetscFunctionReturn(0); 415 } 416 417 PetscErrorCode SNESMSSetType_MS(SNES snes,SNESMSType mstype) 418 { 419 PetscErrorCode ierr; 420 SNES_MS *ms = (SNES_MS*)snes->data; 421 SNESMSTableauLink link; 422 PetscBool match; 423 424 PetscFunctionBegin; 425 if (ms->tableau) { 426 ierr = PetscStrcmp(ms->tableau->name,mstype,&match);CHKERRQ(ierr); 427 if (match) PetscFunctionReturn(0); 428 } 429 for (link = SNESMSTableauList; link; link=link->next) { 430 ierr = PetscStrcmp(link->tab.name,mstype,&match);CHKERRQ(ierr); 431 if (match) { 432 ierr = SNESReset_MS(snes);CHKERRQ(ierr); 433 ms->tableau = &link->tab; 434 PetscFunctionReturn(0); 435 } 436 } 437 SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_UNKNOWN_TYPE,"Could not find '%s'",mstype); 438 PetscFunctionReturn(0); 439 } 440 441 /*@C 442 SNESMSSetType - Set the type of multistage smoother 443 444 Logically collective 445 446 Input Parameter: 447 + snes - nonlinear solver context 448 - mstype - type of multistage method 449 450 Level: beginner 451 452 .seealso: SNESMSGetType(), SNESMS 453 @*/ 454 PetscErrorCode SNESMSSetType(SNES snes,SNESMSType rostype) 455 { 456 PetscErrorCode ierr; 457 458 PetscFunctionBegin; 459 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 460 ierr = PetscTryMethod(snes,"SNESMSSetType_C",(SNES,SNESMSType),(snes,rostype));CHKERRQ(ierr); 461 PetscFunctionReturn(0); 462 } 463 464 /* -------------------------------------------------------------------------- */ 465 /*MC 466 SNESMS - multi-stage smoothers 467 468 Options Database: 469 470 + -snes_ms_type - type of multi-stage smoother 471 - -snes_ms_damping - damping for multi-stage method 472 473 Notes: 474 These multistage methods are explicit Runge-Kutta methods that are often used as smoothers for 475 FAS multigrid for transport problems. In the linear case, these are equivalent to polynomial smoothers (such as Chebyshev). 476 477 Multi-stage smoothers should usually be preconditioned by point-block Jacobi to ensure proper scaling and to normalize the wave speeds. 478 479 The methods are specified in low storage form (Ketcheson 2010). New methods can be registered with SNESMSRegister(). 480 481 References: 482 + 1. - Ketcheson (2010) Runge Kutta methods with minimum storage implementations. 483 . 2. - Jameson (1983) Solution of the Euler equations for two dimensional transonic flow by a multigrid method. 484 - 3. - Pierce and Giles (1997) Preconditioned multigrid methods for compressible flow calculations on stretched meshes. 485 486 Level: beginner 487 488 .seealso: SNESCreate(), SNES, SNESSetType(), SNESMS, SNESFAS, KSPCHEBYSHEV 489 490 M*/ 491 PETSC_EXTERN PetscErrorCode SNESCreate_MS(SNES snes) 492 { 493 PetscErrorCode ierr; 494 SNES_MS *ms; 495 496 PetscFunctionBegin; 497 ierr = SNESMSInitializePackage();CHKERRQ(ierr); 498 499 snes->ops->setup = SNESSetUp_MS; 500 snes->ops->solve = SNESSolve_MS; 501 snes->ops->destroy = SNESDestroy_MS; 502 snes->ops->setfromoptions = SNESSetFromOptions_MS; 503 snes->ops->view = SNESView_MS; 504 snes->ops->reset = SNESReset_MS; 505 506 snes->usesnpc = PETSC_FALSE; 507 snes->usesksp = PETSC_TRUE; 508 509 snes->alwayscomputesfinalresidual = PETSC_FALSE; 510 511 ierr = PetscNewLog(snes,&ms);CHKERRQ(ierr); 512 snes->data = (void*)ms; 513 ms->damping = 0.9; 514 ms->norms = PETSC_FALSE; 515 516 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",SNESMSSetType_MS);CHKERRQ(ierr); 517 PetscFunctionReturn(0); 518 } 519 520