1 /* Defines the basic SNES object */ 2 #include <../src/snes/impls/ngmres/snesngmres.h> 3 #include <petscblaslapack.h> 4 5 6 7 8 #undef __FUNCT__ 9 #define __FUNCT__ "SNESReset_NGMRES" 10 PetscErrorCode SNESReset_NGMRES(SNES snes) 11 { 12 SNES_NGMRES *ngmres = (SNES_NGMRES*) snes->data; 13 PetscErrorCode ierr; 14 15 PetscFunctionBegin; 16 ierr = VecDestroyVecs(ngmres->msize, &ngmres->Fdot);CHKERRQ(ierr); 17 ierr = VecDestroyVecs(ngmres->msize, &ngmres->Xdot);CHKERRQ(ierr); 18 if (snes->work) {ierr = VecDestroyVecs(snes->nwork, &snes->work);CHKERRQ(ierr);} 19 PetscFunctionReturn(0); 20 } 21 22 #undef __FUNCT__ 23 #define __FUNCT__ "SNESDestroy_NGMRES" 24 PetscErrorCode SNESDestroy_NGMRES(SNES snes) 25 { 26 PetscErrorCode ierr; 27 28 PetscFunctionBegin; 29 ierr = SNESReset_NGMRES(snes);CHKERRQ(ierr); 30 if (snes->data) { 31 SNES_NGMRES * ngmres = (SNES_NGMRES *)snes->data; 32 ierr = PetscFree5(ngmres->h, ngmres->beta, ngmres->xi, ngmres->fnorms, ngmres->q);CHKERRQ(ierr); 33 ierr = PetscFree(ngmres->s);CHKERRQ(ierr); 34 #if PETSC_USE_COMPLEX 35 ierr = PetscFree(ngmres->rwork); 36 #endif 37 ierr = PetscFree(ngmres->work); 38 } 39 ierr = PetscFree(snes->data); 40 PetscFunctionReturn(0); 41 } 42 43 #undef __FUNCT__ 44 #define __FUNCT__ "SNESSetUp_NGMRES" 45 PetscErrorCode SNESSetUp_NGMRES(SNES snes) 46 { 47 SNES_NGMRES *ngmres = (SNES_NGMRES *) snes->data; 48 PetscInt msize,hsize; 49 PetscErrorCode ierr; 50 51 PetscFunctionBegin; 52 msize = ngmres->msize; /* restart size */ 53 hsize = msize * msize; 54 55 56 /* explicit least squares minimization solve */ 57 ierr = PetscMalloc5(hsize,PetscScalar,&ngmres->h, 58 msize,PetscScalar,&ngmres->beta, 59 msize,PetscScalar,&ngmres->xi, 60 msize,PetscReal, &ngmres->fnorms, 61 hsize,PetscScalar,&ngmres->q);CHKERRQ(ierr); 62 ngmres->nrhs = 1; 63 ngmres->lda = msize; 64 ngmres->ldb = msize; 65 ierr = PetscMalloc(msize*sizeof(PetscScalar),&ngmres->s);CHKERRQ(ierr); 66 ierr = PetscMemzero(ngmres->h, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 67 ierr = PetscMemzero(ngmres->q, hsize*sizeof(PetscScalar));CHKERRQ(ierr); 68 ierr = PetscMemzero(ngmres->xi, msize*sizeof(PetscScalar));CHKERRQ(ierr); 69 ierr = PetscMemzero(ngmres->beta,msize*sizeof(PetscScalar));CHKERRQ(ierr); 70 ngmres->lwork = 12*msize; 71 #if PETSC_USE_COMPLEX 72 ierr = PetscMalloc(sizeof(PetscReal)*ngmres->lwork,&ngmres->rwork); 73 #endif 74 ierr = PetscMalloc(sizeof(PetscScalar)*ngmres->lwork,&ngmres->work); 75 76 ierr = VecDuplicateVecs(snes->vec_sol, ngmres->msize, &ngmres->Xdot);CHKERRQ(ierr); 77 ierr = VecDuplicateVecs(snes->vec_sol, ngmres->msize, &ngmres->Fdot);CHKERRQ(ierr); 78 ierr = SNESDefaultGetWork(snes, 5);CHKERRQ(ierr); 79 PetscFunctionReturn(0); 80 } 81 82 #undef __FUNCT__ 83 #define __FUNCT__ "SNESSetFromOptions_NGMRES" 84 PetscErrorCode SNESSetFromOptions_NGMRES(SNES snes) 85 { 86 SNES_NGMRES *ngmres = (SNES_NGMRES *) snes->data; 87 PetscErrorCode ierr; 88 PetscBool debug; 89 90 PetscFunctionBegin; 91 ierr = PetscOptionsHead("SNES NGMRES options");CHKERRQ(ierr); 92 ierr = PetscOptionsInt("-snes_ngmres_m", "Number of directions", "SNES", ngmres->msize, &ngmres->msize, PETSC_NULL);CHKERRQ(ierr); 93 ierr = PetscOptionsInt("-snes_ngmres_restart_it","Tolerance iterations before restart","SNES", ngmres->restart_it, &ngmres->restart_it, PETSC_NULL);CHKERRQ(ierr); 94 ierr = PetscOptionsBool("-snes_ngmres_monitor", "Monitor actions of NGMRES", "SNES", ngmres->monitor ? PETSC_TRUE: PETSC_FALSE, &debug, PETSC_NULL);CHKERRQ(ierr); 95 if (debug) { 96 ngmres->monitor = PETSC_VIEWER_STDOUT_(((PetscObject)snes)->comm);CHKERRQ(ierr); 97 } 98 ierr = PetscOptionsBool("-snes_ngmres_ngs", "Use nonlinear Gauss-Seidel", "SNES", ngmres->useGS, &ngmres->useGS, PETSC_NULL);CHKERRQ(ierr); 99 ierr = PetscOptionsReal("-snes_ngmres_gammaA", "Residual selection constant", "SNES", ngmres->gammaA, &ngmres->gammaA, PETSC_NULL);CHKERRQ(ierr); 100 ierr = PetscOptionsReal("-snes_ngmres_gammaC", " Residual restart constant", "SNES", ngmres->gammaC, &ngmres->gammaC, PETSC_NULL);CHKERRQ(ierr); 101 ierr = PetscOptionsReal("-snes_ngmres_epsilonB", "Difference selection constant", "SNES", ngmres->epsilonB, &ngmres->epsilonB, PETSC_NULL);CHKERRQ(ierr); 102 ierr = PetscOptionsReal("-snes_ngmres_deltaB", "Difference residual selection constant", "SNES", ngmres->deltaB, &ngmres->deltaB, PETSC_NULL);CHKERRQ(ierr); 103 ierr = PetscOptionsTail();CHKERRQ(ierr); 104 if ((ngmres->gammaA > ngmres->gammaC) && (ngmres->gammaC > 2.)) ngmres->gammaC = ngmres->gammaA; 105 PetscFunctionReturn(0); 106 } 107 108 #undef __FUNCT__ 109 #define __FUNCT__ "SNESView_NGMRES" 110 PetscErrorCode SNESView_NGMRES(SNES snes, PetscViewer viewer) 111 { 112 SNES_NGMRES *ngmres = (SNES_NGMRES *) snes->data; 113 PetscBool iascii; 114 const char *cstr; 115 PetscErrorCode ierr; 116 117 PetscFunctionBegin; 118 ierr = PetscTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 119 if (iascii) { 120 cstr = SNESLineSearchTypeName(snes->ls_type); 121 ierr = PetscViewerASCIIPrintf(viewer, " line search variant: %s\n",cstr);CHKERRQ(ierr); 122 ierr = PetscViewerASCIIPrintf(viewer, " Number of stored past updates: %d\n", ngmres->msize);CHKERRQ(ierr); 123 ierr = PetscViewerASCIIPrintf(viewer, " Residual selection: gammaA=%1.0e, gammaC=%1.0e\n", ngmres->gammaA, ngmres->gammaC);CHKERRQ(ierr); 124 ierr = PetscViewerASCIIPrintf(viewer, " Difference restart: epsilonB=%1.0e, deltaB=%1.0e\n", ngmres->epsilonB, ngmres->deltaB);CHKERRQ(ierr); 125 } 126 PetscFunctionReturn(0); 127 } 128 129 EXTERN_C_BEGIN 130 #undef __FUNCT__ 131 #define __FUNCT__ "SNESLineSearchSetType_NGMRES" 132 PetscErrorCode SNESLineSearchSetType_NGMRES(SNES snes, SNESLineSearchType type) 133 { 134 PetscErrorCode ierr; 135 PetscFunctionBegin; 136 137 switch (type) { 138 case SNES_LS_BASIC: 139 ierr = SNESLineSearchSet(snes,SNESLineSearchNo,PETSC_NULL);CHKERRQ(ierr); 140 break; 141 case SNES_LS_BASIC_NONORMS: 142 ierr = SNESLineSearchSet(snes,SNESLineSearchNoNorms,PETSC_NULL);CHKERRQ(ierr); 143 break; 144 case SNES_LS_QUADRATIC: 145 ierr = SNESLineSearchSet(snes,SNESLineSearchQuadraticSecant,PETSC_NULL);CHKERRQ(ierr); 146 break; 147 default: 148 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP,"Unknown line search type."); 149 break; 150 } 151 snes->ls_type = type; 152 PetscFunctionReturn(0); 153 } 154 EXTERN_C_END 155 156 157 #undef __FUNCT__ 158 #define __FUNCT__ "SNESSolve_NGMRES" 159 160 PetscErrorCode SNESSolve_NGMRES(SNES snes) 161 { 162 SNES_NGMRES *ngmres = (SNES_NGMRES *) snes->data; 163 /* present solution, residual, and preconditioned residual */ 164 Vec X, F, B, D, G, W, Y; 165 166 /* candidate linear combination answers */ 167 Vec XA, FA; 168 169 /* previous iterations to construct the subspace */ 170 Vec *Fdot = ngmres->Fdot; 171 Vec *Xdot = ngmres->Xdot; 172 173 /* coefficients and RHS to the minimization problem */ 174 PetscScalar *beta = ngmres->beta; 175 PetscScalar *xi = ngmres->xi; 176 PetscReal fnorm, fAnorm, gnorm, ynorm, xnorm = 0.0; 177 PetscReal nu; 178 PetscScalar alph_total = 0.; 179 PetscScalar qentry; 180 PetscInt i, j, k, k_restart, l, ivec, restart_count = 0; 181 182 /* solution selection data */ 183 PetscBool selectA, selectRestart; 184 PetscReal dnorm, dminnorm, dcurnorm; 185 PetscReal fminnorm; 186 187 SNESConvergedReason reason; 188 PetscBool lssucceed; 189 PetscErrorCode ierr; 190 191 PetscFunctionBegin; 192 /* variable initialization */ 193 snes->reason = SNES_CONVERGED_ITERATING; 194 X = snes->vec_sol; 195 F = snes->vec_func; 196 B = snes->vec_rhs; 197 XA = snes->vec_sol_update; 198 FA = snes->work[0]; 199 D = snes->work[1]; 200 201 /* work for the line search */ 202 Y = snes->work[2]; 203 G = snes->work[3]; 204 W = snes->work[4]; 205 206 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 207 snes->iter = 0; 208 snes->norm = 0.; 209 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 210 211 /* initialization */ 212 213 /* r = F(x) */ 214 ierr = SNESComputeFunction(snes, X, F);CHKERRQ(ierr); 215 if (snes->domainerror) { 216 snes->reason = SNES_DIVERGED_FUNCTION_DOMAIN; 217 PetscFunctionReturn(0); 218 } 219 220 /* nu = (r, r) */ 221 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 222 fminnorm = fnorm; 223 nu = fnorm*fnorm; 224 if (PetscIsInfOrNanReal(fnorm)) SETERRQ(PETSC_COMM_SELF, PETSC_ERR_FP, "Infinite or not-a-number generated in function evaluation"); 225 226 /* q_{00} = nu */ 227 Q(0,0) = nu; 228 ngmres->fnorms[0] = fnorm; 229 /* Fdot[0] = F */ 230 ierr = VecCopy(X, Xdot[0]);CHKERRQ(ierr); 231 ierr = VecCopy(F, Fdot[0]);CHKERRQ(ierr); 232 233 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 234 snes->norm = fnorm; 235 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 236 SNESLogConvHistory(snes, fnorm, 0); 237 ierr = SNESMonitor(snes, 0, fnorm);CHKERRQ(ierr); 238 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 239 if (snes->reason) PetscFunctionReturn(0); 240 241 k_restart = 1; 242 l = 1; 243 for (k=1; k<snes->max_its; k++) { 244 /* select which vector of the stored subspace will be updated */ 245 ivec = k_restart % ngmres->msize; /* replace the last used part of the subspace */ 246 247 /* Computation of x^M */ 248 if (snes->pc) { 249 ierr = SNESSolve(snes->pc, B, X);CHKERRQ(ierr); 250 ierr = SNESGetConvergedReason(snes->pc,&reason);CHKERRQ(ierr); 251 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 252 snes->reason = SNES_DIVERGED_INNER; 253 PetscFunctionReturn(0); 254 } 255 ierr = SNESComputeFunction(snes, X, F);CHKERRQ(ierr); 256 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 257 } else if (ngmres->useGS && snes->ops->computegs) { 258 /* compute the update using the supplied Gauss-Seidel routine */ 259 ierr = SNESComputeGS(snes, B, X);CHKERRQ(ierr); 260 ierr = SNESComputeFunction(snes, X, F);CHKERRQ(ierr); 261 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 262 } else { 263 /* no preconditioner -- just take gradient descent with line search */ 264 ierr = VecCopy(F, Y);CHKERRQ(ierr); 265 ierr = VecScale(Y, -1.0);CHKERRQ(ierr); 266 ierr = (*snes->ops->linesearch)(snes,snes->lsP,X,F,Y,fnorm,xnorm,G,W,&ynorm,&gnorm,&lssucceed);CHKERRQ(ierr); 267 if (!lssucceed) { 268 if (++snes->numFailures >= snes->maxFailures) { 269 snes->reason = SNES_DIVERGED_LINE_SEARCH; 270 PetscFunctionReturn(0); 271 } 272 } 273 fnorm = gnorm; 274 ierr = VecCopy(G, F);CHKERRQ(ierr); 275 ierr = VecCopy(W, X);CHKERRQ(ierr); 276 } 277 278 /* r = F(x) */ 279 nu = fnorm*fnorm; 280 if (fminnorm > fnorm) fminnorm = fnorm; /* the minimum norm is now of F^M */ 281 282 /* construct the right hand side and xi factors */ 283 for (i = 0; i < l; i++) { 284 ierr = VecDot(Fdot[i], F, &xi[i]);CHKERRQ(ierr); 285 beta[i] = nu - xi[i]; 286 } 287 288 /* construct h */ 289 for (j = 0; j < l; j++) { 290 for (i = 0; i < l; i++) { 291 H(i, j) = Q(i, j) - xi[i] - xi[j] + nu; 292 } 293 } 294 295 if(l == 1) { 296 /* simply set alpha[0] = beta[0] / H[0, 0] */ 297 beta[0] = beta[0] / H(0, 0); 298 } else { 299 #ifdef PETSC_MISSING_LAPACK_GELSS 300 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "NGMRES with LS requires the LAPACK GELSS routine."); 301 #else 302 ngmres->m = PetscBLASIntCast(l); 303 ngmres->n = PetscBLASIntCast(l); 304 ngmres->info = PetscBLASIntCast(0); 305 ngmres->rcond = -1.; 306 #ifdef PETSC_USE_COMPLEX 307 LAPACKgelss_(&ngmres->m, 308 &ngmres->n, 309 &ngmres->nrhs, 310 ngmres->h, 311 &ngmres->lda, 312 ngmres->beta, 313 &ngmres->ldb, 314 ngmres->s, 315 &ngmres->rcond, 316 &ngmres->rank, 317 ngmres->work, 318 &ngmres->lwork, 319 ngmres->rwork, 320 &ngmres->info); 321 #else 322 LAPACKgelss_(&ngmres->m, 323 &ngmres->n, 324 &ngmres->nrhs, 325 ngmres->h, 326 &ngmres->lda, 327 ngmres->beta, 328 &ngmres->ldb, 329 ngmres->s, 330 &ngmres->rcond, 331 &ngmres->rank, 332 ngmres->work, 333 &ngmres->lwork, 334 &ngmres->info); 335 #endif 336 if (ngmres->info < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Bad argument to GELSS"); 337 if (ngmres->info > 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"SVD failed to converge"); 338 #endif 339 } 340 341 alph_total = 0.; 342 for (i = 0; i < l; i++) { 343 alph_total += beta[i]; 344 } 345 346 ierr = VecCopy(X, XA);CHKERRQ(ierr); 347 ierr = VecScale(XA, 1. - alph_total);CHKERRQ(ierr); 348 349 for(i=0;i<l;i++){ 350 ierr= VecAXPY(XA, beta[i], Xdot[i]);CHKERRQ(ierr); 351 } 352 ierr = SNESComputeFunction(snes, XA, FA);CHKERRQ(ierr); 353 ierr = VecNorm(FA, NORM_2, &fAnorm);CHKERRQ(ierr); 354 355 selectA = PETSC_TRUE; 356 /* Conditions for choosing the accelerated answer */ 357 358 /* Criterion A -- the norm of the function isn't increased above the minimum by too much */ 359 if (fAnorm >= ngmres->gammaA*fminnorm) { 360 selectA = PETSC_FALSE; 361 } 362 363 /* Criterion B -- the choice of x^A isn't too close to some other choice */ 364 ierr=VecCopy(XA, D);CHKERRQ(ierr); 365 ierr=VecAXPY(D, -1.0, X);CHKERRQ(ierr); 366 ierr=VecNorm(D, NORM_2, &dnorm);CHKERRQ(ierr); 367 dminnorm = -1.0; 368 for(i=0;i<l;i++) { 369 ierr=VecCopy(XA, D);CHKERRQ(ierr); 370 ierr=VecAXPY(D, -1.0, Xdot[i]);CHKERRQ(ierr); 371 ierr=VecNorm(D, NORM_2, &dcurnorm);CHKERRQ(ierr); 372 if((dcurnorm < dminnorm) || (dminnorm < 0.0)) dminnorm = dcurnorm; 373 } 374 if (ngmres->epsilonB*dnorm<dminnorm || PetscSqrtReal(fnorm)<ngmres->deltaB*PetscSqrtReal(fminnorm)) { 375 } else { 376 selectA=PETSC_FALSE; 377 } 378 379 380 if (selectA) { 381 if (ngmres->monitor) { 382 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "picked X_A, ||F_A||_2 = %e, ||F_M||_2 = %e\n", fAnorm, fnorm);CHKERRQ(ierr); 383 } 384 /* copy it over */ 385 fnorm = fAnorm; 386 nu = fnorm*fnorm; 387 ierr = VecCopy(FA, F);CHKERRQ(ierr); 388 ierr = VecCopy(XA, X);CHKERRQ(ierr); 389 } else { 390 if (ngmres->monitor) { 391 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "picked X_M, ||F_A||_2 = %e, ||F_M||_2 = %e\n", fAnorm, fnorm);CHKERRQ(ierr); 392 } 393 } 394 395 selectRestart = PETSC_FALSE; 396 397 /* difference stagnation restart */ 398 if((ngmres->epsilonB*dnorm > dminnorm) && (PetscSqrtReal(fAnorm) > ngmres->deltaB*PetscSqrtReal(fminnorm))) { 399 if (ngmres->monitor) { 400 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "difference restart: %e > %e\n", ngmres->epsilonB*dnorm, dminnorm);CHKERRQ(ierr); 401 } 402 selectRestart = PETSC_TRUE; 403 } 404 /* residual stagnation restart */ 405 if (PetscSqrtReal(fAnorm) > ngmres->gammaC*PetscSqrtReal(fminnorm)) { 406 if (ngmres->monitor) { 407 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "residual restart: %e > %e\n", PetscSqrtReal(fAnorm), ngmres->gammaC*PetscSqrtReal(fminnorm));CHKERRQ(ierr); 408 } 409 selectRestart = PETSC_TRUE; 410 } 411 412 /* if the restart conditions persist for more than restart_it iterations, restart. */ 413 if (selectRestart) { 414 restart_count++; 415 } else { 416 restart_count = 0; 417 } 418 419 /* restart after restart conditions have persisted for a fixed number of iterations */ 420 if (restart_count >= ngmres->restart_it) { 421 if (ngmres->monitor){ 422 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "Restarted at iteration %d\n", k_restart);CHKERRQ(ierr); 423 } 424 restart_count = 0; 425 k_restart = 1; 426 l = 1; 427 /* q_{00} = nu */ 428 ngmres->fnorms[0] = fnorm; 429 nu = fnorm*fnorm; 430 Q(0,0) = nu; 431 /* Fdot[0] = F */ 432 ierr = VecCopy(X, Xdot[0]);CHKERRQ(ierr); 433 ierr = VecCopy(F, Fdot[0]);CHKERRQ(ierr); 434 } else { 435 /* select the current size of the subspace */ 436 if (l < ngmres->msize) l++; 437 k_restart++; 438 /* place the current entry in the list of previous entries */ 439 ierr = VecCopy(F, Fdot[ivec]);CHKERRQ(ierr); 440 ierr = VecCopy(X, Xdot[ivec]);CHKERRQ(ierr); 441 ngmres->fnorms[ivec] = fnorm; 442 if (fminnorm > fnorm) fminnorm = fnorm; /* the minimum norm is now of r^A */ 443 for (i = 0; i < l; i++) { 444 ierr = VecDot(F, Fdot[i], &qentry);CHKERRQ(ierr); 445 Q(i, ivec) = qentry; 446 Q(ivec, i) = qentry; 447 } 448 } 449 450 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 451 snes->iter = k; 452 snes->norm = fnorm; 453 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 454 SNESLogConvHistory(snes, snes->norm, snes->iter); 455 ierr = SNESMonitor(snes, snes->iter, snes->norm);CHKERRQ(ierr); 456 457 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 458 if (snes->reason) PetscFunctionReturn(0); 459 } 460 snes->reason = SNES_DIVERGED_MAX_IT; 461 PetscFunctionReturn(0); 462 } 463 464 /*MC 465 SNESNGMRES - The Nonlinear Generalized Minimum Residual (NGMRES) method of Oosterlee and Washio. 466 467 Level: beginner 468 469 "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows", C. W. Oosterlee and T. Washio, 470 SIAM Journal on Scientific Computing, 21(5), 2000. 471 472 This is also the same as the algorithm called Anderson acceleration introduced in "D. G. Anderson. Iterative procedures for nonlinear integral equations. 473 J. Assoc. Comput. Mach., 12:547–560, 1965." 474 475 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 476 M*/ 477 478 EXTERN_C_BEGIN 479 #undef __FUNCT__ 480 #define __FUNCT__ "SNESCreate_NGMRES" 481 PetscErrorCode SNESCreate_NGMRES(SNES snes) 482 { 483 SNES_NGMRES *ngmres; 484 PetscErrorCode ierr; 485 486 PetscFunctionBegin; 487 snes->ops->destroy = SNESDestroy_NGMRES; 488 snes->ops->setup = SNESSetUp_NGMRES; 489 snes->ops->setfromoptions = SNESSetFromOptions_NGMRES; 490 snes->ops->view = SNESView_NGMRES; 491 snes->ops->solve = SNESSolve_NGMRES; 492 snes->ops->reset = SNESReset_NGMRES; 493 494 snes->usespc = PETSC_TRUE; 495 snes->usesksp = PETSC_FALSE; 496 497 ierr = PetscNewLog(snes, SNES_NGMRES, &ngmres);CHKERRQ(ierr); 498 snes->data = (void*) ngmres; 499 ngmres->msize = 10; 500 ngmres->useGS = PETSC_FALSE; 501 502 ngmres->restart_it = 2; 503 ngmres->gammaA = 2.0; 504 ngmres->gammaC = 2.0; 505 ngmres->deltaB = 0.9; 506 ngmres->epsilonB = 0.1; 507 ierr = PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetType_C","SNESLineSearchSetType_NGMRES",SNESLineSearchSetType_NGMRES);CHKERRQ(ierr); 508 ierr = SNESLineSearchSetType(snes, SNES_LS_QUADRATIC);CHKERRQ(ierr); 509 PetscFunctionReturn(0); 510 } 511 EXTERN_C_END 512