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