1 #include <../src/snes/impls/gs/gsimpl.h> /*I "petscsnes.h" I*/ 2 3 /*@ 4 SNESNGSSetTolerances - Sets various parameters used in convergence tests. 5 6 Logically Collective on SNES 7 8 Input Parameters: 9 + snes - the SNES context 10 . abstol - absolute convergence tolerance 11 . rtol - relative convergence tolerance 12 . stol - convergence tolerance in terms of the norm of the change in the solution between steps, || delta x || < stol*|| x || 13 - maxit - maximum number of iterations 14 15 16 Options Database Keys: 17 + -snes_ngs_atol <abstol> - Sets abstol 18 . -snes_ngs_rtol <rtol> - Sets rtol 19 . -snes_ngs_stol <stol> - Sets stol 20 - -snes_max_it <maxit> - Sets maxit 21 22 Level: intermediate 23 24 .seealso: SNESSetTrustRegionTolerance() 25 @*/ 26 PetscErrorCode SNESNGSSetTolerances(SNES snes,PetscReal abstol,PetscReal rtol,PetscReal stol,PetscInt maxit) 27 { 28 SNES_NGS *gs = (SNES_NGS*)snes->data; 29 30 PetscFunctionBegin; 31 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 32 33 if (abstol != PETSC_DEFAULT) { 34 if (abstol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Absolute tolerance %g must be non-negative",(double)abstol); 35 gs->abstol = abstol; 36 } 37 if (rtol != PETSC_DEFAULT) { 38 if (rtol < 0.0 || 1.0 <= rtol) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Relative tolerance %g must be non-negative and less than 1.0",(double)rtol); 39 gs->rtol = rtol; 40 } 41 if (stol != PETSC_DEFAULT) { 42 if (stol < 0.0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Step tolerance %g must be non-negative",(double)stol); 43 gs->stol = stol; 44 } 45 if (maxit != PETSC_DEFAULT) { 46 if (maxit < 0) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_OUTOFRANGE,"Maximum number of iterations %D must be non-negative",maxit); 47 gs->max_its = maxit; 48 } 49 PetscFunctionReturn(0); 50 } 51 52 /*@ 53 SNESNGSGetTolerances - Gets various parameters used in convergence tests. 54 55 Not Collective 56 57 Input Parameters: 58 + snes - the SNES context 59 . atol - absolute convergence tolerance 60 . rtol - relative convergence tolerance 61 . stol - convergence tolerance in terms of the norm 62 of the change in the solution between steps 63 - maxit - maximum number of iterations 64 65 Notes: 66 The user can specify NULL for any parameter that is not needed. 67 68 Level: intermediate 69 70 .seealso: SNESSetTolerances() 71 @*/ 72 PetscErrorCode SNESNGSGetTolerances(SNES snes,PetscReal *atol,PetscReal *rtol,PetscReal *stol,PetscInt *maxit) 73 { 74 SNES_NGS *gs = (SNES_NGS*)snes->data; 75 76 PetscFunctionBegin; 77 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 78 if (atol) *atol = gs->abstol; 79 if (rtol) *rtol = gs->rtol; 80 if (stol) *stol = gs->stol; 81 if (maxit) *maxit = gs->max_its; 82 PetscFunctionReturn(0); 83 } 84 85 /*@ 86 SNESNGSSetSweeps - Sets the number of sweeps of GS to use. 87 88 Input Parameters: 89 + snes - the SNES context 90 - sweeps - the number of sweeps of GS to perform. 91 92 Level: intermediate 93 94 .seealso: SNESSetNGS(), SNESGetNGS(), SNESSetNPC(), SNESNGSGetSweeps() 95 @*/ 96 97 PetscErrorCode SNESNGSSetSweeps(SNES snes, PetscInt sweeps) 98 { 99 SNES_NGS *gs = (SNES_NGS*)snes->data; 100 101 PetscFunctionBegin; 102 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 103 gs->sweeps = sweeps; 104 PetscFunctionReturn(0); 105 } 106 107 /*@ 108 SNESNGSGetSweeps - Gets the number of sweeps GS will use. 109 110 Input Parameters: 111 . snes - the SNES context 112 113 Output Parameters: 114 . sweeps - the number of sweeps of GS to perform. 115 116 Level: intermediate 117 118 .seealso: SNESSetNGS(), SNESGetNGS(), SNESSetNPC(), SNESNGSSetSweeps() 119 @*/ 120 PetscErrorCode SNESNGSGetSweeps(SNES snes, PetscInt * sweeps) 121 { 122 SNES_NGS *gs = (SNES_NGS*)snes->data; 123 124 PetscFunctionBegin; 125 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 126 *sweeps = gs->sweeps; 127 PetscFunctionReturn(0); 128 } 129 130 PetscErrorCode SNESReset_NGS(SNES snes) 131 { 132 SNES_NGS *gs = (SNES_NGS*)snes->data; 133 PetscErrorCode ierr; 134 135 PetscFunctionBegin; 136 ierr = ISColoringDestroy(&gs->coloring);CHKERRQ(ierr); 137 PetscFunctionReturn(0); 138 } 139 140 PetscErrorCode SNESDestroy_NGS(SNES snes) 141 { 142 PetscErrorCode ierr; 143 144 PetscFunctionBegin; 145 ierr = SNESReset_NGS(snes);CHKERRQ(ierr); 146 ierr = PetscFree(snes->data);CHKERRQ(ierr); 147 PetscFunctionReturn(0); 148 } 149 150 PetscErrorCode SNESSetUp_NGS(SNES snes) 151 { 152 PetscErrorCode ierr; 153 PetscErrorCode (*f)(SNES,Vec,Vec,void*); 154 155 PetscFunctionBegin; 156 ierr = SNESGetNGS(snes,&f,NULL);CHKERRQ(ierr); 157 if (!f) { 158 ierr = SNESSetNGS(snes,SNESComputeNGSDefaultSecant,NULL);CHKERRQ(ierr); 159 } 160 PetscFunctionReturn(0); 161 } 162 163 PetscErrorCode SNESSetFromOptions_NGS(PetscOptionItems *PetscOptionsObject,SNES snes) 164 { 165 SNES_NGS *gs = (SNES_NGS*)snes->data; 166 PetscErrorCode ierr; 167 PetscInt sweeps,max_its=PETSC_DEFAULT; 168 PetscReal rtol=PETSC_DEFAULT,atol=PETSC_DEFAULT,stol=PETSC_DEFAULT; 169 PetscBool flg,flg1,flg2,flg3; 170 171 PetscFunctionBegin; 172 ierr = PetscOptionsHead(PetscOptionsObject,"SNES GS options");CHKERRQ(ierr); 173 /* GS Options */ 174 ierr = PetscOptionsInt("-snes_ngs_sweeps","Number of sweeps of GS to apply","SNESComputeGS",gs->sweeps,&sweeps,&flg);CHKERRQ(ierr); 175 if (flg) { 176 ierr = SNESNGSSetSweeps(snes,sweeps);CHKERRQ(ierr); 177 } 178 ierr = PetscOptionsReal("-snes_ngs_atol","Absolute residual tolerance for GS iteration","SNESComputeGS",gs->abstol,&atol,&flg);CHKERRQ(ierr); 179 ierr = PetscOptionsReal("-snes_ngs_rtol","Relative residual tolerance for GS iteration","SNESComputeGS",gs->rtol,&rtol,&flg1);CHKERRQ(ierr); 180 ierr = PetscOptionsReal("-snes_ngs_stol","Absolute update tolerance for GS iteration","SNESComputeGS",gs->stol,&stol,&flg2);CHKERRQ(ierr); 181 ierr = PetscOptionsInt("-snes_ngs_max_it","Maximum number of sweeps of GS to apply","SNESComputeGS",gs->max_its,&max_its,&flg3);CHKERRQ(ierr); 182 if (flg || flg1 || flg2 || flg3) { 183 ierr = SNESNGSSetTolerances(snes,atol,rtol,stol,max_its);CHKERRQ(ierr); 184 } 185 flg = PETSC_FALSE; 186 ierr = PetscOptionsBool("-snes_ngs_secant","Use finite difference secant approximation with coloring","",flg,&flg,NULL);CHKERRQ(ierr); 187 if (flg) { 188 ierr = SNESSetNGS(snes,SNESComputeNGSDefaultSecant,NULL);CHKERRQ(ierr); 189 ierr = PetscInfo(snes,"Setting default finite difference secant approximation with coloring\n");CHKERRQ(ierr); 190 } 191 ierr = PetscOptionsReal("-snes_ngs_secant_h","Differencing parameter for secant search","",gs->h,&gs->h,NULL);CHKERRQ(ierr); 192 ierr = PetscOptionsBool("-snes_ngs_secant_mat_coloring","Use the graph coloring of the Jacobian for the secant GS","",gs->secant_mat,&gs->secant_mat,&flg);CHKERRQ(ierr); 193 194 ierr = PetscOptionsTail();CHKERRQ(ierr); 195 PetscFunctionReturn(0); 196 } 197 198 PetscErrorCode SNESView_NGS(SNES snes, PetscViewer viewer) 199 { 200 PetscErrorCode ierr; 201 PetscErrorCode (*f)(SNES,Vec,Vec,void*); 202 SNES_NGS *gs = (SNES_NGS*)snes->data; 203 PetscBool iascii; 204 205 PetscFunctionBegin; 206 ierr = PetscObjectTypeCompare((PetscObject)viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 207 if (iascii) { 208 ierr = DMSNESGetNGS(snes->dm,&f,NULL);CHKERRQ(ierr); 209 if (f == SNESComputeNGSDefaultSecant) { 210 ierr = PetscViewerASCIIPrintf(viewer," Use finite difference secant approximation with coloring with h = %g \n",(double)gs->h);CHKERRQ(ierr); 211 } 212 } 213 PetscFunctionReturn(0); 214 } 215 216 PetscErrorCode SNESSolve_NGS(SNES snes) 217 { 218 Vec F; 219 Vec X; 220 Vec B; 221 PetscInt i; 222 PetscReal fnorm; 223 PetscErrorCode ierr; 224 SNESNormSchedule normschedule; 225 226 PetscFunctionBegin; 227 228 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); 229 230 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 231 X = snes->vec_sol; 232 F = snes->vec_func; 233 B = snes->vec_rhs; 234 235 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 236 snes->iter = 0; 237 snes->norm = 0.; 238 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 239 snes->reason = SNES_CONVERGED_ITERATING; 240 241 ierr = SNESGetNormSchedule(snes, &normschedule);CHKERRQ(ierr); 242 if (normschedule == SNES_NORM_ALWAYS || normschedule == SNES_NORM_INITIAL_ONLY || normschedule == SNES_NORM_INITIAL_FINAL_ONLY) { 243 /* compute the initial function and preconditioned update delX */ 244 if (!snes->vec_func_init_set) { 245 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 246 } else snes->vec_func_init_set = PETSC_FALSE; 247 248 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 249 SNESCheckFunctionNorm(snes,fnorm); 250 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 251 snes->iter = 0; 252 snes->norm = fnorm; 253 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 254 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 255 ierr = SNESMonitor(snes,0,snes->norm);CHKERRQ(ierr); 256 257 /* test convergence */ 258 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 259 if (snes->reason) PetscFunctionReturn(0); 260 } else { 261 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 262 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 263 } 264 265 /* Call general purpose update function */ 266 if (snes->ops->update) { 267 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 268 } 269 270 for (i = 0; i < snes->max_its; i++) { 271 ierr = SNESComputeNGS(snes, B, X);CHKERRQ(ierr); 272 /* only compute norms if requested or about to exit due to maximum iterations */ 273 if (normschedule == SNES_NORM_ALWAYS || ((i == snes->max_its - 1) && (normschedule == SNES_NORM_INITIAL_FINAL_ONLY || normschedule == SNES_NORM_FINAL_ONLY))) { 274 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 275 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); /* fnorm <- ||F|| */ 276 SNESCheckFunctionNorm(snes,fnorm); 277 /* Monitor convergence */ 278 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 279 snes->iter = i+1; 280 snes->norm = fnorm; 281 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 282 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 283 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 284 } 285 /* Test for convergence */ 286 if (normschedule == SNES_NORM_ALWAYS) ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 287 if (snes->reason) PetscFunctionReturn(0); 288 /* Call general purpose update function */ 289 if (snes->ops->update) { 290 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 291 } 292 } 293 if (normschedule == SNES_NORM_ALWAYS) { 294 if (i == snes->max_its) { 295 ierr = PetscInfo1(snes,"Maximum number of iterations has been reached: %D\n",snes->max_its);CHKERRQ(ierr); 296 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 297 } 298 } else if (!snes->reason) snes->reason = SNES_CONVERGED_ITS; /* GS is meant to be used as a preconditioner */ 299 PetscFunctionReturn(0); 300 } 301 302 /*MC 303 SNESNGS - Either calls the user-provided solution routine provided with SNESSetNGS() or does a finite difference secant approximation 304 using coloring. 305 306 Level: advanced 307 308 Options Database: 309 + -snes_ngs_sweeps <n> - Number of sweeps of GS to apply 310 . -snes_ngs_atol <atol> - Absolute residual tolerance for GS iteration 311 . -snes_ngs_rtol <rtol> - Relative residual tolerance for GS iteration 312 . -snes_ngs_stol <stol> - Absolute update tolerance for GS iteration 313 . -snes_ngs_max_it <maxit> - Maximum number of sweeps of GS to apply 314 . -snes_ngs_secant - Use pointwise secant local Jacobian approximation with coloring instead of user provided Gauss-Seidel routine, this is 315 used by default if no user provided Gauss-Seidel routine is available. Requires either that a DM that can compute a coloring 316 is available or a Jacobian sparse matrix is provided (from which to get the coloring). 317 . -snes_ngs_secant_h <h> - Differencing parameter for secant approximation 318 . -snes_ngs_secant_mat_coloring - Use the graph coloring of the Jacobian for the secant GS even if a DM is available. 319 - -snes_norm_schedule <none, always, initialonly, finalonly, initalfinalonly> - how often the residual norms are computed 320 321 Notes: 322 the Gauss-Seidel smoother is inherited through composition. If a solver has been created with SNESGetNPC(), it will have 323 its parent's Gauss-Seidel routine associated with it. 324 325 By default this routine computes the solution norm at each iteration, this can be time consuming, you can turn this off with SNESSetNormSchedule() or -snes_norm_schedule 326 References: 327 . 1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 328 SIAM Review, 57(4), 2015 329 330 .seealso: SNESCreate(), SNES, SNESSetType(), SNESSetNGS(), SNESType (for list of available types), SNESNGSSetSweeps(), SNESNGSSetTolerances(), 331 SNESSetNormSchedule() 332 M*/ 333 334 PETSC_EXTERN PetscErrorCode SNESCreate_NGS(SNES snes) 335 { 336 SNES_NGS *gs; 337 PetscErrorCode ierr; 338 339 PetscFunctionBegin; 340 snes->ops->destroy = SNESDestroy_NGS; 341 snes->ops->setup = SNESSetUp_NGS; 342 snes->ops->setfromoptions = SNESSetFromOptions_NGS; 343 snes->ops->view = SNESView_NGS; 344 snes->ops->solve = SNESSolve_NGS; 345 snes->ops->reset = SNESReset_NGS; 346 347 snes->usesksp = PETSC_FALSE; 348 snes->usesnpc = PETSC_FALSE; 349 350 snes->alwayscomputesfinalresidual = PETSC_FALSE; 351 352 if (!snes->tolerancesset) { 353 snes->max_its = 10000; 354 snes->max_funcs = 10000; 355 } 356 357 ierr = PetscNewLog(snes,&gs);CHKERRQ(ierr); 358 359 gs->sweeps = 1; 360 gs->rtol = 1e-5; 361 gs->abstol = PETSC_MACHINE_EPSILON; 362 gs->stol = 1000*PETSC_MACHINE_EPSILON; 363 gs->max_its = 50; 364 gs->h = PETSC_SQRT_MACHINE_EPSILON; 365 366 snes->data = (void*) gs; 367 PetscFunctionReturn(0); 368 } 369