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