1 #include <../src/snes/impls/ngmres/snesngmres.h> /*I "petscsnes.h" I*/ 2 #include <petscblaslapack.h> 3 #include <petscdm.h> 4 5 const char *const SNESNGMRESRestartTypes[] = {"NONE","PERIODIC","DIFFERENCE","SNESNGMRESRestartType","SNES_NGMRES_RESTART_",0}; 6 const char *const SNESNGMRESSelectTypes[] = {"NONE","DIFFERENCE","LINESEARCH","SNESNGMRESSelectType","SNES_NGMRES_SELECT_",0}; 7 8 PetscErrorCode SNESReset_NGMRES(SNES snes) 9 { 10 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 11 PetscErrorCode ierr; 12 13 PetscFunctionBegin; 14 ierr = VecDestroyVecs(ngmres->msize,&ngmres->Fdot);CHKERRQ(ierr); 15 ierr = VecDestroyVecs(ngmres->msize,&ngmres->Xdot);CHKERRQ(ierr); 16 ierr = SNESLineSearchDestroy(&ngmres->additive_linesearch);CHKERRQ(ierr); 17 PetscFunctionReturn(0); 18 } 19 20 PetscErrorCode SNESDestroy_NGMRES(SNES snes) 21 { 22 PetscErrorCode ierr; 23 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 24 25 PetscFunctionBegin; 26 ierr = SNESReset_NGMRES(snes);CHKERRQ(ierr); 27 ierr = PetscFree5(ngmres->h,ngmres->beta,ngmres->xi,ngmres->fnorms,ngmres->q);CHKERRQ(ierr); 28 ierr = PetscFree(ngmres->s);CHKERRQ(ierr); 29 ierr = PetscFree(ngmres->xnorms);CHKERRQ(ierr); 30 #if defined(PETSC_USE_COMPLEX) 31 ierr = PetscFree(ngmres->rwork);CHKERRQ(ierr); 32 #endif 33 ierr = PetscFree(ngmres->work);CHKERRQ(ierr); 34 ierr = PetscFree(snes->data);CHKERRQ(ierr); 35 PetscFunctionReturn(0); 36 } 37 38 PetscErrorCode SNESSetUp_NGMRES(SNES snes) 39 { 40 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 41 const char *optionsprefix; 42 PetscInt msize,hsize; 43 PetscErrorCode ierr; 44 DM dm; 45 46 PetscFunctionBegin; 47 if (snes->pc && snes->pcside == PC_LEFT && snes->functype == SNES_FUNCTION_UNPRECONDITIONED) { 48 SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE,"SNESNGMRES does not support left preconditioning with unpreconditioned function"); 49 } 50 if (snes->pcside == PC_LEFT && snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_PRECONDITIONED; 51 ierr = SNESSetWorkVecs(snes,5);CHKERRQ(ierr); 52 53 if (!snes->vec_sol) { 54 ierr = SNESGetDM(snes,&dm);CHKERRQ(ierr); 55 ierr = DMCreateGlobalVector(dm,&snes->vec_sol);CHKERRQ(ierr); 56 } 57 58 if (!ngmres->Xdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Xdot);CHKERRQ(ierr);} 59 if (!ngmres->Fdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Fdot);CHKERRQ(ierr);} 60 if (!ngmres->setup_called) { 61 msize = ngmres->msize; /* restart size */ 62 hsize = msize * msize; 63 64 /* explicit least squares minimization solve */ 65 ierr = PetscMalloc5(hsize,&ngmres->h, msize,&ngmres->beta, msize,&ngmres->xi, msize,&ngmres->fnorms, hsize,&ngmres->q);CHKERRQ(ierr); 66 ierr = PetscMalloc1(msize,&ngmres->xnorms);CHKERRQ(ierr); 67 ngmres->nrhs = 1; 68 ngmres->lda = msize; 69 ngmres->ldb = msize; 70 ierr = PetscMalloc1(msize,&ngmres->s);CHKERRQ(ierr); 71 ierr = PetscMemzero(ngmres->h, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 72 ierr = PetscMemzero(ngmres->q, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 73 ierr = PetscMemzero(ngmres->xi, msize*sizeof(PetscScalar));CHKERRQ(ierr); 74 ierr = PetscMemzero(ngmres->beta,msize*sizeof(PetscScalar));CHKERRQ(ierr); 75 ngmres->lwork = 12*msize; 76 #if defined(PETSC_USE_COMPLEX) 77 ierr = PetscMalloc1(ngmres->lwork,&ngmres->rwork);CHKERRQ(ierr); 78 #endif 79 ierr = PetscMalloc1(ngmres->lwork,&ngmres->work);CHKERRQ(ierr); 80 } 81 82 /* linesearch setup */ 83 ierr = SNESGetOptionsPrefix(snes,&optionsprefix);CHKERRQ(ierr); 84 85 if (ngmres->select_type == SNES_NGMRES_SELECT_LINESEARCH) { 86 ierr = SNESLineSearchCreate(PetscObjectComm((PetscObject)snes),&ngmres->additive_linesearch);CHKERRQ(ierr); 87 ierr = SNESLineSearchSetSNES(ngmres->additive_linesearch,snes);CHKERRQ(ierr); 88 ierr = SNESLineSearchSetType(ngmres->additive_linesearch,SNESLINESEARCHL2);CHKERRQ(ierr); 89 ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,"additive_");CHKERRQ(ierr); 90 ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,optionsprefix);CHKERRQ(ierr); 91 ierr = SNESLineSearchSetFromOptions(ngmres->additive_linesearch);CHKERRQ(ierr); 92 } 93 94 ngmres->setup_called = PETSC_TRUE; 95 PetscFunctionReturn(0); 96 } 97 98 PetscErrorCode SNESSetFromOptions_NGMRES(PetscOptionItems *PetscOptionsObject,SNES snes) 99 { 100 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 101 PetscErrorCode ierr; 102 PetscBool debug = PETSC_FALSE; 103 SNESLineSearch linesearch; 104 105 PetscFunctionBegin; 106 ierr = PetscOptionsHead(PetscOptionsObject,"SNES NGMRES options");CHKERRQ(ierr); 107 ierr = PetscOptionsEnum("-snes_ngmres_select_type","Select type","SNESNGMRESSetSelectType",SNESNGMRESSelectTypes, 108 (PetscEnum)ngmres->select_type,(PetscEnum*)&ngmres->select_type,NULL);CHKERRQ(ierr); 109 ierr = PetscOptionsEnum("-snes_ngmres_restart_type","Restart type","SNESNGMRESSetRestartType",SNESNGMRESRestartTypes, 110 (PetscEnum)ngmres->restart_type,(PetscEnum*)&ngmres->restart_type,NULL);CHKERRQ(ierr); 111 ierr = PetscOptionsBool("-snes_ngmres_candidate", "Use candidate storage", "SNES",ngmres->candidate,&ngmres->candidate,NULL);CHKERRQ(ierr); 112 ierr = PetscOptionsBool("-snes_ngmres_approxfunc","Linearly approximate the function", "SNES",ngmres->approxfunc,&ngmres->approxfunc,NULL);CHKERRQ(ierr); 113 ierr = PetscOptionsInt("-snes_ngmres_m", "Number of directions", "SNES",ngmres->msize,&ngmres->msize,NULL);CHKERRQ(ierr); 114 ierr = PetscOptionsInt("-snes_ngmres_restart", "Iterations before forced restart", "SNES",ngmres->restart_periodic,&ngmres->restart_periodic,NULL);CHKERRQ(ierr); 115 ierr = PetscOptionsInt("-snes_ngmres_restart_it", "Tolerance iterations before restart","SNES",ngmres->restart_it,&ngmres->restart_it,NULL);CHKERRQ(ierr); 116 ierr = PetscOptionsBool("-snes_ngmres_monitor", "Monitor actions of NGMRES", "SNES",ngmres->monitor ? PETSC_TRUE : PETSC_FALSE,&debug,NULL);CHKERRQ(ierr); 117 if (debug) { 118 ngmres->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));CHKERRQ(ierr); 119 } 120 ierr = PetscOptionsReal("-snes_ngmres_gammaA", "Residual selection constant", "SNES",ngmres->gammaA,&ngmres->gammaA,NULL);CHKERRQ(ierr); 121 ierr = PetscOptionsReal("-snes_ngmres_gammaC", "Residual restart constant", "SNES",ngmres->gammaC,&ngmres->gammaC,NULL);CHKERRQ(ierr); 122 ierr = PetscOptionsReal("-snes_ngmres_epsilonB", "Difference selection constant", "SNES",ngmres->epsilonB,&ngmres->epsilonB,NULL);CHKERRQ(ierr); 123 ierr = PetscOptionsReal("-snes_ngmres_deltaB", "Difference residual selection constant", "SNES",ngmres->deltaB,&ngmres->deltaB,NULL);CHKERRQ(ierr); 124 ierr = PetscOptionsBool("-snes_ngmres_single_reduction", "Aggregate reductions", "SNES",ngmres->singlereduction,&ngmres->singlereduction,NULL);CHKERRQ(ierr); 125 ierr = PetscOptionsBool("-snes_ngmres_restart_fm_rise", "Restart on F_M residual rise", "SNESNGMRESSetRestartFmRise",ngmres->restart_fm_rise,&ngmres->restart_fm_rise,NULL);CHKERRQ(ierr); 126 ierr = PetscOptionsTail();CHKERRQ(ierr); 127 if ((ngmres->gammaA > ngmres->gammaC) && (ngmres->gammaC > 2.)) ngmres->gammaC = ngmres->gammaA; 128 129 /* set the default type of the line search if the user hasn't already. */ 130 if (!snes->linesearch) { 131 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 132 ierr = SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);CHKERRQ(ierr); 133 } 134 PetscFunctionReturn(0); 135 } 136 137 PetscErrorCode SNESView_NGMRES(SNES snes,PetscViewer viewer) 138 { 139 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 140 PetscBool iascii; 141 PetscErrorCode ierr; 142 143 PetscFunctionBegin; 144 ierr = PetscObjectTypeCompare((PetscObject) viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 145 if (iascii) { 146 ierr = PetscViewerASCIIPrintf(viewer," Number of stored past updates: %d\n", ngmres->msize);CHKERRQ(ierr); 147 ierr = PetscViewerASCIIPrintf(viewer," Residual selection: gammaA=%1.0e, gammaC=%1.0e\n",ngmres->gammaA,ngmres->gammaC);CHKERRQ(ierr); 148 ierr = PetscViewerASCIIPrintf(viewer," Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n",ngmres->epsilonB,ngmres->deltaB);CHKERRQ(ierr); 149 ierr = PetscViewerASCIIPrintf(viewer," Restart on F_M residual increase: %s\n",ngmres->restart_fm_rise?"TRUE":"FALSE");CHKERRQ(ierr); 150 } 151 PetscFunctionReturn(0); 152 } 153 154 PetscErrorCode SNESSolve_NGMRES(SNES snes) 155 { 156 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 157 /* present solution, residual, and preconditioned residual */ 158 Vec X,F,B,D,Y; 159 160 /* candidate linear combination answers */ 161 Vec XA,FA,XM,FM; 162 163 /* coefficients and RHS to the minimization problem */ 164 PetscReal fnorm,fMnorm,fAnorm; 165 PetscReal xnorm,xMnorm,xAnorm; 166 PetscReal ynorm,yMnorm,yAnorm; 167 PetscInt k,k_restart,l,ivec,restart_count = 0; 168 169 /* solution selection data */ 170 PetscBool selectRestart; 171 /* 172 These two variables are initialized to prevent compilers/analyzers from producing false warnings about these variables being passed 173 to SNESNGMRESSelect_Private() without being set when SNES_NGMRES_RESTART_DIFFERENCE, the values are not used in the subroutines in that case 174 so the code is correct as written. 175 */ 176 PetscReal dnorm = 0.0,dminnorm = 0.0; 177 PetscReal fminnorm; 178 179 SNESConvergedReason reason; 180 SNESLineSearchReason lssucceed; 181 PetscErrorCode ierr; 182 183 PetscFunctionBegin; 184 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); 185 186 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 187 /* variable initialization */ 188 snes->reason = SNES_CONVERGED_ITERATING; 189 X = snes->vec_sol; 190 F = snes->vec_func; 191 B = snes->vec_rhs; 192 XA = snes->vec_sol_update; 193 FA = snes->work[0]; 194 D = snes->work[1]; 195 196 /* work for the line search */ 197 Y = snes->work[2]; 198 XM = snes->work[3]; 199 FM = snes->work[4]; 200 201 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 202 snes->iter = 0; 203 snes->norm = 0.; 204 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 205 206 /* initialization */ 207 208 if (snes->pc && snes->pcside == PC_LEFT) { 209 ierr = SNESApplyNPC(snes,X,NULL,F);CHKERRQ(ierr); 210 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 211 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 212 snes->reason = SNES_DIVERGED_INNER; 213 PetscFunctionReturn(0); 214 } 215 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 216 } else { 217 if (!snes->vec_func_init_set) { 218 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 219 } else snes->vec_func_init_set = PETSC_FALSE; 220 221 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 222 SNESCheckFunctionNorm(snes,fnorm); 223 } 224 fminnorm = fnorm; 225 226 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 227 snes->norm = fnorm; 228 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 229 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 230 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 231 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 232 if (snes->reason) PetscFunctionReturn(0); 233 SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X); 234 235 k_restart = 1; 236 l = 1; 237 ivec = 0; 238 for (k=1; k < snes->max_its+1; k++) { 239 /* Computation of x^M */ 240 if (snes->pc && snes->pcside == PC_RIGHT) { 241 ierr = VecCopy(X,XM);CHKERRQ(ierr); 242 ierr = SNESSetInitialFunction(snes->pc,F);CHKERRQ(ierr); 243 244 ierr = PetscLogEventBegin(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 245 ierr = SNESSolve(snes->pc,B,XM);CHKERRQ(ierr); 246 ierr = PetscLogEventEnd(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 247 248 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 249 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 250 snes->reason = SNES_DIVERGED_INNER; 251 PetscFunctionReturn(0); 252 } 253 ierr = SNESGetNPCFunction(snes,FM,&fMnorm);CHKERRQ(ierr); 254 } else { 255 /* no preconditioner -- just take gradient descent with line search */ 256 ierr = VecCopy(F,Y);CHKERRQ(ierr); 257 ierr = VecCopy(F,FM);CHKERRQ(ierr); 258 ierr = VecCopy(X,XM);CHKERRQ(ierr); 259 260 fMnorm = fnorm; 261 262 ierr = SNESLineSearchApply(snes->linesearch,XM,FM,&fMnorm,Y);CHKERRQ(ierr); 263 ierr = SNESLineSearchGetReason(snes->linesearch,&lssucceed);CHKERRQ(ierr); 264 if (lssucceed) { 265 if (++snes->numFailures >= snes->maxFailures) { 266 snes->reason = SNES_DIVERGED_LINE_SEARCH; 267 PetscFunctionReturn(0); 268 } 269 } 270 } 271 272 ierr = SNESNGMRESFormCombinedSolution_Private(snes,ivec,l,XM,FM,fMnorm,X,XA,FA);CHKERRQ(ierr); 273 /* r = F(x) */ 274 if (fminnorm > fMnorm) fminnorm = fMnorm; /* the minimum norm is now of F^M */ 275 276 /* differences for selection and restart */ 277 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE || ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) { 278 ierr = SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,&dnorm,&dminnorm,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);CHKERRQ(ierr); 279 } else { 280 ierr = SNESNGMRESNorms_Private(snes,l,X,F,XM,FM,XA,FA,D,NULL,NULL,&xMnorm,NULL,&yMnorm,&xAnorm,&fAnorm,&yAnorm);CHKERRQ(ierr); 281 } 282 SNESCheckFunctionNorm(snes,fnorm); 283 284 /* combination (additive) or selection (multiplicative) of the N-GMRES solution */ 285 ierr = SNESNGMRESSelect_Private(snes,k_restart,XM,FM,xMnorm,fMnorm,yMnorm,XA,FA,xAnorm,fAnorm,yAnorm,dnorm,fminnorm,dminnorm,X,F,Y,&xnorm,&fnorm,&ynorm);CHKERRQ(ierr); 286 selectRestart = PETSC_FALSE; 287 288 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE) { 289 ierr = SNESNGMRESSelectRestart_Private(snes,l,fMnorm,fAnorm,dnorm,fminnorm,dminnorm,&selectRestart);CHKERRQ(ierr); 290 291 /* if the restart conditions persist for more than restart_it iterations, restart. */ 292 if (selectRestart) restart_count++; 293 else restart_count = 0; 294 } else if (ngmres->restart_type == SNES_NGMRES_RESTART_PERIODIC) { 295 if (k_restart > ngmres->restart_periodic) { 296 if (ngmres->monitor) ierr = PetscViewerASCIIPrintf(ngmres->monitor,"periodic restart after %D iterations\n",k_restart);CHKERRQ(ierr); 297 restart_count = ngmres->restart_it; 298 } 299 } 300 301 ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */ 302 303 /* restart after restart conditions have persisted for a fixed number of iterations */ 304 if (restart_count >= ngmres->restart_it) { 305 if (ngmres->monitor) { 306 ierr = PetscViewerASCIIPrintf(ngmres->monitor,"Restarted at iteration %d\n",k_restart);CHKERRQ(ierr); 307 } 308 restart_count = 0; 309 k_restart = 1; 310 l = 1; 311 ivec = 0; 312 /* q_{00} = nu */ 313 ierr = SNESNGMRESUpdateSubspace_Private(snes,0,0,FM,fMnorm,XM);CHKERRQ(ierr); 314 } else { 315 /* select the current size of the subspace */ 316 if (l < ngmres->msize) l++; 317 k_restart++; 318 /* place the current entry in the list of previous entries */ 319 if (ngmres->candidate) { 320 if (fminnorm > fMnorm) fminnorm = fMnorm; 321 ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,FM,fMnorm,XM);CHKERRQ(ierr); 322 } else { 323 if (fminnorm > fnorm) fminnorm = fnorm; 324 ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,F,fnorm,X);CHKERRQ(ierr); 325 } 326 } 327 328 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 329 snes->iter = k; 330 snes->norm = fnorm; 331 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 332 ierr = SNESLogConvergenceHistory(snes,snes->norm,snes->iter);CHKERRQ(ierr); 333 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 334 ierr = (*snes->ops->converged)(snes,snes->iter,0,0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 335 if (snes->reason) PetscFunctionReturn(0); 336 } 337 snes->reason = SNES_DIVERGED_MAX_IT; 338 PetscFunctionReturn(0); 339 } 340 341 /*@ 342 SNESNGMRESSetRestartFmRise - Increase the restart count if the step x_M increases the residual F_M 343 344 Input Parameters: 345 + snes - the SNES context. 346 - flg - boolean value deciding whether to use the option or not 347 348 Options Database: 349 + -snes_ngmres_restart_fm_rise - Increase the restart count if the step x_M increases the residual F_M 350 351 Level: intermediate 352 353 Notes: 354 If the proposed step x_M increases the residual F_M, it might be trying to get out of a stagnation area. 355 To help the solver do that, reset the Krylov subspace whenever F_M increases. 356 357 This option must be used with SNES_NGMRES_RESTART_DIFFERENCE 358 359 The default is FALSE. 360 .seealso: SNES_NGMRES_RESTART_DIFFERENCE 361 @*/ 362 PetscErrorCode SNESNGMRESSetRestartFmRise(SNES snes,PetscBool flg) 363 { 364 PetscErrorCode (*f)(SNES,PetscBool); 365 PetscErrorCode ierr; 366 367 PetscFunctionBegin; 368 ierr = PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",&f);CHKERRQ(ierr); 369 if (f) {ierr = (f)(snes,flg);CHKERRQ(ierr);} 370 PetscFunctionReturn(0); 371 } 372 373 PetscErrorCode SNESNGMRESSetRestartFmRise_NGMRES(SNES snes,PetscBool flg) 374 { 375 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 376 377 PetscFunctionBegin; 378 ngmres->restart_fm_rise = flg; 379 PetscFunctionReturn(0); 380 } 381 382 PetscErrorCode SNESNGMRESGetRestartFmRise(SNES snes,PetscBool *flg) 383 { 384 PetscErrorCode (*f)(SNES,PetscBool*); 385 PetscErrorCode ierr; 386 387 PetscFunctionBegin; 388 ierr = PetscObjectQueryFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",&f);CHKERRQ(ierr); 389 if (f) {ierr = (f)(snes,flg);CHKERRQ(ierr);} 390 PetscFunctionReturn(0); 391 } 392 393 PetscErrorCode SNESNGMRESGetRestartFmRise_NGMRES(SNES snes,PetscBool *flg) 394 { 395 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 396 397 PetscFunctionBegin; 398 *flg = ngmres->restart_fm_rise; 399 PetscFunctionReturn(0); 400 } 401 402 403 /*@ 404 SNESNGMRESSetRestartType - Sets the restart type for SNESNGMRES. 405 406 Logically Collective on SNES 407 408 Input Parameters: 409 + snes - the iterative context 410 - rtype - restart type 411 412 Options Database: 413 + -snes_ngmres_restart_type<difference,periodic,none> - set the restart type 414 - -snes_ngmres_restart[30] - sets the number of iterations before restart for periodic 415 416 Level: intermediate 417 418 SNESNGMRESRestartTypes: 419 + SNES_NGMRES_RESTART_NONE - never restart 420 . SNES_NGMRES_RESTART_DIFFERENCE - restart based upon difference criteria 421 - SNES_NGMRES_RESTART_PERIODIC - restart after a fixed number of iterations 422 423 Notes: 424 The default line search used is the L2 line search and it requires two additional function evaluations. 425 426 .keywords: SNES, SNESNGMRES, restart, type, set SNESLineSearch 427 @*/ 428 PetscErrorCode SNESNGMRESSetRestartType(SNES snes,SNESNGMRESRestartType rtype) 429 { 430 PetscErrorCode ierr; 431 432 PetscFunctionBegin; 433 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 434 ierr = PetscTryMethod(snes,"SNESNGMRESSetRestartType_C",(SNES,SNESNGMRESRestartType),(snes,rtype));CHKERRQ(ierr); 435 PetscFunctionReturn(0); 436 } 437 438 /*@ 439 SNESNGMRESSetSelectType - Sets the selection type for SNESNGMRES. This determines how the candidate solution and 440 combined solution are used to create the next iterate. 441 442 Logically Collective on SNES 443 444 Input Parameters: 445 + snes - the iterative context 446 - stype - selection type 447 448 Options Database: 449 . -snes_ngmres_select_type<difference,none,linesearch> 450 451 Level: intermediate 452 453 SNESNGMRESSelectTypes: 454 + SNES_NGMRES_SELECT_NONE - choose the combined solution all the time 455 . SNES_NGMRES_SELECT_DIFFERENCE - choose based upon the selection criteria 456 - SNES_NGMRES_SELECT_LINESEARCH - choose based upon line search combination 457 458 Notes: 459 The default line search used is the L2 line search and it requires two additional function evaluations. 460 461 .keywords: SNES, SNESNGMRES, selection, type, set SNESLineSearch 462 @*/ 463 PetscErrorCode SNESNGMRESSetSelectType(SNES snes,SNESNGMRESSelectType stype) 464 { 465 PetscErrorCode ierr; 466 467 PetscFunctionBegin; 468 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 469 ierr = PetscTryMethod(snes,"SNESNGMRESSetSelectType_C",(SNES,SNESNGMRESSelectType),(snes,stype));CHKERRQ(ierr); 470 PetscFunctionReturn(0); 471 } 472 473 PetscErrorCode SNESNGMRESSetSelectType_NGMRES(SNES snes,SNESNGMRESSelectType stype) 474 { 475 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 476 477 PetscFunctionBegin; 478 ngmres->select_type = stype; 479 PetscFunctionReturn(0); 480 } 481 482 PetscErrorCode SNESNGMRESSetRestartType_NGMRES(SNES snes,SNESNGMRESRestartType rtype) 483 { 484 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 485 486 PetscFunctionBegin; 487 ngmres->restart_type = rtype; 488 PetscFunctionReturn(0); 489 } 490 491 /*MC 492 SNESNGMRES - The Nonlinear Generalized Minimum Residual method. 493 494 Level: beginner 495 496 Options Database: 497 + -snes_ngmres_select_type<difference,none,linesearch> - choose the select between candidate and combined solution 498 . -snes_ngmres_restart_type<difference,none,periodic> - choose the restart conditions 499 . -snes_ngmres_candidate - Use NGMRES variant which combines candidate solutions instead of actual solutions 500 . -snes_ngmres_m - Number of stored previous solutions and residuals 501 . -snes_ngmres_restart_it - Number of iterations the restart conditions hold before restart 502 . -snes_ngmres_gammaA - Residual tolerance for solution select between the candidate and combination 503 . -snes_ngmres_gammaC - Residual tolerance for restart 504 . -snes_ngmres_epsilonB - Difference tolerance between subsequent solutions triggering restart 505 . -snes_ngmres_deltaB - Difference tolerance between residuals triggering restart 506 . -snes_ngmres_restart_fm_rise - Restart on residual rise from x_M step 507 . -snes_ngmres_monitor - Prints relevant information about the ngmres iteration 508 . -snes_linesearch_type <basic,l2,cp> - Line search type used for the default smoother 509 - -additive_snes_linesearch_type - linesearch type used to select between the candidate and combined solution with additive select type 510 511 Notes: 512 513 The N-GMRES method combines m previous solutions into a minimum-residual solution by solving a small linearized 514 optimization problem at each iteration. 515 516 Very similar to the SNESANDERSON algorithm. 517 518 References: 519 + 1. - C. W. Oosterlee and T. Washio, "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows", 520 SIAM Journal on Scientific Computing, 21(5), 2000. 521 - 2. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu, "Composing Scalable Nonlinear Algebraic Solvers", 522 SIAM Review, 57(4), 2015 523 524 525 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 526 M*/ 527 528 PETSC_EXTERN PetscErrorCode SNESCreate_NGMRES(SNES snes) 529 { 530 SNES_NGMRES *ngmres; 531 PetscErrorCode ierr; 532 533 PetscFunctionBegin; 534 snes->ops->destroy = SNESDestroy_NGMRES; 535 snes->ops->setup = SNESSetUp_NGMRES; 536 snes->ops->setfromoptions = SNESSetFromOptions_NGMRES; 537 snes->ops->view = SNESView_NGMRES; 538 snes->ops->solve = SNESSolve_NGMRES; 539 snes->ops->reset = SNESReset_NGMRES; 540 541 snes->usespc = PETSC_TRUE; 542 snes->usesksp = PETSC_FALSE; 543 snes->pcside = PC_RIGHT; 544 545 snes->alwayscomputesfinalresidual = PETSC_TRUE; 546 547 ierr = PetscNewLog(snes,&ngmres);CHKERRQ(ierr); 548 snes->data = (void*) ngmres; 549 ngmres->msize = 30; 550 551 if (!snes->tolerancesset) { 552 snes->max_funcs = 30000; 553 snes->max_its = 10000; 554 } 555 556 ngmres->candidate = PETSC_FALSE; 557 558 ngmres->additive_linesearch = NULL; 559 ngmres->approxfunc = PETSC_FALSE; 560 ngmres->restart_it = 2; 561 ngmres->restart_periodic = 30; 562 ngmres->gammaA = 2.0; 563 ngmres->gammaC = 2.0; 564 ngmres->deltaB = 0.9; 565 ngmres->epsilonB = 0.1; 566 ngmres->restart_fm_rise = PETSC_FALSE; 567 568 ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE; 569 ngmres->select_type = SNES_NGMRES_SELECT_DIFFERENCE; 570 571 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);CHKERRQ(ierr); 572 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);CHKERRQ(ierr); 573 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartFmRise_C",SNESNGMRESSetRestartFmRise_NGMRES);CHKERRQ(ierr); 574 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESGetRestartFmRise_C",SNESNGMRESGetRestartFmRise_NGMRES);CHKERRQ(ierr); 575 PetscFunctionReturn(0); 576 } 577