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 - Number of sweeps of GS to apply 310 . -snes_ngs_atol - Absolute residual tolerance for GS iteration 311 . -snes_ngs_rtol - Relative residual tolerance for GS iteration 312 . -snes_ngs_stol - Absolute update tolerance for GS iteration 313 . -snes_ngs_max_it - 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 - 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 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 References: 326 . 1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 327 SIAM Review, 57(4), 2015 328 329 .seealso: SNESCreate(), SNES, SNESSetType(), SNESSetNGS(), SNESType (for list of available types), SNESNGSSetSweeps(), SNESNGSSetTolerances() 330 M*/ 331 332 PETSC_EXTERN PetscErrorCode SNESCreate_NGS(SNES snes) 333 { 334 SNES_NGS *gs; 335 PetscErrorCode ierr; 336 337 PetscFunctionBegin; 338 snes->ops->destroy = SNESDestroy_NGS; 339 snes->ops->setup = SNESSetUp_NGS; 340 snes->ops->setfromoptions = SNESSetFromOptions_NGS; 341 snes->ops->view = SNESView_NGS; 342 snes->ops->solve = SNESSolve_NGS; 343 snes->ops->reset = SNESReset_NGS; 344 345 snes->usesksp = PETSC_FALSE; 346 snes->usesnpc = PETSC_FALSE; 347 348 snes->alwayscomputesfinalresidual = PETSC_FALSE; 349 350 if (!snes->tolerancesset) { 351 snes->max_its = 10000; 352 snes->max_funcs = 10000; 353 } 354 355 ierr = PetscNewLog(snes,&gs);CHKERRQ(ierr); 356 357 gs->sweeps = 1; 358 gs->rtol = 1e-5; 359 gs->abstol = 1e-15; 360 gs->stol = 1e-12; 361 gs->max_its = 50; 362 gs->h = 1e-8; 363 364 snes->data = (void*) gs; 365 PetscFunctionReturn(0); 366 } 367