1 #include <../src/snes/impls/ncg/snesncgimpl.h> /*I "petscsnes.h" I*/ 2 const char *const SNESNCGTypes[] = {"FR","PRP","HS","DY","CD","SNESNCGType","SNES_NCG_",0}; 3 4 PetscErrorCode SNESReset_NCG(SNES snes) 5 { 6 PetscFunctionBegin; 7 PetscFunctionReturn(0); 8 } 9 10 #define SNESLINESEARCHNCGLINEAR "ncglinear" 11 12 /* 13 SNESDestroy_NCG - Destroys the private SNES_NCG context that was created with SNESCreate_NCG(). 14 15 Input Parameter: 16 . snes - the SNES context 17 18 Application Interface Routine: SNESDestroy() 19 */ 20 PetscErrorCode SNESDestroy_NCG(SNES snes) 21 { 22 PetscErrorCode ierr; 23 24 PetscFunctionBegin; 25 ierr = PetscFree(snes->data);CHKERRQ(ierr); 26 PetscFunctionReturn(0); 27 } 28 29 /* 30 SNESSetUp_NCG - Sets up the internal data structures for the later use 31 of the SNESNCG nonlinear solver. 32 33 Input Parameters: 34 + snes - the SNES context 35 - x - the solution vector 36 37 Application Interface Routine: SNESSetUp() 38 */ 39 40 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_NCGLinear(SNESLineSearch); 41 42 PetscErrorCode SNESSetUp_NCG(SNES snes) 43 { 44 PetscErrorCode ierr; 45 46 PetscFunctionBegin; 47 ierr = SNESSetWorkVecs(snes,2);CHKERRQ(ierr); 48 if (snes->npcside== PC_RIGHT) SETERRQ(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNESNCG only supports left preconditioning"); 49 if (snes->functype == SNES_FUNCTION_DEFAULT) snes->functype = SNES_FUNCTION_UNPRECONDITIONED; 50 PetscFunctionReturn(0); 51 } 52 /* 53 SNESSetFromOptions_NCG - Sets various parameters for the SNESNCG method. 54 55 Input Parameter: 56 . snes - the SNES context 57 58 Application Interface Routine: SNESSetFromOptions() 59 */ 60 static PetscErrorCode SNESSetFromOptions_NCG(PetscOptionItems *PetscOptionsObject,SNES snes) 61 { 62 SNES_NCG *ncg = (SNES_NCG*)snes->data; 63 PetscErrorCode ierr; 64 PetscBool debug = PETSC_FALSE; 65 SNESNCGType ncgtype=ncg->type; 66 SNESLineSearch linesearch; 67 68 PetscFunctionBegin; 69 ierr = PetscOptionsHead(PetscOptionsObject,"SNES NCG options");CHKERRQ(ierr); 70 ierr = PetscOptionsBool("-snes_ncg_monitor","Monitor NCG iterations","SNES",ncg->monitor ? PETSC_TRUE : PETSC_FALSE, &debug, NULL);CHKERRQ(ierr); 71 if (debug) { 72 ncg->monitor = PETSC_VIEWER_STDOUT_(PetscObjectComm((PetscObject)snes));CHKERRQ(ierr); 73 } 74 ierr = PetscOptionsEnum("-snes_ncg_type","NCG Beta type used","SNESNCGSetType",SNESNCGTypes,(PetscEnum)ncg->type,(PetscEnum*)&ncgtype,NULL);CHKERRQ(ierr); 75 ierr = SNESNCGSetType(snes, ncgtype);CHKERRQ(ierr); 76 ierr = PetscOptionsTail();CHKERRQ(ierr); 77 if (!snes->linesearch) { 78 ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr); 79 if (!snes->npc) { 80 ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHCP);CHKERRQ(ierr); 81 } else { 82 ierr = SNESLineSearchSetType(linesearch, SNESLINESEARCHL2);CHKERRQ(ierr); 83 } 84 } 85 ierr = SNESLineSearchRegister(SNESLINESEARCHNCGLINEAR, SNESLineSearchCreate_NCGLinear);CHKERRQ(ierr); 86 PetscFunctionReturn(0); 87 } 88 89 /* 90 SNESView_NCG - Prints info from the SNESNCG data structure. 91 92 Input Parameters: 93 + SNES - the SNES context 94 - viewer - visualization context 95 96 Application Interface Routine: SNESView() 97 */ 98 static PetscErrorCode SNESView_NCG(SNES snes, PetscViewer viewer) 99 { 100 SNES_NCG *ncg = (SNES_NCG *) snes->data; 101 PetscBool iascii; 102 PetscErrorCode ierr; 103 104 PetscFunctionBegin; 105 ierr = PetscObjectTypeCompare((PetscObject) viewer, PETSCVIEWERASCII, &iascii);CHKERRQ(ierr); 106 if (iascii) { 107 ierr = PetscViewerASCIIPrintf(viewer, " type: %s\n", SNESNCGTypes[ncg->type]);CHKERRQ(ierr); 108 } 109 PetscFunctionReturn(0); 110 } 111 112 PetscErrorCode SNESLineSearchApply_NCGLinear(SNESLineSearch linesearch) 113 { 114 PetscScalar alpha, ptAp; 115 Vec X, Y, F, W; 116 SNES snes; 117 PetscErrorCode ierr; 118 PetscReal *fnorm, *xnorm, *ynorm; 119 120 PetscFunctionBegin; 121 ierr = SNESLineSearchGetSNES(linesearch, &snes);CHKERRQ(ierr); 122 X = linesearch->vec_sol; 123 W = linesearch->vec_sol_new; 124 F = linesearch->vec_func; 125 Y = linesearch->vec_update; 126 fnorm = &linesearch->fnorm; 127 xnorm = &linesearch->xnorm; 128 ynorm = &linesearch->ynorm; 129 130 if (!snes->jacobian) { 131 ierr = SNESSetUpMatrices(snes);CHKERRQ(ierr); 132 } 133 134 /* 135 136 The exact step size for unpreconditioned linear CG is just: 137 alpha = (r, r) / (p, Ap) = (f, f) / (y, Jy) 138 */ 139 ierr = SNESComputeJacobian(snes, X, snes->jacobian, snes->jacobian_pre);CHKERRQ(ierr); 140 ierr = VecDot(F, F, &alpha);CHKERRQ(ierr); 141 ierr = MatMult(snes->jacobian, Y, W);CHKERRQ(ierr); 142 ierr = VecDot(Y, W, &ptAp);CHKERRQ(ierr); 143 alpha = alpha / ptAp; 144 ierr = VecAXPY(X,-alpha,Y);CHKERRQ(ierr); 145 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 146 147 ierr = VecNorm(F,NORM_2,fnorm);CHKERRQ(ierr); 148 ierr = VecNorm(X,NORM_2,xnorm);CHKERRQ(ierr); 149 ierr = VecNorm(Y,NORM_2,ynorm);CHKERRQ(ierr); 150 PetscFunctionReturn(0); 151 } 152 153 PETSC_EXTERN PetscErrorCode SNESLineSearchCreate_NCGLinear(SNESLineSearch linesearch) 154 { 155 PetscFunctionBegin; 156 linesearch->ops->apply = SNESLineSearchApply_NCGLinear; 157 linesearch->ops->destroy = NULL; 158 linesearch->ops->setfromoptions = NULL; 159 linesearch->ops->reset = NULL; 160 linesearch->ops->view = NULL; 161 linesearch->ops->setup = NULL; 162 PetscFunctionReturn(0); 163 } 164 165 /* 166 167 Assuming F = SNESComputeFunction(X) compute Y^tJ^tF using a simple secant approximation of the jacobian. 168 169 */ 170 PetscErrorCode SNESNCGComputeYtJtF_Private(SNES snes, Vec X, Vec F, Vec Y, Vec W, Vec G, PetscScalar * ytJtf) 171 { 172 PetscErrorCode ierr; 173 PetscScalar ftf, ftg, fty, h; 174 175 PetscFunctionBegin; 176 ierr = VecDot(F, F, &ftf);CHKERRQ(ierr); 177 ierr = VecDot(F, Y, &fty);CHKERRQ(ierr); 178 h = 1e-5*fty / fty; 179 ierr = VecCopy(X, W);CHKERRQ(ierr); 180 ierr = VecAXPY(W, -h, Y);CHKERRQ(ierr); /* this is arbitrary */ 181 ierr = SNESComputeFunction(snes, W, G);CHKERRQ(ierr); 182 ierr = VecDot(G, F, &ftg);CHKERRQ(ierr); 183 *ytJtf = (ftg - ftf) / h; 184 PetscFunctionReturn(0); 185 } 186 187 /*@ 188 SNESNCGSetType - Sets the conjugate update type for SNESNCG. 189 190 Logically Collective on SNES 191 192 Input Parameters: 193 + snes - the iterative context 194 - btype - update type 195 196 Options Database: 197 . -snes_ncg_type<prp,fr,hs,dy,cd> 198 199 Level: intermediate 200 201 SNESNCGSelectTypes: 202 + SNES_NCG_FR - Fletcher-Reeves update 203 . SNES_NCG_PRP - Polak-Ribiere-Polyak update 204 . SNES_NCG_HS - Hestenes-Steifel update 205 . SNES_NCG_DY - Dai-Yuan update 206 - SNES_NCG_CD - Conjugate Descent update 207 208 Notes: 209 PRP is the default, and the only one that tolerates generalized search directions. 210 211 @*/ 212 PetscErrorCode SNESNCGSetType(SNES snes, SNESNCGType btype) 213 { 214 PetscErrorCode ierr; 215 216 PetscFunctionBegin; 217 PetscValidHeaderSpecific(snes,SNES_CLASSID,1); 218 ierr = PetscTryMethod(snes,"SNESNCGSetType_C",(SNES,SNESNCGType),(snes,btype));CHKERRQ(ierr); 219 PetscFunctionReturn(0); 220 } 221 222 PetscErrorCode SNESNCGSetType_NCG(SNES snes, SNESNCGType btype) 223 { 224 SNES_NCG *ncg = (SNES_NCG*)snes->data; 225 226 PetscFunctionBegin; 227 ncg->type = btype; 228 PetscFunctionReturn(0); 229 } 230 231 /* 232 SNESSolve_NCG - Solves a nonlinear system with the Nonlinear Conjugate Gradient method. 233 234 Input Parameters: 235 . snes - the SNES context 236 237 Output Parameter: 238 . outits - number of iterations until termination 239 240 Application Interface Routine: SNESSolve() 241 */ 242 PetscErrorCode SNESSolve_NCG(SNES snes) 243 { 244 SNES_NCG *ncg = (SNES_NCG*)snes->data; 245 Vec X,dX,lX,F,dXold; 246 PetscReal fnorm, ynorm, xnorm, beta = 0.0; 247 PetscScalar dXdotdX, dXolddotdXold, dXdotdXold, lXdotdX, lXdotdXold; 248 PetscInt maxits, i; 249 PetscErrorCode ierr; 250 SNESLineSearchReason lsresult = SNES_LINESEARCH_SUCCEEDED; 251 SNESLineSearch linesearch; 252 SNESConvergedReason reason; 253 254 PetscFunctionBegin; 255 if (snes->xl || snes->xu || snes->ops->computevariablebounds) SETERRQ1(PetscObjectComm((PetscObject)snes),PETSC_ERR_ARG_WRONGSTATE, "SNES solver %s does not support bounds", ((PetscObject)snes)->type_name); 256 257 ierr = PetscCitationsRegister(SNESCitation,&SNEScite);CHKERRQ(ierr); 258 snes->reason = SNES_CONVERGED_ITERATING; 259 260 maxits = snes->max_its; /* maximum number of iterations */ 261 X = snes->vec_sol; /* X^n */ 262 dXold = snes->work[0]; /* The previous iterate of X */ 263 dX = snes->work[1]; /* the preconditioned direction */ 264 lX = snes->vec_sol_update; /* the conjugate direction */ 265 F = snes->vec_func; /* residual vector */ 266 267 ierr = SNESGetLineSearch(snes, &linesearch);CHKERRQ(ierr); 268 269 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 270 snes->iter = 0; 271 snes->norm = 0.; 272 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 273 274 /* compute the initial function and preconditioned update dX */ 275 276 if (snes->npc && snes->functype == SNES_FUNCTION_PRECONDITIONED) { 277 ierr = SNESApplyNPC(snes,X,NULL,dX);CHKERRQ(ierr); 278 ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr); 279 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 280 snes->reason = SNES_DIVERGED_INNER; 281 PetscFunctionReturn(0); 282 } 283 ierr = VecCopy(dX,F);CHKERRQ(ierr); 284 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 285 } else { 286 if (!snes->vec_func_init_set) { 287 ierr = SNESComputeFunction(snes,X,F);CHKERRQ(ierr); 288 } else snes->vec_func_init_set = PETSC_FALSE; 289 290 /* convergence test */ 291 ierr = VecNorm(F,NORM_2,&fnorm);CHKERRQ(ierr); 292 SNESCheckFunctionNorm(snes,fnorm); 293 ierr = VecCopy(F,dX);CHKERRQ(ierr); 294 } 295 if (snes->npc) { 296 if (snes->functype == SNES_FUNCTION_UNPRECONDITIONED) { 297 ierr = SNESApplyNPC(snes,X,F,dX);CHKERRQ(ierr); 298 ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr); 299 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 300 snes->reason = SNES_DIVERGED_INNER; 301 PetscFunctionReturn(0); 302 } 303 } 304 } 305 ierr = VecCopy(dX,lX);CHKERRQ(ierr); 306 ierr = VecDot(dX, dX, &dXdotdX);CHKERRQ(ierr); 307 308 309 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 310 snes->norm = fnorm; 311 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 312 ierr = SNESLogConvergenceHistory(snes,fnorm,0);CHKERRQ(ierr); 313 ierr = SNESMonitor(snes,0,fnorm);CHKERRQ(ierr); 314 315 /* test convergence */ 316 ierr = (*snes->ops->converged)(snes,0,0.0,0.0,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 317 if (snes->reason) PetscFunctionReturn(0); 318 319 /* Call general purpose update function */ 320 if (snes->ops->update) { 321 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 322 } 323 324 /* first update -- just use the (preconditioned) residual direction for the initial conjugate direction */ 325 326 for (i = 1; i < maxits + 1; i++) { 327 /* some update types require the old update direction or conjugate direction */ 328 if (ncg->type != SNES_NCG_FR) { 329 ierr = VecCopy(dX, dXold);CHKERRQ(ierr); 330 } 331 ierr = SNESLineSearchApply(linesearch,X,F,&fnorm,lX);CHKERRQ(ierr); 332 ierr = SNESLineSearchGetReason(linesearch, &lsresult);CHKERRQ(ierr); 333 ierr = SNESLineSearchGetNorms(linesearch, &xnorm, &fnorm, &ynorm);CHKERRQ(ierr); 334 if (lsresult) { 335 if (++snes->numFailures >= snes->maxFailures) { 336 snes->reason = SNES_DIVERGED_LINE_SEARCH; 337 PetscFunctionReturn(0); 338 } 339 } 340 if (snes->nfuncs >= snes->max_funcs && snes->max_funcs >= 0) { 341 snes->reason = SNES_DIVERGED_FUNCTION_COUNT; 342 PetscFunctionReturn(0); 343 } 344 /* Monitor convergence */ 345 ierr = PetscObjectSAWsTakeAccess((PetscObject)snes);CHKERRQ(ierr); 346 snes->iter = i; 347 snes->norm = fnorm; 348 snes->xnorm = xnorm; 349 snes->ynorm = ynorm; 350 ierr = PetscObjectSAWsGrantAccess((PetscObject)snes);CHKERRQ(ierr); 351 ierr = SNESLogConvergenceHistory(snes,snes->norm,0);CHKERRQ(ierr); 352 ierr = SNESMonitor(snes,snes->iter,snes->norm);CHKERRQ(ierr); 353 354 /* Test for convergence */ 355 ierr = (*snes->ops->converged)(snes,snes->iter,xnorm,ynorm,fnorm,&snes->reason,snes->cnvP);CHKERRQ(ierr); 356 if (snes->reason) PetscFunctionReturn(0); 357 358 /* Call general purpose update function */ 359 if (snes->ops->update) { 360 ierr = (*snes->ops->update)(snes, snes->iter);CHKERRQ(ierr); 361 } 362 if (snes->npc) { 363 if (snes->functype == SNES_FUNCTION_PRECONDITIONED) { 364 ierr = SNESApplyNPC(snes,X,NULL,dX);CHKERRQ(ierr); 365 ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr); 366 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 367 snes->reason = SNES_DIVERGED_INNER; 368 PetscFunctionReturn(0); 369 } 370 ierr = VecCopy(dX,F);CHKERRQ(ierr); 371 } else { 372 ierr = SNESApplyNPC(snes,X,F,dX);CHKERRQ(ierr); 373 ierr = SNESGetConvergedReason(snes->npc,&reason);CHKERRQ(ierr); 374 if (reason < 0 && reason != SNES_DIVERGED_MAX_IT) { 375 snes->reason = SNES_DIVERGED_INNER; 376 PetscFunctionReturn(0); 377 } 378 } 379 } else { 380 ierr = VecCopy(F,dX);CHKERRQ(ierr); 381 } 382 383 /* compute the conjugate direction lX = dX + beta*lX with beta = ((dX, dX) / (dX_old, dX_old) (Fletcher-Reeves update)*/ 384 switch (ncg->type) { 385 case SNES_NCG_FR: /* Fletcher-Reeves */ 386 dXolddotdXold= dXdotdX; 387 ierr = VecDot(dX, dX, &dXdotdX);CHKERRQ(ierr); 388 beta = PetscRealPart(dXdotdX / dXolddotdXold); 389 break; 390 case SNES_NCG_PRP: /* Polak-Ribiere-Poylak */ 391 dXolddotdXold= dXdotdX; 392 ierr = VecDotBegin(dX, dX, &dXdotdXold);CHKERRQ(ierr); 393 ierr = VecDotBegin(dXold, dX, &dXdotdXold);CHKERRQ(ierr); 394 ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr); 395 ierr = VecDotEnd(dXold, dX, &dXdotdXold);CHKERRQ(ierr); 396 beta = PetscRealPart(((dXdotdX - dXdotdXold) / dXolddotdXold)); 397 if (beta < 0.0) beta = 0.0; /* restart */ 398 break; 399 case SNES_NCG_HS: /* Hestenes-Stiefel */ 400 ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr); 401 ierr = VecDotBegin(dX, dXold, &dXdotdXold);CHKERRQ(ierr); 402 ierr = VecDotBegin(lX, dX, &lXdotdX);CHKERRQ(ierr); 403 ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 404 ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr); 405 ierr = VecDotEnd(dX, dXold, &dXdotdXold);CHKERRQ(ierr); 406 ierr = VecDotEnd(lX, dX, &lXdotdX);CHKERRQ(ierr); 407 ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 408 beta = PetscRealPart((dXdotdX - dXdotdXold) / (lXdotdX - lXdotdXold)); 409 break; 410 case SNES_NCG_DY: /* Dai-Yuan */ 411 ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr); 412 ierr = VecDotBegin(lX, dX, &lXdotdX);CHKERRQ(ierr); 413 ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 414 ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr); 415 ierr = VecDotEnd(lX, dX, &lXdotdX);CHKERRQ(ierr); 416 ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 417 beta = PetscRealPart(dXdotdX / (lXdotdXold - lXdotdX));CHKERRQ(ierr); 418 break; 419 case SNES_NCG_CD: /* Conjugate Descent */ 420 ierr = VecDotBegin(dX, dX, &dXdotdX);CHKERRQ(ierr); 421 ierr = VecDotBegin(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 422 ierr = VecDotEnd(dX, dX, &dXdotdX);CHKERRQ(ierr); 423 ierr = VecDotEnd(lX, dXold, &lXdotdXold);CHKERRQ(ierr); 424 beta = PetscRealPart(dXdotdX / lXdotdXold);CHKERRQ(ierr); 425 break; 426 } 427 if (ncg->monitor) { 428 ierr = PetscViewerASCIIPrintf(ncg->monitor, "beta = %e\n", (double)beta);CHKERRQ(ierr); 429 } 430 ierr = VecAYPX(lX, beta, dX);CHKERRQ(ierr); 431 } 432 ierr = PetscInfo1(snes, "Maximum number of iterations has been reached: %D\n", maxits);CHKERRQ(ierr); 433 if (!snes->reason) snes->reason = SNES_DIVERGED_MAX_IT; 434 PetscFunctionReturn(0); 435 } 436 437 438 439 /*MC 440 SNESNCG - Nonlinear Conjugate-Gradient method for the solution of general nonlinear systems. 441 442 Level: beginner 443 444 Options Database: 445 + -snes_ncg_type <fr, prp, dy, hs, cd> - Choice of conjugate-gradient update parameter, default is prp. 446 . -snes_linesearch_type <cp,l2,basic> - Line search type. 447 - -snes_ncg_monitor - Print relevant information about the ncg iteration. 448 449 Notes: 450 This solves the nonlinear system of equations F(x) = 0 using the nonlinear generalization of the conjugate 451 gradient method. This may be used with a nonlinear preconditioner used to pick the new search directions, but otherwise 452 chooses the initial search direction as F(x) for the initial guess x. 453 454 Only supports left non-linear preconditioning. 455 456 References: 457 . 1. - Peter R. Brune, Matthew G. Knepley, Barry F. Smith, and Xuemin Tu,"Composing Scalable Nonlinear Algebraic Solvers", 458 SIAM Review, 57(4), 2015 459 460 461 .seealso: SNESCreate(), SNES, SNESSetType(), SNESNEWTONLS, SNESNEWTONTR, SNESNGMRES, SNESQN 462 M*/ 463 PETSC_EXTERN PetscErrorCode SNESCreate_NCG(SNES snes) 464 { 465 PetscErrorCode ierr; 466 SNES_NCG * neP; 467 468 PetscFunctionBegin; 469 snes->ops->destroy = SNESDestroy_NCG; 470 snes->ops->setup = SNESSetUp_NCG; 471 snes->ops->setfromoptions = SNESSetFromOptions_NCG; 472 snes->ops->view = SNESView_NCG; 473 snes->ops->solve = SNESSolve_NCG; 474 snes->ops->reset = SNESReset_NCG; 475 476 snes->usesksp = PETSC_FALSE; 477 snes->usesnpc = PETSC_TRUE; 478 snes->npcside = PC_LEFT; 479 480 snes->alwayscomputesfinalresidual = PETSC_TRUE; 481 482 if (!snes->tolerancesset) { 483 snes->max_funcs = 30000; 484 snes->max_its = 10000; 485 snes->stol = 1e-20; 486 } 487 488 ierr = PetscNewLog(snes,&neP);CHKERRQ(ierr); 489 snes->data = (void*) neP; 490 neP->monitor = NULL; 491 neP->type = SNES_NCG_PRP; 492 ierr = PetscObjectComposeFunction((PetscObject)snes,"SNESNCGSetType_C", SNESNCGSetType_NCG);CHKERRQ(ierr); 493 PetscFunctionReturn(0); 494 } 495