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 ierr = SNESSetWorkVecs(snes,5);CHKERRQ(ierr); 55 if (!ngmres->Xdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Xdot);CHKERRQ(ierr);} 56 if (!ngmres->Fdot) {ierr = VecDuplicateVecs(snes->vec_sol,ngmres->msize,&ngmres->Fdot);CHKERRQ(ierr);} 57 if (!ngmres->setup_called) { 58 msize = ngmres->msize; /* restart size */ 59 hsize = msize * msize; 60 61 /* explicit least squares minimization solve */ 62 ierr = PetscMalloc5(hsize,PetscScalar,&ngmres->h, 63 msize,PetscScalar,&ngmres->beta, 64 msize,PetscScalar,&ngmres->xi, 65 msize,PetscReal, &ngmres->fnorms, 66 hsize,PetscScalar,&ngmres->q);CHKERRQ(ierr); 67 if (ngmres->singlereduction) { 68 ierr = PetscMalloc(msize*sizeof(PetscReal),&ngmres->xnorms);CHKERRQ(ierr); 69 } 70 ngmres->nrhs = 1; 71 ngmres->lda = msize; 72 ngmres->ldb = msize; 73 ierr = PetscMalloc(msize*sizeof(PetscScalar),&ngmres->s);CHKERRQ(ierr); 74 ierr = PetscMemzero(ngmres->h, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 75 ierr = PetscMemzero(ngmres->q, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 76 ierr = PetscMemzero(ngmres->xi, msize*sizeof(PetscScalar));CHKERRQ(ierr); 77 ierr = PetscMemzero(ngmres->beta,msize*sizeof(PetscScalar));CHKERRQ(ierr); 78 ngmres->lwork = 12*msize; 79 #if PETSC_USE_COMPLEX 80 ierr = PetscMalloc(sizeof(PetscReal)*ngmres->lwork,&ngmres->rwork);CHKERRQ(ierr); 81 #endif 82 ierr = PetscMalloc(sizeof(PetscScalar)*ngmres->lwork,&ngmres->work);CHKERRQ(ierr); 83 } 84 85 /* linesearch setup */ 86 ierr = SNESGetOptionsPrefix(snes,&optionsprefix);CHKERRQ(ierr); 87 88 if (ngmres->select_type == SNES_NGMRES_SELECT_LINESEARCH) { 89 ierr = SNESLineSearchCreate(PetscObjectComm((PetscObject)snes),&ngmres->additive_linesearch);CHKERRQ(ierr); 90 ierr = SNESLineSearchSetSNES(ngmres->additive_linesearch,snes);CHKERRQ(ierr); 91 ierr = SNESLineSearchSetType(ngmres->additive_linesearch,SNESLINESEARCHL2);CHKERRQ(ierr); 92 ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,"additive_");CHKERRQ(ierr); 93 ierr = SNESLineSearchAppendOptionsPrefix(ngmres->additive_linesearch,optionsprefix);CHKERRQ(ierr); 94 ierr = SNESLineSearchSetFromOptions(ngmres->additive_linesearch);CHKERRQ(ierr); 95 } 96 97 ngmres->setup_called = PETSC_TRUE; 98 PetscFunctionReturn(0); 99 } 100 101 #undef __FUNCT__ 102 #define __FUNCT__ "SNESSetFromOptions_NGMRES" 103 PetscErrorCode SNESSetFromOptions_NGMRES(SNES snes) 104 { 105 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 106 PetscErrorCode ierr; 107 PetscBool debug; 108 SNESLineSearch linesearch; 109 110 PetscFunctionBegin; 111 ierr = PetscOptionsHead("SNES NGMRES options");CHKERRQ(ierr); 112 ierr = PetscOptionsEnum("-snes_ngmres_select_type","Select type","SNESNGMRESSetSelectType",SNESNGMRESSelectTypes, 113 (PetscEnum)ngmres->select_type,(PetscEnum*)&ngmres->select_type,NULL);CHKERRQ(ierr); 114 ierr = PetscOptionsEnum("-snes_ngmres_restart_type","Restart type","SNESNGMRESSetRestartType",SNESNGMRESRestartTypes, 115 (PetscEnum)ngmres->restart_type,(PetscEnum*)&ngmres->restart_type,NULL);CHKERRQ(ierr); 116 ierr = PetscOptionsBool("-snes_ngmres_candidate", "Use candidate storage", "SNES",ngmres->candidate,&ngmres->candidate,NULL);CHKERRQ(ierr); 117 ierr = PetscOptionsBool("-snes_ngmres_approxfunc","Linearly approximate the function", "SNES",ngmres->approxfunc,&ngmres->approxfunc,NULL);CHKERRQ(ierr); 118 ierr = PetscOptionsInt("-snes_ngmres_m", "Number of directions", "SNES",ngmres->msize,&ngmres->msize,NULL);CHKERRQ(ierr); 119 ierr = PetscOptionsInt("-snes_ngmres_restart", "Iterations before forced restart", "SNES",ngmres->restart_periodic,&ngmres->restart_periodic,NULL);CHKERRQ(ierr); 120 ierr = PetscOptionsInt("-snes_ngmres_restart_it", "Tolerance iterations before restart","SNES",ngmres->restart_it,&ngmres->restart_it,NULL);CHKERRQ(ierr); 121 ierr = PetscOptionsBool("-snes_ngmres_monitor", "Monitor actions of NGMRES", "SNES",ngmres->monitor ? PETSC_TRUE : PETSC_FALSE,&debug,NULL);CHKERRQ(ierr); 122 if (debug) { 123 ngmres->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));CHKERRQ(ierr); 124 } 125 ierr = PetscOptionsReal("-snes_ngmres_gammaA", "Residual selection constant", "SNES",ngmres->gammaA,&ngmres->gammaA,NULL);CHKERRQ(ierr); 126 ierr = PetscOptionsReal("-snes_ngmres_gammaC", "Residual restart constant", "SNES",ngmres->gammaC,&ngmres->gammaC,NULL);CHKERRQ(ierr); 127 ierr = PetscOptionsReal("-snes_ngmres_epsilonB", "Difference selection constant", "SNES",ngmres->epsilonB,&ngmres->epsilonB,NULL);CHKERRQ(ierr); 128 ierr = PetscOptionsReal("-snes_ngmres_deltaB", "Difference residual selection constant", "SNES",ngmres->deltaB,&ngmres->deltaB,NULL);CHKERRQ(ierr); 129 ierr = PetscOptionsBool("-snes_ngmres_single_reduction", "Aggregate reductions", "SNES",ngmres->singlereduction,&ngmres->singlereduction,NULL);CHKERRQ(ierr); 130 ierr = PetscOptionsTail();CHKERRQ(ierr); 131 if ((ngmres->gammaA > ngmres->gammaC) && (ngmres->gammaC > 2.)) ngmres->gammaC = ngmres->gammaA; 132 133 /* set the default type of the line search if the user hasn't already. */ 134 if (!snes->linesearch) { 135 ierr = SNESGetLineSearch(snes,&linesearch);CHKERRQ(ierr); 136 ierr = SNESLineSearchSetType(linesearch,SNESLINESEARCHBASIC);CHKERRQ(ierr); 137 } 138 PetscFunctionReturn(0); 139 } 140 141 #undef __FUNCT__ 142 #define __FUNCT__ "SNESView_NGMRES" 143 PetscErrorCode SNESView_NGMRES(SNES snes,PetscViewer viewer) 144 { 145 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 146 PetscBool iascii; 147 PetscErrorCode ierr; 148 149 PetscFunctionBegin; 150 ierr = PetscObjectTypeCompare((PetscObject) viewer,PETSCVIEWERASCII,&iascii);CHKERRQ(ierr); 151 if (iascii) { 152 ierr = PetscViewerASCIIPrintf(viewer," Number of stored past updates: %d\n", ngmres->msize);CHKERRQ(ierr); 153 ierr = PetscViewerASCIIPrintf(viewer," Residual selection: gammaA=%1.0e, gammaC=%1.0e\n",ngmres->gammaA,ngmres->gammaC);CHKERRQ(ierr); 154 ierr = PetscViewerASCIIPrintf(viewer," Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n",ngmres->epsilonB,ngmres->deltaB);CHKERRQ(ierr); 155 } 156 PetscFunctionReturn(0); 157 } 158 159 #undef __FUNCT__ 160 #define __FUNCT__ "SNESSolve_NGMRES" 161 PetscErrorCode SNESSolve_NGMRES(SNES snes) 162 { 163 164 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 165 /* present solution, residual, and preconditioned residual */ 166 Vec X,F,B,D,Y; 167 168 /* candidate linear combination answers */ 169 Vec XA,FA,XM,FM,FPC; 170 171 /* coefficients and RHS to the minimization problem */ 172 PetscReal fnorm,fMnorm,fAnorm; 173 PetscInt k,k_restart,l,ivec,restart_count = 0; 174 175 /* solution selection data */ 176 PetscBool selectRestart; 177 PetscReal dnorm,dminnorm = 0.0; 178 PetscReal fminnorm,xnorm,ynorm; 179 180 SNESConvergedReason reason; 181 PetscBool lssucceed; 182 PetscErrorCode ierr; 183 184 PetscFunctionBegin; 185 /* variable initialization */ 186 snes->reason = SNES_CONVERGED_ITERATING; 187 X = snes->vec_sol; 188 F = snes->vec_func; 189 B = snes->vec_rhs; 190 XA = snes->vec_sol_update; 191 FA = snes->work[0]; 192 D = snes->work[1]; 193 194 /* work for the line search */ 195 Y = snes->work[2]; 196 XM = snes->work[3]; 197 FM = snes->work[4]; 198 199 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 200 snes->iter = 0; 201 snes->norm = 0.; 202 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 203 204 /* initialization */ 205 206 if (snes->pc && snes->pcside == PC_LEFT) { 207 ierr = SNESApplyPC(snes,X,NULL,NULL,F);CHKERRQ(ierr); 208 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 209 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 210 snes->reason = SNES_DIVERGED_INNER; 211 PetscFunctionReturn(0); 212 } 213 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 214 } else { 215 if (!snes->vec_func_init_set) { 216 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 217 if (snes->domainerror) { 218 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 219 PetscFunctionReturn(0); 220 } 221 } else snes->vec_func_init_set = PETSC_FALSE; 222 if (!snes->norm_init_set) { 223 /* convergence test */ 224 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 225 if (PetscIsInfOrNanReal(fnorm)) { 226 snes->reason = SNES_DIVERGED_FNORM_NAN; 227 PetscFunctionReturn(0); 228 } 229 } else { 230 fnorm = snes->norm_init; 231 snes->norm_init_set = PETSC_FALSE; 232 } 233 } 234 fminnorm = fnorm; 235 236 /* q_{00} = nu */ 237 SNESNGMRESUpdateSubspace_Private(snes,0,0,F,fnorm,X); 238 239 ierr = PetscObjectAMSTakeAccess((PetscObject)snes);CHKERRQ(ierr); 240 snes->norm = fnorm; 241 ierr = PetscObjectAMSGrantAccess((PetscObject)snes);CHKERRQ(ierr); 242 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 243 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 244 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 245 if (snes->reason) PetscFunctionReturn(0); 246 247 k_restart = 1; 248 l = 1; 249 for (k=1; k < snes->max_its+1; k++) { 250 /* select which vector of the stored subspace will be updated */ 251 ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */ 252 253 /* Computation of x^M */ 254 if (snes->pc && snes->pcside == PC_RIGHT) { 255 ierr = VecCopy(X,XM);CHKERRQ(ierr); 256 ierr = SNESSetInitialFunction(snes->pc,F);CHKERRQ(ierr); 257 ierr = SNESSetInitialFunctionNorm(snes->pc,fnorm);CHKERRQ(ierr); 258 259 ierr = PetscLogEventBegin(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 260 ierr = SNESSolve(snes->pc,B,XM);CHKERRQ(ierr); 261 ierr = PetscLogEventEnd(SNES_NPCSolve,snes->pc,XM,B,0);CHKERRQ(ierr); 262 263 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 264 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 265 snes->reason = SNES_DIVERGED_INNER; 266 PetscFunctionReturn(0); 267 } 268 ierr = SNESGetFunction(snes->pc,&FPC,NULL,NULL);CHKERRQ(ierr); 269 ierr = VecCopy(FPC,FM);CHKERRQ(ierr); 270 ierr = SNESGetFunctionNorm(snes->pc,&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 snes->functype = SNES_FUNCTION_PRECONDITIONED; 509 510 ierr = PetscNewLog(snes,SNES_NGMRES,&ngmres);CHKERRQ(ierr); 511 snes->data = (void*) ngmres; 512 ngmres->msize = 30; 513 514 if (!snes->tolerancesset) { 515 snes->max_funcs = 30000; 516 snes->max_its = 10000; 517 } 518 519 ngmres->candidate = PETSC_FALSE; 520 521 ngmres->additive_linesearch = NULL; 522 ngmres->approxfunc = PETSC_FALSE; 523 ngmres->restart_it = 2; 524 ngmres->restart_periodic = 30; 525 ngmres->gammaA = 2.0; 526 ngmres->gammaC = 2.0; 527 ngmres->deltaB = 0.9; 528 ngmres->epsilonB = 0.1; 529 530 ngmres->restart_type = SNES_NGMRES_RESTART_DIFFERENCE; 531 ngmres->select_type = SNES_NGMRES_SELECT_DIFFERENCE; 532 533 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetSelectType_C",SNESNGMRESSetSelectType_NGMRES);CHKERRQ(ierr); 534 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNGMRESSetRestartType_C",SNESNGMRESSetRestartType_NGMRES);CHKERRQ(ierr); 535 PetscFunctionReturn(0); 536 } 537 538