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 SNESMSRegister(). 106 107 Not Collective 108 109 Level: advanced 110 111 .seealso: SNESMSRegister(), SNESMSRegisterAll() 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 PetscErrorCode ierr; 347 348 PetscFunctionBegin; 349 ierr = SNESSetWorkVecs(snes,3);CHKERRQ(ierr); 350 ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 351 PetscFunctionReturn(0); 352 } 353 354 static PetscErrorCode SNESReset_MS(SNES snes) 355 { 356 PetscFunctionBegin; 357 PetscFunctionReturn(0); 358 } 359 360 static PetscErrorCode SNESDestroy_MS(SNES snes) 361 { 362 PetscErrorCode ierr; 363 364 PetscFunctionBegin; 365 ierr = SNESReset_MS(snes);CHKERRQ(ierr); 366 ierr = PetscFree(snes->data);CHKERRQ(ierr); 367 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetType_C",NULL);CHKERRQ(ierr); 368 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",NULL);CHKERRQ(ierr); 369 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetDamping_C",NULL);CHKERRQ(ierr); 370 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetDamping_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 SNESMSTableau tab = ms->tableau; 380 381 PetscFunctionBegin; 382 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 383 if (iascii) { 384 ierr = PetscViewerASCIIPrintf(viewer," multi-stage method type: %s\n",tab->name);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 SNESMSType mstype; 402 PetscReal damping; 403 404 ierr = SNESMSGetType(snes,&mstype);CHKERRQ(ierr); 405 for (link=SNESMSTableauList,count=0; link; link=link->next,count++) ; 406 ierr = PetscMalloc1(count,(char***)&namelist);CHKERRQ(ierr); 407 for (link=SNESMSTableauList,count=0; link; link=link->next,count++) namelist[count] = link->tab.name; 408 ierr = PetscOptionsEList("-snes_ms_type","Multistage smoother type","SNESMSSetType",(const char*const*)namelist,count,mstype,&choice,&flg);CHKERRQ(ierr); 409 if (flg) {ierr = SNESMSSetType(snes,namelist[choice]);CHKERRQ(ierr);} 410 ierr = PetscFree(namelist);CHKERRQ(ierr); 411 ierr = SNESMSGetDamping(snes,&damping);CHKERRQ(ierr); 412 ierr = PetscOptionsReal("-snes_ms_damping","Damping for multistage method","SNESMSSetDamping",damping,&damping,&flg);CHKERRQ(ierr); 413 if (flg) {ierr = SNESMSSetDamping(snes,damping);CHKERRQ(ierr);} 414 ierr = PetscOptionsBool("-snes_ms_norms","Compute norms for monitoring","none",ms->norms,&ms->norms,NULL);CHKERRQ(ierr); 415 } 416 ierr = PetscOptionsTail();CHKERRQ(ierr); 417 PetscFunctionReturn(0); 418 } 419 420 static PetscErrorCode SNESMSGetType_MS(SNES snes,SNESMSType *mstype) 421 { 422 SNES_MS *ms = (SNES_MS*)snes->data; 423 SNESMSTableau tab = ms->tableau; 424 425 PetscFunctionBegin; 426 *mstype = tab->name; 427 PetscFunctionReturn(0); 428 } 429 430 static PetscErrorCode SNESMSSetType_MS(SNES snes,SNESMSType mstype) 431 { 432 SNES_MS *ms = (SNES_MS*)snes->data; 433 SNESMSTableauLink link; 434 PetscBool match; 435 PetscErrorCode ierr; 436 437 PetscFunctionBegin; 438 if (ms->tableau) { 439 ierr = PetscStrcmp(ms->tableau->name,mstype,&match);CHKERRQ(ierr); 440 if (match) PetscFunctionReturn(0); 441 } 442 for (link = SNESMSTableauList; link; link=link->next) { 443 ierr = PetscStrcmp(link->tab.name,mstype,&match);CHKERRQ(ierr); 444 if (match) { 445 if (snes->setupcalled) {ierr = SNESReset_MS(snes);CHKERRQ(ierr);} 446 ms->tableau = &link->tab; 447 if (snes->setupcalled) {ierr = SNESSetUp_MS(snes);CHKERRQ(ierr);} 448 PetscFunctionReturn(0); 449 } 450 } 451 SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_UNKNOWN_TYPE,"Could not find '%s'",mstype); 452 PetscFunctionReturn(0); 453 } 454 455 /*@C 456 SNESMSGetType - Get the type of multistage smoother 457 458 Not collective 459 460 Input Parameter: 461 . snes - nonlinear solver context 462 463 Output Parameter: 464 . mstype - type of multistage method 465 466 Level: beginner 467 468 .seealso: SNESMSSetType(), SNESMSType, SNESMS 469 @*/ 470 PetscErrorCode SNESMSGetType(SNES snes,SNESMSType *mstype) 471 { 472 PetscErrorCode ierr; 473 474 PetscFunctionBegin; 475 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 476 PetscValidPointer(mstype,2); 477 ierr = PetscUseMethod(snes,"SNESMSGetType_C",(SNES,SNESMSType*),(snes,mstype));CHKERRQ(ierr); 478 PetscFunctionReturn(0); 479 } 480 481 /*@C 482 SNESMSSetType - Set the type of multistage smoother 483 484 Logically collective 485 486 Input Parameter: 487 + snes - nonlinear solver context 488 - mstype - type of multistage method 489 490 Level: beginner 491 492 .seealso: SNESMSGetType(), SNESMSType, SNESMS 493 @*/ 494 PetscErrorCode SNESMSSetType(SNES snes,SNESMSType mstype) 495 { 496 PetscErrorCode ierr; 497 498 PetscFunctionBegin; 499 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 500 PetscValidPointer(mstype,2); 501 ierr = PetscTryMethod(snes,"SNESMSSetType_C",(SNES,SNESMSType),(snes,mstype));CHKERRQ(ierr); 502 PetscFunctionReturn(0); 503 } 504 505 static PetscErrorCode SNESMSGetDamping_MS(SNES snes,PetscReal *damping) 506 { 507 SNES_MS *ms = (SNES_MS*)snes->data; 508 509 PetscFunctionBegin; 510 *damping = ms->damping; 511 PetscFunctionReturn(0); 512 } 513 514 static PetscErrorCode SNESMSSetDamping_MS(SNES snes,PetscReal damping) 515 { 516 SNES_MS *ms = (SNES_MS*)snes->data; 517 518 PetscFunctionBegin; 519 ms->damping = damping; 520 PetscFunctionReturn(0); 521 } 522 523 /*@C 524 SNESMSGetDamping - Get the damping parameter 525 526 Not collective 527 528 Input Parameter: 529 . snes - nonlinear solver context 530 531 Output Parameter: 532 . damping - damping parameter 533 534 Level: advanced 535 536 .seealso: SNESMSSetDamping(), SNESMS 537 @*/ 538 PetscErrorCode SNESMSGetDamping(SNES snes,PetscReal *damping) 539 { 540 PetscErrorCode ierr; 541 542 PetscFunctionBegin; 543 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 544 PetscValidPointer(damping,2); 545 ierr = PetscUseMethod(snes,"SNESMSGetDamping_C",(SNES,PetscReal*),(snes,damping));CHKERRQ(ierr); 546 PetscFunctionReturn(0); 547 } 548 549 /*@C 550 SNESMSSetDamping - Set the damping parameter 551 552 Logically collective 553 554 Input Parameter: 555 + snes - nonlinear solver context 556 - damping - damping parameter 557 558 Level: advanced 559 560 .seealso: SNESMSGetDamping(), SNESMS 561 @*/ 562 PetscErrorCode SNESMSSetDamping(SNES snes,PetscReal damping) 563 { 564 PetscErrorCode ierr; 565 566 PetscFunctionBegin; 567 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 568 PetscValidLogicalCollectiveReal(snes,damping,2); 569 ierr = PetscTryMethod(snes,"SNESMSSetDamping_C",(SNES,PetscReal),(snes,damping));CHKERRQ(ierr); 570 PetscFunctionReturn(0); 571 } 572 573 /* -------------------------------------------------------------------------- */ 574 /*MC 575 SNESMS - multi-stage smoothers 576 577 Options Database: 578 579 + -snes_ms_type - type of multi-stage smoother 580 - -snes_ms_damping - damping for multi-stage method 581 582 Notes: 583 These multistage methods are explicit Runge-Kutta methods that are often used as smoothers for 584 FAS multigrid for transport problems. In the linear case, these are equivalent to polynomial smoothers (such as Chebyshev). 585 586 Multi-stage smoothers should usually be preconditioned by point-block Jacobi to ensure proper scaling and to normalize the wave speeds. 587 588 The methods are specified in low storage form (Ketcheson 2010). New methods can be registered with SNESMSRegister(). 589 590 References: 591 + 1. - Ketcheson (2010) Runge Kutta methods with minimum storage implementations (https://doi.org/10.1016/j.jcp.2009.11.006). 592 . 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). 593 . 3. - Pierce and Giles (1997) Preconditioned multigrid methods for compressible flow calculations on stretched meshes (https://doi.org/10.1006/jcph.1997.5772). 594 - 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). 595 596 Level: beginner 597 598 .seealso: SNESCreate(), SNES, SNESSetType(), SNESMS, SNESFAS, KSPCHEBYSHEV 599 600 M*/ 601 PETSC_EXTERN PetscErrorCode SNESCreate_MS(SNES snes) 602 { 603 PetscErrorCode ierr; 604 SNES_MS *ms; 605 606 PetscFunctionBegin; 607 ierr = SNESMSInitializePackage();CHKERRQ(ierr); 608 609 snes->ops->setup = SNESSetUp_MS; 610 snes->ops->solve = SNESSolve_MS; 611 snes->ops->destroy = SNESDestroy_MS; 612 snes->ops->setfromoptions = SNESSetFromOptions_MS; 613 snes->ops->view = SNESView_MS; 614 snes->ops->reset = SNESReset_MS; 615 616 snes->usesnpc = PETSC_FALSE; 617 snes->usesksp = PETSC_TRUE; 618 619 snes->alwayscomputesfinalresidual = PETSC_FALSE; 620 621 ierr = PetscNewLog(snes,&ms);CHKERRQ(ierr); 622 snes->data = (void*)ms; 623 ms->damping = 0.9; 624 ms->norms = PETSC_FALSE; 625 626 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetType_C",SNESMSGetType_MS);CHKERRQ(ierr); 627 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetType_C",SNESMSSetType_MS);CHKERRQ(ierr); 628 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSGetDamping_C",SNESMSGetDamping_MS);CHKERRQ(ierr); 629 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESMSSetDamping_C",SNESMSSetDamping_MS);CHKERRQ(ierr); 630 631 ierr = SNESMSSetType(snes,SNESMSDefault);CHKERRQ(ierr); 632 PetscFunctionReturn(0); 633 } 634 635