Actual source code: nash.c
1: #define PETSCKSP_DLL
3: #include private/kspimpl.h
4: #include ../src/ksp/ksp/impls/cg/nash/nash.h
6: #define NASH_PRECONDITIONED_DIRECTION 0
7: #define NASH_UNPRECONDITIONED_DIRECTION 1
8: #define NASH_DIRECTION_TYPES 2
10: static const char *DType_Table[64] = {
11: "preconditioned", "unpreconditioned"
12: };
16: /*@
17: KSPNASHSetRadius - Sets the radius of the trust region.
19: Collective on KSP
21: Input Parameters:
22: + ksp - the iterative context
23: - radius - the trust region radius (Infinity is the default)
25: Options Database Key:
26: . -ksp_nash_radius <r>
28: Level: advanced
30: .keywords: KSP, NASH, set, trust region radius
31: @*/
32: PetscErrorCode KSPNASHSetRadius(KSP ksp, PetscReal radius)
33: {
34: PetscErrorCode ierr, (*f)(KSP, PetscReal);
38: if (radius < 0.0) SETERRQ(PETSC_ERR_ARG_OUTOFRANGE, "Radius negative");
39: PetscObjectQueryFunction((PetscObject)ksp, "KSPNASHSetRadius_C", (void (**)(void))&f);
40: if (f) {
41: (*f)(ksp, radius);
42: }
43: return(0);
44: }
48: /*@
49: KSPNASHGetNormD - Got norm of the direction.
51: Collective on KSP
53: Input Parameters:
54: + ksp - the iterative context
55: - norm_d - the norm of the direction
57: Level: advanced
59: .keywords: KSP, NASH, get, norm direction
60: @*/
61: PetscErrorCode KSPNASHGetNormD(KSP ksp, PetscReal *norm_d)
62: {
63: PetscErrorCode ierr, (*f)(KSP, PetscReal *);
67: PetscObjectQueryFunction((PetscObject)ksp, "KSPNASHGetNormD_C", (void (**)(void))&f);
68: if (f) {
69: (*f)(ksp, norm_d);
70: }
71: return(0);
72: }
76: /*@
77: KSPNASHGetObjFcn - Get objective function value.
79: Collective on KSP
81: Input Parameters:
82: + ksp - the iterative context
83: - o_fcn - the objective function value
85: Level: advanced
87: .keywords: KSP, NASH, get, objective function
88: @*/
89: PetscErrorCode KSPNASHGetObjFcn(KSP ksp, PetscReal *o_fcn)
90: {
91: PetscErrorCode ierr, (*f)(KSP, PetscReal *);
95: PetscObjectQueryFunction((PetscObject)ksp, "KSPNASHGetObjFcn_C", (void (**)(void))&f);
96: if (f) {
97: (*f)(ksp, o_fcn);
98: }
99: return(0);
100: }
104: /*
105: KSPSolve_NASH - Use preconditioned conjugate gradient to compute
106: an approximate minimizer of the quadratic function
108: q(s) = g^T * s + 0.5 * s^T * H * s
110: subject to the trust region constraint
112: || s || <= delta,
114: where
116: delta is the trust region radius,
117: g is the gradient vector,
118: H is the Hessian approximation, and
119: M is the positive definite preconditioner matrix.
121: KSPConvergedReason may be
122: $ KSP_CONVERGED_CG_NEG_CURVE if convergence is reached along a negative curvature direction,
123: $ KSP_CONVERGED_CG_CONSTRAINED if convergence is reached along a constrained step,
124: $ other KSP converged/diverged reasons
126: Notes:
127: The preconditioner supplied should be symmetric and positive definite.
128: */
129: PetscErrorCode KSPSolve_NASH(KSP ksp)
130: {
131: #ifdef PETSC_USE_COMPLEX
132: SETERRQ(PETSC_ERR_SUP, "NASH is not available for complex systems");
133: #else
134: KSP_NASH *cg = (KSP_NASH *)ksp->data;
137: MatStructure pflag;
138: Mat Qmat, Mmat;
139: Vec r, z, p, d;
140: PC pc;
142: PetscReal norm_r, norm_d, norm_dp1, norm_p, dMp;
143: PetscReal alpha, beta, kappa, rz, rzm1;
144: PetscReal rr, r2, step;
146: PetscInt max_cg_its;
148: PetscTruth diagonalscale;
152: /***************************************************************************/
153: /* Check the arguments and parameters. */
154: /***************************************************************************/
156: PCDiagonalScale(ksp->pc, &diagonalscale);
157: if (diagonalscale) {
158: SETERRQ1(PETSC_ERR_SUP, "Krylov method %s does not support diagonal scaling", ((PetscObject)ksp)->type_name);
159: }
161: if (cg->radius < 0.0) {
162: SETERRQ(PETSC_ERR_ARG_OUTOFRANGE, "Input error: radius < 0");
163: }
165: /***************************************************************************/
166: /* Get the workspace vectors and initialize variables */
167: /***************************************************************************/
169: r2 = cg->radius * cg->radius;
170: r = ksp->work[0];
171: z = ksp->work[1];
172: p = ksp->work[2];
173: d = ksp->vec_sol;
174: pc = ksp->pc;
176: PCGetOperators(pc, &Qmat, &Mmat, &pflag);
178: VecGetSize(d, &max_cg_its);
179: max_cg_its = PetscMin(max_cg_its, ksp->max_it);
180: ksp->its = 0;
182: /***************************************************************************/
183: /* Initialize objective function and direction. */
184: /***************************************************************************/
186: cg->o_fcn = 0.0;
188: VecSet(d, 0.0); /* d = 0 */
189: cg->norm_d = 0.0;
191: /***************************************************************************/
192: /* Begin the conjugate gradient method. Check the right-hand side for */
193: /* numerical problems. The check for not-a-number and infinite values */
194: /* need be performed only once. */
195: /***************************************************************************/
197: VecCopy(ksp->vec_rhs, r); /* r = -grad */
198: VecDot(r, r, &rr); /* rr = r^T r */
199: if (PetscIsInfOrNanScalar(rr)) {
200: /*************************************************************************/
201: /* The right-hand side contains not-a-number or an infinite value. */
202: /* The gradient step does not work; return a zero value for the step. */
203: /*************************************************************************/
205: ksp->reason = KSP_DIVERGED_NAN;
206: PetscInfo1(ksp, "KSPSolve_NASH: bad right-hand side: rr=%g\n", rr);
207: return(0);
208: }
210: /***************************************************************************/
211: /* Check the preconditioner for numerical problems and for positive */
212: /* definiteness. The check for not-a-number and infinite values need be */
213: /* performed only once. */
214: /***************************************************************************/
216: KSP_PCApply(ksp, r, z); /* z = inv(M) r */
217: VecDot(r, z, &rz); /* rz = r^T inv(M) r */
218: if (PetscIsInfOrNanScalar(rz)) {
219: /*************************************************************************/
220: /* The preconditioner contains not-a-number or an infinite value. */
221: /* Return the gradient direction intersected with the trust region. */
222: /*************************************************************************/
224: ksp->reason = KSP_DIVERGED_NAN;
225: PetscInfo1(ksp, "KSPSolve_NASH: bad preconditioner: rz=%g\n", rz);
227: if (cg->radius) {
228: if (r2 >= rr) {
229: alpha = 1.0;
230: cg->norm_d = sqrt(rr);
231: }
232: else {
233: alpha = sqrt(r2 / rr);
234: cg->norm_d = cg->radius;
235: }
237: VecAXPY(d, alpha, r); /* d = d + alpha r */
239: /***********************************************************************/
240: /* Compute objective function. */
241: /***********************************************************************/
243: KSP_MatMult(ksp, Qmat, d, z); CHKERRQ(ierr)
244: VecAYPX(z, -0.5, ksp->vec_rhs);
245: VecDot(d, z, &cg->o_fcn);
246: cg->o_fcn = -cg->o_fcn;
247: ++ksp->its;
248: }
249: return(0);
250: }
252: if (rz < 0.0) {
253: /*************************************************************************/
254: /* The preconditioner is indefinite. Because this is the first */
255: /* and we do not have a direction yet, we use the gradient step. Note */
256: /* that we cannot use the preconditioned norm when computing the step */
257: /* because the matrix is indefinite. */
258: /*************************************************************************/
260: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
261: PetscInfo1(ksp, "KSPSolve_NASH: indefinite preconditioner: rz=%g\n", rz);
263: if (cg->radius) {
264: if (r2 >= rr) {
265: alpha = 1.0;
266: cg->norm_d = sqrt(rr);
267: }
268: else {
269: alpha = sqrt(r2 / rr);
270: cg->norm_d = cg->radius;
271: }
273: VecAXPY(d, alpha, r); /* d = d + alpha r */
275: /***********************************************************************/
276: /* Compute objective function. */
277: /***********************************************************************/
279: KSP_MatMult(ksp, Qmat, d, z); CHKERRQ(ierr)
280: VecAYPX(z, -0.5, ksp->vec_rhs);
281: VecDot(d, z, &cg->o_fcn);
282: cg->o_fcn = -cg->o_fcn;
283: ++ksp->its;
284: }
285: return(0);
286: }
288: /***************************************************************************/
289: /* As far as we know, the preconditioner is positive semidefinite. */
290: /* Compute and log the residual. Check convergence because this */
291: /* initializes things, but do not terminate until at least one conjugate */
292: /* gradient iteration has been performed. */
293: /***************************************************************************/
295: switch(ksp->normtype) {
296: case KSP_NORM_PRECONDITIONED:
297: VecNorm(z, NORM_2, &norm_r); /* norm_r = |z| */
298: break;
300: case KSP_NORM_UNPRECONDITIONED:
301: norm_r = sqrt(rr); /* norm_r = |r| */
302: break;
304: case KSP_NORM_NATURAL:
305: norm_r = sqrt(rz); /* norm_r = |r|_M */
306: break;
308: default:
309: norm_r = 0.0;
310: break;
311: }
313: KSPLogResidualHistory(ksp, norm_r);
314: KSPMonitor(ksp, ksp->its, norm_r);
315: ksp->rnorm = norm_r;
317: (*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP);
319: /***************************************************************************/
320: /* Compute the first direction and update the iteration. */
321: /***************************************************************************/
323: VecCopy(z, p); /* p = z */
324: KSP_MatMult(ksp, Qmat, p, z); /* z = Q * p */
325: ++ksp->its;
327: /***************************************************************************/
328: /* Check the matrix for numerical problems. */
329: /***************************************************************************/
331: VecDot(p, z, &kappa); /* kappa = p^T Q p */
332: if (PetscIsInfOrNanScalar(kappa)) {
333: /*************************************************************************/
334: /* The matrix produced not-a-number or an infinite value. In this case, */
335: /* we must stop and use the gradient direction. This condition need */
336: /* only be checked once. */
337: /*************************************************************************/
339: ksp->reason = KSP_DIVERGED_NAN;
340: PetscInfo1(ksp, "KSPSolve_NASH: bad matrix: kappa=%g\n", kappa);
342: if (cg->radius) {
343: if (r2 >= rr) {
344: alpha = 1.0;
345: cg->norm_d = sqrt(rr);
346: }
347: else {
348: alpha = sqrt(r2 / rr);
349: cg->norm_d = cg->radius;
350: }
352: VecAXPY(d, alpha, r); /* d = d + alpha r */
354: /***********************************************************************/
355: /* Compute objective function. */
356: /***********************************************************************/
358: KSP_MatMult(ksp, Qmat, d, z); CHKERRQ(ierr)
359: VecAYPX(z, -0.5, ksp->vec_rhs);
360: VecDot(d, z, &cg->o_fcn);
361: cg->o_fcn = -cg->o_fcn;
362: ++ksp->its;
363: }
364: return(0);
365: }
367: /***************************************************************************/
368: /* Initialize variables for calculating the norm of the direction. */
369: /***************************************************************************/
371: dMp = 0.0;
372: norm_d = 0.0;
373: switch(cg->dtype) {
374: case NASH_PRECONDITIONED_DIRECTION:
375: norm_p = rz;
376: break;
378: default:
379: VecDot(p, p, &norm_p);
380: break;
381: }
383: /***************************************************************************/
384: /* Check for negative curvature. */
385: /***************************************************************************/
387: if (kappa <= 0.0) {
388: /*************************************************************************/
389: /* In this case, the matrix is indefinite and we have encountered a */
390: /* direction of negative curvature. Because negative curvature occurs */
391: /* during the first step, we must follow a direction. */
392: /*************************************************************************/
394: ksp->reason = KSP_CONVERGED_CG_NEG_CURVE;
395: PetscInfo1(ksp, "KSPSolve_NASH: negative curvature: kappa=%g\n", kappa);
397: if (cg->radius && norm_p > 0.0) {
398: /***********************************************************************/
399: /* Follow direction of negative curvature to the boundary of the */
400: /* trust region. */
401: /***********************************************************************/
403: step = sqrt(r2 / norm_p);
404: cg->norm_d = cg->radius;
406: VecAXPY(d, step, p); /* d = d + step p */
408: /***********************************************************************/
409: /* Update objective function. */
410: /***********************************************************************/
412: cg->o_fcn += step * (0.5 * step * kappa - rz);
413: }
414: else if (cg->radius) {
415: /***********************************************************************/
416: /* The norm of the preconditioned direction is zero; use the gradient */
417: /* step. */
418: /***********************************************************************/
420: if (r2 >= rr) {
421: alpha = 1.0;
422: cg->norm_d = sqrt(rr);
423: }
424: else {
425: alpha = sqrt(r2 / rr);
426: cg->norm_d = cg->radius;
427: }
429: VecAXPY(d, alpha, r); /* d = d + alpha r */
431: /***********************************************************************/
432: /* Compute objective function. */
433: /***********************************************************************/
435: KSP_MatMult(ksp, Qmat, d, z); CHKERRQ(ierr)
436: VecAYPX(z, -0.5, ksp->vec_rhs);
437: VecDot(d, z, &cg->o_fcn);
438: cg->o_fcn = -cg->o_fcn;
439: ++ksp->its;
440: }
441: return(0);
442: }
444: /***************************************************************************/
445: /* Run the conjugate gradient method until either the problem is solved, */
446: /* we encounter the boundary of the trust region, or the conjugate */
447: /* gradient method breaks down. */
448: /***************************************************************************/
450: while(1) {
451: /*************************************************************************/
452: /* Know that kappa is nonzero, because we have not broken down, so we */
453: /* can compute the steplength. */
454: /*************************************************************************/
456: alpha = rz / kappa;
458: /*************************************************************************/
459: /* Compute the steplength and check for intersection with the trust */
460: /* region. */
461: /*************************************************************************/
463: norm_dp1 = norm_d + alpha*(2.0*dMp + alpha*norm_p);
464: if (cg->radius && norm_dp1 >= r2) {
465: /***********************************************************************/
466: /* In this case, the matrix is positive definite as far as we know. */
467: /* However, the full step goes beyond the trust region. */
468: /***********************************************************************/
470: ksp->reason = KSP_CONVERGED_CG_CONSTRAINED;
471: PetscInfo1(ksp, "KSPSolve_NASH: constrained step: radius=%g\n", cg->radius);
473: if (norm_p > 0.0) {
474: /*********************************************************************/
475: /* Follow the direction to the boundary of the trust region. */
476: /*********************************************************************/
478: step = (sqrt(dMp*dMp+norm_p*(r2-norm_d))-dMp)/norm_p;
479: cg->norm_d = cg->radius;
481: VecAXPY(d, step, p); /* d = d + step p */
483: /*********************************************************************/
484: /* Update objective function. */
485: /*********************************************************************/
487: cg->o_fcn += step * (0.5 * step * kappa - rz);
488: }
489: else {
490: /*********************************************************************/
491: /* The norm of the direction is zero; there is nothing to follow. */
492: /*********************************************************************/
493: }
494: break;
495: }
497: /*************************************************************************/
498: /* Now we can update the direction and residual. */
499: /*************************************************************************/
501: VecAXPY(d, alpha, p); /* d = d + alpha p */
502: VecAXPY(r, -alpha, z); /* r = r - alpha Q p */
503: KSP_PCApply(ksp, r, z); /* z = inv(M) r */
505: switch(cg->dtype) {
506: case NASH_PRECONDITIONED_DIRECTION:
507: norm_d = norm_dp1;
508: break;
510: default:
511: VecDot(d, d, &norm_d);
512: break;
513: }
514: cg->norm_d = sqrt(norm_d);
516: /*************************************************************************/
517: /* Update objective function. */
518: /*************************************************************************/
520: cg->o_fcn -= 0.5 * alpha * rz;
522: /*************************************************************************/
523: /* Check that the preconditioner appears positive semidefinite. */
524: /*************************************************************************/
526: rzm1 = rz;
527: VecDot(r, z, &rz); /* rz = r^T z */
528: if (rz < 0.0) {
529: /***********************************************************************/
530: /* The preconditioner is indefinite. */
531: /***********************************************************************/
533: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
534: PetscInfo1(ksp, "KSPSolve_NASH: cg indefinite preconditioner: rz=%g\n", rz);
535: break;
536: }
538: /*************************************************************************/
539: /* As far as we know, the preconditioner is positive semidefinite. */
540: /* Compute the residual and check for convergence. */
541: /*************************************************************************/
543: switch(ksp->normtype) {
544: case KSP_NORM_PRECONDITIONED:
545: VecNorm(z, NORM_2, &norm_r); /* norm_r = |z| */
546: break;
548: case KSP_NORM_UNPRECONDITIONED:
549: VecNorm(r, NORM_2, &norm_r); /* norm_r = |r| */
550: break;
552: case KSP_NORM_NATURAL:
553: norm_r = sqrt(rz); /* norm_r = |r|_M */
554: break;
556: default:
557: norm_r = 0;
558: break;
559: }
561: KSPLogResidualHistory(ksp, norm_r);
562: KSPMonitor(ksp, ksp->its, norm_r);
563: ksp->rnorm = norm_r;
564:
565: (*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP);
566: if (ksp->reason) {
567: /***********************************************************************/
568: /* The method has converged. */
569: /***********************************************************************/
571: PetscInfo2(ksp, "KSPSolve_NASH: truncated step: rnorm=%g, radius=%g\n", norm_r, cg->radius);
572: break;
573: }
575: /*************************************************************************/
576: /* We have not converged yet. Check for breakdown. */
577: /*************************************************************************/
579: beta = rz / rzm1;
580: if (fabs(beta) <= 0.0) {
581: /***********************************************************************/
582: /* Conjugate gradients has broken down. */
583: /***********************************************************************/
585: ksp->reason = KSP_DIVERGED_BREAKDOWN;
586: PetscInfo1(ksp, "KSPSolve_NASH: breakdown: beta=%g\n", beta);
587: break;
588: }
590: /*************************************************************************/
591: /* Check iteration limit. */
592: /*************************************************************************/
594: if (ksp->its >= max_cg_its) {
595: ksp->reason = KSP_DIVERGED_ITS;
596: PetscInfo1(ksp, "KSPSolve_NASH: iterlim: its=%d\n", ksp->its);
597: break;
598: }
600: /*************************************************************************/
601: /* Update p and the norms. */
602: /*************************************************************************/
604: VecAYPX(p, beta, z); /* p = z + beta p */
606: switch(cg->dtype) {
607: case NASH_PRECONDITIONED_DIRECTION:
608: dMp = beta*(dMp + alpha*norm_p);
609: norm_p = beta*(rzm1 + beta*norm_p);
610: break;
612: default:
613: VecDot(d, p, &dMp);
614: VecDot(p, p, &norm_p);
615: break;
616: }
618: /*************************************************************************/
619: /* Compute the new direction and update the iteration. */
620: /*************************************************************************/
622: KSP_MatMult(ksp, Qmat, p, z); /* z = Q * p */
623: VecDot(p, z, &kappa); /* kappa = p^T Q p */
624: ++ksp->its;
626: /*************************************************************************/
627: /* Check for negative curvature. */
628: /*************************************************************************/
630: if (kappa <= 0.0) {
631: /***********************************************************************/
632: /* In this case, the matrix is indefinite and we have encountered */
633: /* a direction of negative curvature. Stop at the base. */
634: /***********************************************************************/
636: ksp->reason = KSP_CONVERGED_CG_NEG_CURVE;
637: PetscInfo1(ksp, "KSPSolve_NASH: negative curvature: kappa=%g\n", kappa);
638: break;
639: }
640: }
642: return(0);
643: #endif
644: }
648: PetscErrorCode KSPSetUp_NASH(KSP ksp)
649: {
654: /***************************************************************************/
655: /* This implementation of CG only handles left preconditioning so generate */
656: /* an error otherwise. */
657: /***************************************************************************/
659: if (ksp->pc_side == PC_RIGHT) {
660: SETERRQ(PETSC_ERR_SUP, "No right preconditioning for KSPNASH");
661: } else if (ksp->pc_side == PC_SYMMETRIC) {
662: SETERRQ(PETSC_ERR_SUP, "No symmetric preconditioning for KSPNASH");
663: }
665: /***************************************************************************/
666: /* Set work vectors needed by conjugate gradient method and allocate */
667: /***************************************************************************/
669: KSPDefaultGetWork(ksp, 3);
670: return(0);
671: }
675: PetscErrorCode KSPDestroy_NASH(KSP ksp)
676: {
681: /***************************************************************************/
682: /* Clear composed functions */
683: /***************************************************************************/
685: PetscObjectComposeFunctionDynamic((PetscObject)ksp,"KSPNASHSetRadius_C","",PETSC_NULL);
686: PetscObjectComposeFunctionDynamic((PetscObject)ksp,"KSPNASHGetNormD_C","",PETSC_NULL);
687: PetscObjectComposeFunctionDynamic((PetscObject)ksp,"KSPNASHGetObjFcn_C","",PETSC_NULL);
689: /***************************************************************************/
690: /* Destroy KSP object. */
691: /***************************************************************************/
693: KSPDefaultDestroy(ksp);
694: return(0);
695: }
700: PetscErrorCode KSPNASHSetRadius_NASH(KSP ksp, PetscReal radius)
701: {
702: KSP_NASH *cg = (KSP_NASH *)ksp->data;
705: cg->radius = radius;
706: return(0);
707: }
711: PetscErrorCode KSPNASHGetNormD_NASH(KSP ksp, PetscReal *norm_d)
712: {
713: KSP_NASH *cg = (KSP_NASH *)ksp->data;
716: *norm_d = cg->norm_d;
717: return(0);
718: }
722: PetscErrorCode KSPNASHGetObjFcn_NASH(KSP ksp, PetscReal *o_fcn){
723: KSP_NASH *cg = (KSP_NASH *)ksp->data;
726: *o_fcn = cg->o_fcn;
727: return(0);
728: }
733: PetscErrorCode KSPSetFromOptions_NASH(KSP ksp)
734: {
736: KSP_NASH *cg = (KSP_NASH *)ksp->data;
739: PetscOptionsHead("KSP NASH options");
741: PetscOptionsReal("-ksp_nash_radius", "Trust Region Radius", "KSPNASHSetRadius", cg->radius, &cg->radius, PETSC_NULL);
743: PetscOptionsEList("-ksp_nash_dtype", "Norm used for direction", "", DType_Table, NASH_DIRECTION_TYPES, DType_Table[cg->dtype], &cg->dtype, PETSC_NULL);
745: PetscOptionsTail();
746: return(0);
747: }
749: /*MC
750: KSPNASH - Code to run conjugate gradient method subject to a constraint
751: on the solution norm. This is used in Trust Region methods for
752: nonlinear equations, SNESTR
754: Options Database Keys:
755: . -ksp_nash_radius <r> - Trust Region Radius
757: Notes: This is rarely used directly
759: Level: developer
761: .seealso: KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP, KSPNASHSetRadius(), KSPNASHGetNormD(), KSPNASHGetObjFcn()
762: M*/
767: PetscErrorCode KSPCreate_NASH(KSP ksp)
768: {
770: KSP_NASH *cg;
774: PetscNewLog(ksp,KSP_NASH, &cg);
776: cg->radius = 0.0;
777: cg->dtype = NASH_UNPRECONDITIONED_DIRECTION;
779: ksp->data = (void *) cg;
780: ksp->pc_side = PC_LEFT;
781: ksp->normtype = KSP_NORM_UNPRECONDITIONED;
783: /***************************************************************************/
784: /* Sets the functions that are associated with this data structure */
785: /* (in C++ this is the same as defining virtual functions). */
786: /***************************************************************************/
788: ksp->ops->setup = KSPSetUp_NASH;
789: ksp->ops->solve = KSPSolve_NASH;
790: ksp->ops->destroy = KSPDestroy_NASH;
791: ksp->ops->setfromoptions = KSPSetFromOptions_NASH;
792: ksp->ops->buildsolution = KSPDefaultBuildSolution;
793: ksp->ops->buildresidual = KSPDefaultBuildResidual;
794: ksp->ops->view = 0;
796: PetscObjectComposeFunctionDynamic((PetscObject)ksp,
797: "KSPNASHSetRadius_C",
798: "KSPNASHSetRadius_NASH",
799: KSPNASHSetRadius_NASH);
800: PetscObjectComposeFunctionDynamic((PetscObject)ksp,
801: "KSPNASHGetNormD_C",
802: "KSPNASHGetNormD_NASH",
803: KSPNASHGetNormD_NASH);
804: PetscObjectComposeFunctionDynamic((PetscObject)ksp,
805: "KSPNASHGetObjFcn_C",
806: "KSPNASHGetObjFcn_NASH",
807: KSPNASHGetObjFcn_NASH);
808: return(0);
809: }