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", "Maximum iterations before restart.", "SNES", ngmres->k_rmax, &ngmres->k_rmax, 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, " Maximum iterations before total restart: %d\n", ngmres->k_rmax);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; 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 /* no preconditioner -- just take gradient descent with line search */ 250 ierr = VecCopy(F, Y);CHKERRQ(ierr); 251 ierr = VecScale(Y, -1.0);CHKERRQ(ierr); 252 ierr = (*snes->ops->linesearch)(snes,snes->lsP,X,F,Y,fnorm,xnorm,G,W,&ynorm,&gnorm,&lssucceed);CHKERRQ(ierr); 253 if (!lssucceed) { 254 if (++snes->numFailures >= snes->maxFailures) { 255 snes->reason = SNES_DIVERGED_LINE_SEARCH; 256 PetscFunctionReturn(0); 257 } 258 } 259 fnorm = gnorm; 260 ierr = VecCopy(G, F);CHKERRQ(ierr); 261 ierr = VecCopy(W, X);CHKERRQ(ierr); 262 } else { 263 ierr = SNESSolve(snes->pc, B, X);CHKERRQ(ierr); 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 = SNESComputeFunction(snes, X, F);CHKERRQ(ierr); 270 ierr = VecNorm(F, NORM_2, &fnorm);CHKERRQ(ierr); 271 } 272 273 /* r = F(x) */ 274 nu = fnorm*fnorm; 275 if (fminnorm > fnorm) fminnorm = fnorm; /* the minimum norm is now of F^M */ 276 277 /* construct the right hand side and xi factors */ 278 for (i = 0; i < l; i++) { 279 ierr = VecDot(Fdot[i], F, &xi[i]);CHKERRQ(ierr); 280 beta[i] = nu - xi[i]; 281 } 282 283 /* construct h */ 284 for (j = 0; j < l; j++) { 285 for (i = 0; i < l; i++) { 286 H(i, j) = Q(i, j) - xi[i] - xi[j] + nu; 287 } 288 } 289 290 if(l == 1) { 291 /* simply set alpha[0] = beta[0] / H[0, 0] */ 292 beta[0] = beta[0] / H(0, 0); 293 } else { 294 #ifdef PETSC_MISSING_LAPACK_GELSS 295 SETERRQ(PETSC_COMM_SELF, PETSC_ERR_SUP, "NGMRES with LS requires the LAPACK GELSS routine."); 296 #else 297 ngmres->m = PetscBLASIntCast(l); 298 ngmres->n = PetscBLASIntCast(l); 299 ngmres->info = PetscBLASIntCast(0); 300 ngmres->rcond = -1.; 301 #ifdef PETSC_USE_COMPLEX 302 LAPACKgelss_(&ngmres->m, 303 &ngmres->n, 304 &ngmres->nrhs, 305 ngmres->h, 306 &ngmres->lda, 307 ngmres->beta, 308 &ngmres->ldb, 309 ngmres->s, 310 &ngmres->rcond, 311 &ngmres->rank, 312 ngmres->work, 313 &ngmres->lwork, 314 ngmres->rwork, 315 &ngmres->info); 316 #else 317 LAPACKgelss_(&ngmres->m, 318 &ngmres->n, 319 &ngmres->nrhs, 320 ngmres->h, 321 &ngmres->lda, 322 ngmres->beta, 323 &ngmres->ldb, 324 ngmres->s, 325 &ngmres->rcond, 326 &ngmres->rank, 327 ngmres->work, 328 &ngmres->lwork, 329 &ngmres->info); 330 #endif 331 if (ngmres->info < 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"Bad argument to GELSS"); 332 if (ngmres->info > 0) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_LIB,"SVD failed to converge"); 333 #endif 334 } 335 336 alph_total = 0.; 337 for (i = 0; i < l; i++) { 338 alph_total += beta[i]; 339 } 340 341 ierr = VecCopy(X, XA);CHKERRQ(ierr); 342 ierr = VecScale(XA, 1. - alph_total);CHKERRQ(ierr); 343 344 for(i=0;i<l;i++){ 345 ierr= VecAXPY(XA, beta[i], Xdot[i]);CHKERRQ(ierr); 346 } 347 ierr = SNESComputeFunction(snes, XA, FA);CHKERRQ(ierr); 348 ierr = VecNorm(FA, NORM_2, &fAnorm);CHKERRQ(ierr); 349 350 selectA = PETSC_TRUE; 351 /* Conditions for choosing the accelerated answer */ 352 353 /* Criterion A -- the norm of the function isn't increased above the minimum by too much */ 354 if (fAnorm >= ngmres->gammaA*fminnorm) { 355 selectA = PETSC_FALSE; 356 } 357 358 /* Criterion B -- the choice of x^A isn't too close to some other choice */ 359 ierr=VecCopy(XA, D);CHKERRQ(ierr); 360 ierr=VecAXPY(D, -1.0, X);CHKERRQ(ierr); 361 ierr=VecNorm(D, NORM_2, &dnorm);CHKERRQ(ierr); 362 dminnorm = -1.0; 363 for(i=0;i<l;i++) { 364 ierr=VecCopy(XA, D);CHKERRQ(ierr); 365 ierr=VecAXPY(D, -1.0, Xdot[i]);CHKERRQ(ierr); 366 ierr=VecNorm(D, NORM_2, &dcurnorm);CHKERRQ(ierr); 367 if((dcurnorm < dminnorm) || (dminnorm < 0.0)) dminnorm = dcurnorm; 368 } 369 if (ngmres->epsilonB*dnorm<dminnorm || sqrt(fnorm)<ngmres->deltaB*sqrt(fminnorm)) { 370 } else { 371 selectA=PETSC_FALSE; 372 } 373 374 375 if (selectA) { 376 if (ngmres->monitor) { 377 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "picked X_A, ||F_A||_2 = %e, ||F_M||_2 = %e\n", fAnorm, fnorm);CHKERRQ(ierr); 378 } 379 /* copy it over */ 380 fnorm = fAnorm; 381 nu = fnorm*fnorm; 382 ierr = VecCopy(FA, F);CHKERRQ(ierr); 383 ierr = VecCopy(XA, X);CHKERRQ(ierr); 384 } else { 385 if (ngmres->monitor) { 386 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "picked X_M, ||F_A||_2 = %e, ||F_M||_2 = %e\n", fAnorm, fnorm);CHKERRQ(ierr); 387 } 388 } 389 390 selectRestart = PETSC_FALSE; 391 392 /* maximum iteration criterion */ 393 if (k_restart > ngmres->k_rmax) { 394 selectRestart = PETSC_TRUE; 395 } 396 397 /* difference stagnation restart */ 398 if((ngmres->epsilonB*dnorm > dminnorm) && (sqrt(fAnorm) > ngmres->deltaB*sqrt(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 (sqrt(fAnorm) > ngmres->gammaC*sqrt(fminnorm)) { 406 if (ngmres->monitor) { 407 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "residual restart: %e > %e\n", sqrt(fAnorm), ngmres->gammaC*sqrt(fminnorm));CHKERRQ(ierr); 408 } 409 selectRestart = PETSC_TRUE; 410 } 411 412 if (selectRestart) { 413 if (ngmres->monitor){ 414 ierr = PetscViewerASCIIPrintf(ngmres->monitor, "Restarted at iteration %d\n", k_restart);CHKERRQ(ierr); 415 } 416 k_restart = 1; 417 l = 1; 418 /* q_{00} = nu */ 419 ngmres->fnorms[0] = fnorm; 420 nu = fnorm*fnorm; 421 Q(0,0) = nu; 422 /* Fdot[0] = F */ 423 ierr = VecCopy(X, Xdot[0]);CHKERRQ(ierr); 424 ierr = VecCopy(F, Fdot[0]);CHKERRQ(ierr); 425 } else { 426 /* select the current size of the subspace */ 427 if (l < ngmres->msize) l++; 428 k_restart++; 429 /* place the current entry in the list of previous entries */ 430 ierr = VecCopy(F, Fdot[ivec]);CHKERRQ(ierr); 431 ierr = VecCopy(X, Xdot[ivec]);CHKERRQ(ierr); 432 ngmres->fnorms[ivec] = fnorm; 433 if (fminnorm > fnorm) fminnorm = fnorm; /* the minimum norm is now of r^A */ 434 for (i = 0; i < l; i++) { 435 ierr = VecDot(F, Fdot[i], &qentry);CHKERRQ(ierr); 436 Q(i, ivec) = qentry; 437 Q(ivec, i) = qentry; 438 } 439 } 440 441 ierr = PetscObjectTakeAccess(snes);CHKERRQ(ierr); 442 snes->iter = k; 443 snes->norm = fnorm; 444 ierr = PetscObjectGrantAccess(snes);CHKERRQ(ierr); 445 SNESLogConvHistory(snes, snes->norm, snes->iter); 446 ierr = SNESMonitor(snes, snes->iter, snes->norm);CHKERRQ(ierr); 447 448 ierr = (*snes->ops->converged)(snes,snes->iter,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 449 if (snes->reason) PetscFunctionReturn(0); 450 } 451 snes->reason = SNES_DIVERGED_MAX_IT; 452 PetscFunctionReturn(0); 453 } 454 455 /*MC 456 SNESNGMRES - The Nonlinear Generalized Minimum Residual (NGMRES) method of Oosterlee and Washio. 457 458 Level: beginner 459 460 "Krylov Subspace Acceleration of Nonlinear Multigrid with Application to Recirculating Flows", C. W. Oosterlee and T. Washio, 461 SIAM Journal on Scientific Computing, 21(5), 2000. 462 463 This is also the same as the algorithm called Anderson acceleration introduced in "D. G. Anderson. Iterative procedures for nonlinear integral equations. 464 J. Assoc. Comput. Mach., 12:547–560, 1965." 465 466 .seealso: SNESCreate(), SNES, SNESSetType(), SNESType (for list of available types) 467 M*/ 468 469 EXTERN_C_BEGIN 470 #undef __FUNCT__ 471 #define __FUNCT__ "SNESCreate_NGMRES" 472 PetscErrorCode SNESCreate_NGMRES(SNES snes) 473 { 474 SNES_NGMRES *ngmres; 475 PetscErrorCode ierr; 476 477 PetscFunctionBegin; 478 snes->ops->destroy = SNESDestroy_NGMRES; 479 snes->ops->setup = SNESSetUp_NGMRES; 480 snes->ops->setfromoptions = SNESSetFromOptions_NGMRES; 481 snes->ops->view = SNESView_NGMRES; 482 snes->ops->solve = SNESSolve_NGMRES; 483 snes->ops->reset = SNESReset_NGMRES; 484 485 snes->usespc = PETSC_TRUE; 486 snes->usesksp = PETSC_FALSE; 487 488 ierr = PetscNewLog(snes, SNES_NGMRES, &ngmres);CHKERRQ(ierr); 489 snes->data = (void*) ngmres; 490 ngmres->msize = 10; 491 492 ngmres->gammaA = 2.0; 493 ngmres->gammaC = 2.0; 494 ngmres->deltaB = 0.9; 495 ngmres->epsilonB = 0.1; 496 ngmres->k_rmax = 200; 497 ierr = PetscObjectComposeFunctionDynamic((PetscObject)snes,"SNESLineSearchSetType_C","SNESLineSearchSetType_NGMRES",SNESLineSearchSetType_NGMRES);CHKERRQ(ierr); 498 ierr = SNESLineSearchSetType(snes, SNES_LS_QUADRATIC);CHKERRQ(ierr); 499 PetscFunctionReturn(0); 500 } 501 EXTERN_C_END 502