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 PCSide npcside; 184 PetscErrorCode ierr; 185 186 PetscFunctionBegin; 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 = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 202 snes->iter = 0; 203 snes->norm = 0.; 204 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 205 206 /* initialization */ 207 208 if (snes->pc && snes->pcside == PC_LEFT) { 209 ierr = SNESApplyPC(snes,X,NULL,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 if (snes->domainerror) { 220 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 221 PetscFunctionReturn(0); 222 } 223 } else snes->vec_func_init_set = PETSC_FALSE; 224 if (!snes->norm_init_set) { 225 /* convergence test */ 226 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 227 if (PetscIsInfOrNanReal(fnorm)) { 228 snes->reason = SNES_DIVERGED_FNORM_NAN; 229 PetscFunctionReturn(0); 230 } 231 } else { 232 fnorm = snes->norm_init; 233 snes->norm_init_set = PETSC_FALSE; 234 } 235 } 236 fminnorm = fnorm; 237 238 /* q_{00} = nu */ 239 SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X); 240 241 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 242 snes->norm = fnorm; 243 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 244 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 245 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 246 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 247 if (snes->reason) PetscFunctionReturn(0); 248 249 k_restart = 1; 250 l = 1; 251 for (k=1; k < snes->max_its+1; k++) { 252 /* select which vector of the stored subspace will be updated */ 253 ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */ 254 255 /* Computation of x^M */ 256 if (snes->pc && snes->pcside == PC_RIGHT) { 257 ierr = VecCopy(X,XM);CHKERRQ(ierr); 258 ierr = SNESSetInitialFunction(snes->pc,F);CHKERRQ(ierr); 259 ierr = SNESSetInitialFunctionNorm(snes->pc,fnorm);CHKERRQ(ierr); 260 261 ierr = PetscLogEventBegin(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 262 ierr = SNESSolve(snes->pc,B,XM);CHKERRQ(ierr); 263 ierr = PetscLogEventEnd(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 264 265 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 266 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 267 snes->reason = SNES_DIVERGED_INNER; 268 PetscFunctionReturn(0); 269 } 270 ierr = SNESGetPCFunction(snes,FM,&fMnorm);CHKERRQ(ierr); 271 } else { 272 /* no preconditioner -- just take gradient descent with line search */ 273 ierr = VecCopy(F,Y);CHKERRQ(ierr); 274 ierr = VecCopy(F,FM);CHKERRQ(ierr); 275 ierr = VecCopy(X,XM);CHKERRQ(ierr); 276 277 fMnorm = fnorm; 278 279 ierr = SNESLineSearchApply(snes->linesearch,XM,FM,&fMnorm,Y);CHKERRQ(ierr); 280 ierr = SNESLineSearchGetSuccess(snes->linesearch,&lssucceed);CHKERRQ(ierr); 281 if (!lssucceed) { 282 if (++snes->numFailures >= snes->maxFailures) { 283 snes->reason = SNES_DIVERGED_LINE_SEARCH; 284 PetscFunctionReturn(0); 285 } 286 } 287 } 288 ierr = SNESNGMRESFormCombinedSolution_Private(snes,l,XM,FM,fMnorm,X,XA,FA);CHKERRQ(ierr); 289 /* r = F(x) */ 290 if (fminnorm > fMnorm) fminnorm = fMnorm; /* the minimum norm is now of F^M */ 291 292 /* differences for selection and restart */ 293 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE || ngmres->select_type == SNES_NGMRES_SELECT_DIFFERENCE) { 294 ierr = SNESNGMRESCalculateDifferences_Private(snes,l,X,F,XM,FM,XA,FA,D,&dnorm,&dminnorm,&fAnorm);CHKERRQ(ierr); 295 } else { 296 ierr = VecNorm(FA,NORM_2,&fAnorm);CHKERRQ(ierr); 297 } 298 if (PetscIsInfOrNanReal(fAnorm)) { 299 snes->reason = SNES_DIVERGED_FNORM_NAN; 300 PetscFunctionReturn(0); 301 } 302 303 /* combination (additive) or selection (multiplicative) of the N-GMRES solution */ 304 ierr = SNESNGMRESSelect_Private(snes,k_restart,XM,FM,fMnorm,XA,FA,fAnorm,dnorm,fminnorm,dminnorm,X,F,Y,&fnorm);CHKERRQ(ierr); 305 selectRestart = PETSC_FALSE; 306 if (ngmres->restart_type == SNES_NGMRES_RESTART_DIFFERENCE) { 307 ierr = SNESNGMRESSelectRestart_Private(snes,l,fAnorm,dnorm,fminnorm,dminnorm,&selectRestart);CHKERRQ(ierr); 308 /* if the restart conditions persist for more than restart_it iterations, restart. */ 309 if (selectRestart) restart_count++; 310 else restart_count = 0; 311 } else if (ngmres->restart_type == SNES_NGMRES_RESTART_PERIODIC) { 312 if (k_restart > ngmres->restart_periodic) { 313 if (ngmres->monitor) ierr = PetscViewerASCIIPrintf(ngmres->monitor,"periodic restart after %D iterations\n",k_restart);CHKERRQ(ierr); 314 restart_count = ngmres->restart_it; 315 } 316 } 317 /* restart after restart conditions have persisted for a fixed number of iterations */ 318 if (restart_count >= ngmres->restart_it) { 319 if (ngmres->monitor) { 320 ierr = PetscViewerASCIIPrintf(ngmres->monitor,"Restarted at iteration %d\n",k_restart);CHKERRQ(ierr); 321 } 322 restart_count = 0; 323 k_restart = 1; 324 l = 1; 325 /* q_{00} = nu */ 326 if (ngmres->candidate) { 327 ierr = SNESNGMRESUpdateSubspace_Private(snes,0,0,FM,fMnorm,XM);CHKERRQ(ierr); 328 } else { 329 ierr = SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fMnorm,X);CHKERRQ(ierr); 330 } 331 } else { 332 /* select the current size of the subspace */ 333 if (l < ngmres->msize) l++; 334 k_restart++; 335 /* place the current entry in the list of previous entries */ 336 if (ngmres->candidate) { 337 if (fminnorm > fMnorm) fminnorm = fMnorm; 338 ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,FM,fMnorm,XM);CHKERRQ(ierr); 339 } else { 340 if (fminnorm > fnorm) fminnorm = fnorm; 341 ierr = SNESNGMRESUpdateSubspace_Private(snes,ivec,l,F,fnorm,X);CHKERRQ(ierr); 342 } 343 } 344 345 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 346 snes->iter = k; 347 snes->norm = fnorm; 348 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 349 ierr = SNESLogConvergenceHistory(snes,snes->norm,snes->iter);CHKERRQ(ierr); 350 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 351 ierr = VecNormBegin(Y,NORM_2,&ynorm);CHKERRQ(ierr); 352 ierr = VecNormBegin(X,NORM_2,&xnorm);CHKERRQ(ierr); 353 ierr = VecNormEnd(Y,NORM_2,&ynorm);CHKERRQ(ierr); 354 ierr = VecNormEnd(X,NORM_2,&xnorm);CHKERRQ(ierr); 355 ierr = (*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 356 if (snes->reason) PetscFunctionReturn(0); 357 } 358 snes->reason = SNES_DIVERGED_MAX_IT; 359 PetscFunctionReturn(0); 360 } 361 362 #undef __FUNCT__ 363 #define __FUNCT__ "SNESNGMRESSetRestartType" 364 /*@ 365 SNESNGMRESSetRestartType - Sets the restart type for SNESNGMRES. 366 367 Logically Collective on SNES 368 369 Input Parameters: 370 + snes - the iterative context 371 - rtype - restart type 372 373 Options Database: 374 + -snes_ngmres_restart_type<difference,periodic,none> - set the restart type 375 - -snes_ngmres_restart[30] - sets the number of iterations before restart for periodic 376 377 Level: intermediate 378 379 SNESNGMRESRestartTypes: 380 + SNES_NGMRES_RESTART_NONE - never restart 381 . SNES_NGMRES_RESTART_DIFFERENCE - restart based upon difference criteria 382 - SNES_NGMRES_RESTART_PERIODIC - restart after a fixed number of iterations 383 384 Notes: 385 The default line search used is the L2 line search and it requires two additional function evaluations. 386 387 .keywords: SNES, SNESNGMRES, restart, type, set SNESLineSearch 388 @*/ 389 PetscErrorCode SNESNGMRESSetRestartType(SNES snes,SNESNGMRESRestartType rtype) 390 { 391 PetscErrorCode ierr; 392 393 PetscFunctionBegin; 394 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 395 ierr = PetscTryMethod(snes,"SNESNGMRESSetRestartType_C",(SNES,SNESNGMRESRestartType),(snes,rtype));CHKERRQ(ierr); 396 PetscFunctionReturn(0); 397 } 398 399 #undef __FUNCT__ 400 #define __FUNCT__ "SNESNGMRESSetSelectType" 401 /*@ 402 SNESNGMRESSetSelectType - Sets the selection type for SNESNGMRES. This determines how the candidate solution and 403 combined solution are used to create the next iterate. 404 405 Logically Collective on SNES 406 407 Input Parameters: 408 + snes - the iterative context 409 - stype - selection type 410 411 Options Database: 412 . -snes_ngmres_select_type<difference,none,linesearch> 413 414 Level: intermediate 415 416 SNESNGMRESSelectTypes: 417 + SNES_NGMRES_SELECT_NONE - choose the combined solution all the time 418 . SNES_NGMRES_SELECT_DIFFERENCE - choose based upon the selection criteria 419 - SNES_NGMRES_SELECT_LINESEARCH - choose based upon line search combination 420 421 Notes: 422 The default line search used is the L2 line search and it requires two additional function evaluations. 423 424 .keywords: SNES, SNESNGMRES, selection, type, set SNESLineSearch 425 @*/ 426 PetscErrorCode SNESNGMRESSetSelectType(SNES snes,SNESNGMRESSelectType stype) 427 { 428 PetscErrorCode ierr; 429 430 PetscFunctionBegin; 431 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 432 ierr = PetscTryMethod(snes,"SNESNGMRESSetSelectType_C",(SNES,SNESNGMRESSelectType),(snes,stype));CHKERRQ(ierr); 433 PetscFunctionReturn(0); 434 } 435 436 #undef __FUNCT__ 437 #define __FUNCT__ "SNESNGMRESSetSelectType_NGMRES" 438 PetscErrorCode SNESNGMRESSetSelectType_NGMRES(SNES snes,SNESNGMRESSelectType stype) 439 { 440 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 441 442 PetscFunctionBegin; 443 ngmres->select_type = stype; 444 PetscFunctionReturn(0); 445 } 446 447 #undef __FUNCT__ 448 #define __FUNCT__ "SNESNGMRESSetRestartType_NGMRES" 449 PetscErrorCode SNESNGMRESSetRestartType_NGMRES(SNES snes,SNESNGMRESRestartType rtype) 450 { 451 SNES_NGMRES *ngmres = (SNES_NGMRES*)snes->data; 452 453 PetscFunctionBegin; 454 ngmres->restart_type = rtype; 455 PetscFunctionReturn(0); 456 } 457 458 /*MC 459 SNESNGMRES - The Nonlinear Generalized Minimum Residual method. 460 461 Level: beginner 462 463 Options Database: 464 + -snes_ngmres_select_type<difference,none,linesearch> - choose the select between candidate and combined solution 465 . -snes_ngmres_restart_type<difference,none,periodic> - choose the restart conditions 466 . -snes_ngmres_candidate - Use NGMRES variant which combines candidate solutions instead of actual solutions 467 . -snes_ngmres_m - Number of stored previous solutions and residuals 468 . -snes_ngmres_restart_it - Number of iterations the restart conditions hold before restart 469 . -snes_ngmres_gammaA - Residual tolerance for solution select between the candidate and combination 470 . -snes_ngmres_gammaC - Residual tolerance for restart 471 . -snes_ngmres_epsilonB - Difference tolerance between subsequent solutions triggering restart 472 . -snes_ngmres_deltaB - Difference tolerance between residuals triggering restart 473 . -snes_ngmres_monitor - Prints relevant information about the ngmres iteration 474 . -snes_linesearch_type <basic,l2,cp> - Line search type used for the default smoother 475 - -additive_snes_linesearch_type - linesearch type used to select between the candidate and combined solution with additive select type 476 477 Notes: 478 479 The N-GMRES method combines m previous solutions into a minimum-residual solution by solving a small linearized 480 optimization problem at each iteration. 481 482 References: 483 484 "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows", C. W. Oosterlee and T. Washio, 485 SIAM Journal on Scientific Computing, 21(5), 2000. 486 487 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 488 M*/ 489 490 #undef __FUNCT__ 491 #define __FUNCT__ "SNESCreate_NGMRES" 492 PETSC_EXTERN PetscErrorCode SNESCreate_NGMRES(SNES snes) 493 { 494 SNES_NGMRES *ngmres; 495 PetscErrorCode ierr; 496 497 PetscFunctionBegin; 498 snes->ops->destroy = SNESDestroy_NGMRES; 499 snes->ops->setup = SNESSetUp_NGMRES; 500 snes->ops->setfromoptions = SNESSetFromOptions_NGMRES; 501 snes->ops->view = SNESView_NGMRES; 502 snes->ops->solve = SNESSolve_NGMRES; 503 snes->ops->reset = SNESReset_NGMRES; 504 505 snes->usespc = PETSC_TRUE; 506 snes->usesksp = PETSC_FALSE; 507 snes->pcside = PC_RIGHT; 508 509 ierr = PetscNewLog(snes,SNES_NGMRES,&ngmres);CHKERRQ(ierr); 510 snes->data = (void*) ngmres; 511 ngmres->msize = 30; 512 513 if (!snes->tolerancesset) { 514 snes->max_funcs = 30000; 515 snes->max_its = 10000; 516 } 517 518 ngmres->candidate = PETSC_FALSE; 519 520 ngmres->additive_linesearch = NULL; 521 ngmres->approxfunc = PETSC_FALSE; 522 ngmres->restart_it = 2; 523 ngmres->restart_periodic = 30; 524 ngmres->gammaA = 2.0; 525 ngmres->gammaC = 2.0; 526 ngmres->deltaB = 0.9; 527 ngmres->epsilonB = 0.1; 528 529 ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE; 530 ngmres->select_type = SNES_NGMRES_SELECT_DIFFERENCE; 531 532 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);CHKERRQ(ierr); 533 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);CHKERRQ(ierr); 534 PetscFunctionReturn(0); 535 } 536 537