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