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