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