Actual source code: nash.c
petsc-3.5.4 2015-05-23
2: #include <petsc-private/kspimpl.h> /*I "petscksp.h" I*/
3: #include <../src/ksp/ksp/impls/cg/nash/nashimpl.h>
5: #define NASH_PRECONDITIONED_DIRECTION 0
6: #define NASH_UNPRECONDITIONED_DIRECTION 1
7: #define NASH_DIRECTION_TYPES 2
9: static const char *DType_Table[64] = { "preconditioned", "unpreconditioned"};
13: /*@
14: KSPNASHSetRadius - Sets the radius of the trust region.
16: Logically Collective on KSP
18: Input Parameters:
19: + ksp - the iterative context
20: - radius - the trust region radius (Infinity is the default)
22: Options Database Key:
23: . -ksp_nash_radius <r>
25: Level: advanced
27: .keywords: KSP, NASH, set, trust region radius
28: @*/
29: PetscErrorCode KSPNASHSetRadius(KSP ksp, PetscReal radius)
30: {
35: if (radius < 0.0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE, "Radius negative");
37: PetscTryMethod(ksp,"KSPNASHSetRadius_C",(KSP,PetscReal),(ksp,radius));
38: return(0);
39: }
43: /*@
44: KSPNASHGetNormD - Got norm of the direction.
46: Collective on KSP
48: Input Parameters:
49: + ksp - the iterative context
50: - norm_d - the norm of the direction
52: Level: advanced
54: .keywords: KSP, NASH, get, norm direction
55: @*/
56: PetscErrorCode KSPNASHGetNormD(KSP ksp, PetscReal *norm_d)
57: {
62: PetscUseMethod(ksp,"KSPNASHGetNormD_C",(KSP,PetscReal*),(ksp,norm_d));
63: return(0);
64: }
68: /*@
69: KSPNASHGetObjFcn - Get objective function value.
71: Collective on KSP
73: Input Parameters:
74: + ksp - the iterative context
75: - o_fcn - the objective function value
77: Level: advanced
79: .keywords: KSP, NASH, get, objective function
80: @*/
81: PetscErrorCode KSPNASHGetObjFcn(KSP ksp, PetscReal *o_fcn)
82: {
87: PetscUseMethod(ksp,"KSPNASHGetObjFcn_C",(KSP,PetscReal*),(ksp,o_fcn));
88: return(0);
89: }
94: PetscErrorCode KSPSolve_NASH(KSP ksp)
95: {
96: #if defined(PETSC_USE_COMPLEX)
97: SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP, "NASH is not available for complex systems");
98: #else
99: KSP_NASH *cg = (KSP_NASH*)ksp->data;
101: Mat Qmat, Mmat;
102: Vec r, z, p, d;
103: PC pc;
105: PetscReal norm_r, norm_d, norm_dp1, norm_p, dMp;
106: PetscReal alpha, beta, kappa, rz, rzm1;
107: PetscReal rr, r2, step;
109: PetscInt max_cg_its;
111: PetscBool diagonalscale;
114: /***************************************************************************/
115: /* Check the arguments and parameters. */
116: /***************************************************************************/
118: PCGetDiagonalScale(ksp->pc, &diagonalscale);
119: if (diagonalscale) SETERRQ1(PetscObjectComm((PetscObject)ksp),PETSC_ERR_SUP, "Krylov method %s does not support diagonal scaling", ((PetscObject)ksp)->type_name);
120: if (cg->radius < 0.0) SETERRQ(PetscObjectComm((PetscObject)ksp),PETSC_ERR_ARG_OUTOFRANGE, "Input error: radius < 0");
122: /***************************************************************************/
123: /* Get the workspace vectors and initialize variables */
124: /***************************************************************************/
126: r2 = cg->radius * cg->radius;
127: r = ksp->work[0];
128: z = ksp->work[1];
129: p = ksp->work[2];
130: d = ksp->vec_sol;
131: pc = ksp->pc;
133: PCGetOperators(pc, &Qmat, &Mmat);
135: VecGetSize(d, &max_cg_its);
136: max_cg_its = PetscMin(max_cg_its, ksp->max_it);
137: ksp->its = 0;
139: /***************************************************************************/
140: /* Initialize objective function and direction. */
141: /***************************************************************************/
143: cg->o_fcn = 0.0;
145: VecSet(d, 0.0); /* d = 0 */
146: cg->norm_d = 0.0;
148: /***************************************************************************/
149: /* Begin the conjugate gradient method. Check the right-hand side for */
150: /* numerical problems. The check for not-a-number and infinite values */
151: /* need be performed only once. */
152: /***************************************************************************/
154: VecCopy(ksp->vec_rhs, r); /* r = -grad */
155: VecDot(r, r, &rr); /* rr = r^T r */
156: if (PetscIsInfOrNanScalar(rr)) {
157: /*************************************************************************/
158: /* The right-hand side contains not-a-number or an infinite value. */
159: /* The gradient step does not work; return a zero value for the step. */
160: /*************************************************************************/
162: ksp->reason = KSP_DIVERGED_NANORINF;
163: PetscInfo1(ksp, "KSPSolve_NASH: bad right-hand side: rr=%g\n", (double)rr);
164: return(0);
165: }
167: /***************************************************************************/
168: /* Check the preconditioner for numerical problems and for positive */
169: /* definiteness. The check for not-a-number and infinite values need be */
170: /* performed only once. */
171: /***************************************************************************/
173: KSP_PCApply(ksp, r, z); /* z = inv(M) r */
174: VecDot(r, z, &rz); /* rz = r^T inv(M) r */
175: if (PetscIsInfOrNanScalar(rz)) {
176: /*************************************************************************/
177: /* The preconditioner contains not-a-number or an infinite value. */
178: /* Return the gradient direction intersected with the trust region. */
179: /*************************************************************************/
181: ksp->reason = KSP_DIVERGED_NANORINF;
182: PetscInfo1(ksp, "KSPSolve_NASH: bad preconditioner: rz=%g\n", (double)rz);
184: if (cg->radius) {
185: if (r2 >= rr) {
186: alpha = 1.0;
187: cg->norm_d = PetscSqrtReal(rr);
188: } else {
189: alpha = PetscSqrtReal(r2 / rr);
190: cg->norm_d = cg->radius;
191: }
193: VecAXPY(d, alpha, r); /* d = d + alpha r */
195: /***********************************************************************/
196: /* Compute objective function. */
197: /***********************************************************************/
199: KSP_MatMult(ksp, Qmat, d, z);
200: VecAYPX(z, -0.5, ksp->vec_rhs);
201: VecDot(d, z, &cg->o_fcn);
202: cg->o_fcn = -cg->o_fcn;
203: ++ksp->its;
204: }
205: return(0);
206: }
208: if (rz < 0.0) {
209: /*************************************************************************/
210: /* The preconditioner is indefinite. Because this is the first */
211: /* and we do not have a direction yet, we use the gradient step. Note */
212: /* that we cannot use the preconditioned norm when computing the step */
213: /* because the matrix is indefinite. */
214: /*************************************************************************/
216: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
217: PetscInfo1(ksp, "KSPSolve_NASH: indefinite preconditioner: rz=%g\n", (double)rz);
219: if (cg->radius) {
220: if (r2 >= rr) {
221: alpha = 1.0;
222: cg->norm_d = PetscSqrtReal(rr);
223: } else {
224: alpha = PetscSqrtReal(r2 / rr);
225: cg->norm_d = cg->radius;
226: }
228: VecAXPY(d, alpha, r); /* d = d + alpha r */
230: /***********************************************************************/
231: /* Compute objective function. */
232: /***********************************************************************/
234: KSP_MatMult(ksp, Qmat, d, z);
235: VecAYPX(z, -0.5, ksp->vec_rhs);
236: VecDot(d, z, &cg->o_fcn);
237: cg->o_fcn = -cg->o_fcn;
238: ++ksp->its;
239: }
240: return(0);
241: }
243: /***************************************************************************/
244: /* As far as we know, the preconditioner is positive semidefinite. */
245: /* Compute and log the residual. Check convergence because this */
246: /* initializes things, but do not terminate until at least one conjugate */
247: /* gradient iteration has been performed. */
248: /***************************************************************************/
250: switch (ksp->normtype) {
251: case KSP_NORM_PRECONDITIONED:
252: VecNorm(z, NORM_2, &norm_r); /* norm_r = |z| */
253: break;
255: case KSP_NORM_UNPRECONDITIONED:
256: norm_r = PetscSqrtReal(rr); /* norm_r = |r| */
257: break;
259: case KSP_NORM_NATURAL:
260: norm_r = PetscSqrtReal(rz); /* norm_r = |r|_M */
261: break;
263: default:
264: norm_r = 0.0;
265: break;
266: }
268: KSPLogResidualHistory(ksp, norm_r);
269: KSPMonitor(ksp, ksp->its, norm_r);
270: ksp->rnorm = norm_r;
272: (*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP);
274: /***************************************************************************/
275: /* Compute the first direction and update the iteration. */
276: /***************************************************************************/
278: VecCopy(z, p); /* p = z */
279: KSP_MatMult(ksp, Qmat, p, z); /* z = Q * p */
280: ++ksp->its;
282: /***************************************************************************/
283: /* Check the matrix for numerical problems. */
284: /***************************************************************************/
286: VecDot(p, z, &kappa); /* kappa = p^T Q p */
287: if (PetscIsInfOrNanScalar(kappa)) {
288: /*************************************************************************/
289: /* The matrix produced not-a-number or an infinite value. In this case, */
290: /* we must stop and use the gradient direction. This condition need */
291: /* only be checked once. */
292: /*************************************************************************/
294: ksp->reason = KSP_DIVERGED_NANORINF;
295: PetscInfo1(ksp, "KSPSolve_NASH: bad matrix: kappa=%g\n", (double)kappa);
297: if (cg->radius) {
298: if (r2 >= rr) {
299: alpha = 1.0;
300: cg->norm_d = PetscSqrtReal(rr);
301: } else {
302: alpha = PetscSqrtReal(r2 / rr);
303: cg->norm_d = cg->radius;
304: }
306: VecAXPY(d, alpha, r); /* d = d + alpha r */
308: /***********************************************************************/
309: /* Compute objective function. */
310: /***********************************************************************/
312: KSP_MatMult(ksp, Qmat, d, z);
313: VecAYPX(z, -0.5, ksp->vec_rhs);
314: VecDot(d, z, &cg->o_fcn);
315: cg->o_fcn = -cg->o_fcn;
316: ++ksp->its;
317: }
318: return(0);
319: }
321: /***************************************************************************/
322: /* Initialize variables for calculating the norm of the direction. */
323: /***************************************************************************/
325: dMp = 0.0;
326: norm_d = 0.0;
327: switch (cg->dtype) {
328: case NASH_PRECONDITIONED_DIRECTION:
329: norm_p = rz;
330: break;
332: default:
333: VecDot(p, p, &norm_p);
334: break;
335: }
337: /***************************************************************************/
338: /* Check for negative curvature. */
339: /***************************************************************************/
341: if (kappa <= 0.0) {
342: /*************************************************************************/
343: /* In this case, the matrix is indefinite and we have encountered a */
344: /* direction of negative curvature. Because negative curvature occurs */
345: /* during the first step, we must follow a direction. */
346: /*************************************************************************/
348: ksp->reason = KSP_CONVERGED_CG_NEG_CURVE;
349: PetscInfo1(ksp, "KSPSolve_NASH: negative curvature: kappa=%g\n", (double)kappa);
351: if (cg->radius && norm_p > 0.0) {
352: /***********************************************************************/
353: /* Follow direction of negative curvature to the boundary of the */
354: /* trust region. */
355: /***********************************************************************/
357: step = PetscSqrtReal(r2 / norm_p);
358: cg->norm_d = cg->radius;
360: VecAXPY(d, step, p); /* d = d + step p */
362: /***********************************************************************/
363: /* Update objective function. */
364: /***********************************************************************/
366: cg->o_fcn += step * (0.5 * step * kappa - rz);
367: } else if (cg->radius) {
368: /***********************************************************************/
369: /* The norm of the preconditioned direction is zero; use the gradient */
370: /* step. */
371: /***********************************************************************/
373: if (r2 >= rr) {
374: alpha = 1.0;
375: cg->norm_d = PetscSqrtReal(rr);
376: } else {
377: alpha = PetscSqrtReal(r2 / rr);
378: cg->norm_d = cg->radius;
379: }
381: VecAXPY(d, alpha, r); /* d = d + alpha r */
383: /***********************************************************************/
384: /* Compute objective function. */
385: /***********************************************************************/
387: KSP_MatMult(ksp, Qmat, d, z);
388: VecAYPX(z, -0.5, ksp->vec_rhs);
389: VecDot(d, z, &cg->o_fcn);
390: cg->o_fcn = -cg->o_fcn;
391: ++ksp->its;
392: }
393: return(0);
394: }
396: /***************************************************************************/
397: /* Run the conjugate gradient method until either the problem is solved, */
398: /* we encounter the boundary of the trust region, or the conjugate */
399: /* gradient method breaks down. */
400: /***************************************************************************/
402: while (1) {
403: /*************************************************************************/
404: /* Know that kappa is nonzero, because we have not broken down, so we */
405: /* can compute the steplength. */
406: /*************************************************************************/
408: alpha = rz / kappa;
410: /*************************************************************************/
411: /* Compute the steplength and check for intersection with the trust */
412: /* region. */
413: /*************************************************************************/
415: norm_dp1 = norm_d + alpha*(2.0*dMp + alpha*norm_p);
416: if (cg->radius && norm_dp1 >= r2) {
417: /***********************************************************************/
418: /* In this case, the matrix is positive definite as far as we know. */
419: /* However, the full step goes beyond the trust region. */
420: /***********************************************************************/
422: ksp->reason = KSP_CONVERGED_CG_CONSTRAINED;
423: PetscInfo1(ksp, "KSPSolve_NASH: constrained step: radius=%g\n", (double)cg->radius);
425: if (norm_p > 0.0) {
426: /*********************************************************************/
427: /* Follow the direction to the boundary of the trust region. */
428: /*********************************************************************/
430: step = (PetscSqrtReal(dMp*dMp+norm_p*(r2-norm_d))-dMp)/norm_p;
431: cg->norm_d = cg->radius;
433: VecAXPY(d, step, p); /* d = d + step p */
435: /*********************************************************************/
436: /* Update objective function. */
437: /*********************************************************************/
439: cg->o_fcn += step * (0.5 * step * kappa - rz);
440: } else {
441: /*********************************************************************/
442: /* The norm of the direction is zero; there is nothing to follow. */
443: /*********************************************************************/
444: }
445: break;
446: }
448: /*************************************************************************/
449: /* Now we can update the direction and residual. */
450: /*************************************************************************/
452: VecAXPY(d, alpha, p); /* d = d + alpha p */
453: VecAXPY(r, -alpha, z); /* r = r - alpha Q p */
454: KSP_PCApply(ksp, r, z); /* z = inv(M) r */
456: switch (cg->dtype) {
457: case NASH_PRECONDITIONED_DIRECTION:
458: norm_d = norm_dp1;
459: break;
461: default:
462: VecDot(d, d, &norm_d);
463: break;
464: }
465: cg->norm_d = PetscSqrtReal(norm_d);
467: /*************************************************************************/
468: /* Update objective function. */
469: /*************************************************************************/
471: cg->o_fcn -= 0.5 * alpha * rz;
473: /*************************************************************************/
474: /* Check that the preconditioner appears positive semidefinite. */
475: /*************************************************************************/
477: rzm1 = rz;
478: VecDot(r, z, &rz); /* rz = r^T z */
479: if (rz < 0.0) {
480: /***********************************************************************/
481: /* The preconditioner is indefinite. */
482: /***********************************************************************/
484: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
485: PetscInfo1(ksp, "KSPSolve_NASH: cg indefinite preconditioner: rz=%g\n", (double)rz);
486: break;
487: }
489: /*************************************************************************/
490: /* As far as we know, the preconditioner is positive semidefinite. */
491: /* Compute the residual and check for convergence. */
492: /*************************************************************************/
494: switch (ksp->normtype) {
495: case KSP_NORM_PRECONDITIONED:
496: VecNorm(z, NORM_2, &norm_r); /* norm_r = |z| */
497: break;
499: case KSP_NORM_UNPRECONDITIONED:
500: VecNorm(r, NORM_2, &norm_r); /* norm_r = |r| */
501: break;
503: case KSP_NORM_NATURAL:
504: norm_r = PetscSqrtReal(rz); /* norm_r = |r|_M */
505: break;
507: default:
508: norm_r = 0.;
509: break;
510: }
512: KSPLogResidualHistory(ksp, norm_r);
513: KSPMonitor(ksp, ksp->its, norm_r);
514: ksp->rnorm = norm_r;
516: (*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP);
517: if (ksp->reason) {
518: /***********************************************************************/
519: /* The method has converged. */
520: /***********************************************************************/
522: PetscInfo2(ksp, "KSPSolve_NASH: truncated step: rnorm=%g, radius=%g\n", (double)norm_r, (double)cg->radius);
523: break;
524: }
526: /*************************************************************************/
527: /* We have not converged yet. Check for breakdown. */
528: /*************************************************************************/
530: beta = rz / rzm1;
531: if (PetscAbsReal(beta) <= 0.0) {
532: /***********************************************************************/
533: /* Conjugate gradients has broken down. */
534: /***********************************************************************/
536: ksp->reason = KSP_DIVERGED_BREAKDOWN;
537: PetscInfo1(ksp, "KSPSolve_NASH: breakdown: beta=%g\n", (double)beta);
538: break;
539: }
541: /*************************************************************************/
542: /* Check iteration limit. */
543: /*************************************************************************/
545: if (ksp->its >= max_cg_its) {
546: ksp->reason = KSP_DIVERGED_ITS;
547: PetscInfo1(ksp, "KSPSolve_NASH: iterlim: its=%D\n", ksp->its);
548: break;
549: }
551: /*************************************************************************/
552: /* Update p and the norms. */
553: /*************************************************************************/
555: VecAYPX(p, beta, z); /* p = z + beta p */
557: switch (cg->dtype) {
558: case NASH_PRECONDITIONED_DIRECTION:
559: dMp = beta*(dMp + alpha*norm_p);
560: norm_p = beta*(rzm1 + beta*norm_p);
561: break;
563: default:
564: VecDot(d, p, &dMp);
565: VecDot(p, p, &norm_p);
566: break;
567: }
569: /*************************************************************************/
570: /* Compute the new direction and update the iteration. */
571: /*************************************************************************/
573: KSP_MatMult(ksp, Qmat, p, z); /* z = Q * p */
574: VecDot(p, z, &kappa); /* kappa = p^T Q p */
575: ++ksp->its;
577: /*************************************************************************/
578: /* Check for negative curvature. */
579: /*************************************************************************/
581: if (kappa <= 0.0) {
582: /***********************************************************************/
583: /* In this case, the matrix is indefinite and we have encountered */
584: /* a direction of negative curvature. Stop at the base. */
585: /***********************************************************************/
587: ksp->reason = KSP_CONVERGED_CG_NEG_CURVE;
588: PetscInfo1(ksp, "KSPSolve_NASH: negative curvature: kappa=%g\n", (double)kappa);
589: break;
590: }
591: }
592: return(0);
593: #endif
594: }
598: PetscErrorCode KSPSetUp_NASH(KSP ksp)
599: {
603: /***************************************************************************/
604: /* Set work vectors needed by conjugate gradient method and allocate */
605: /***************************************************************************/
607: KSPSetWorkVecs(ksp, 3);
608: return(0);
609: }
613: PetscErrorCode KSPDestroy_NASH(KSP ksp)
614: {
618: /***************************************************************************/
619: /* Clear composed functions */
620: /***************************************************************************/
622: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHSetRadius_C",NULL);
623: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHGetNormD_C",NULL);
624: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHGetObjFcn_C",NULL);
626: /***************************************************************************/
627: /* Destroy KSP object. */
628: /***************************************************************************/
630: KSPDestroyDefault(ksp);
631: return(0);
632: }
636: static PetscErrorCode KSPNASHSetRadius_NASH(KSP ksp, PetscReal radius)
637: {
638: KSP_NASH *cg = (KSP_NASH*)ksp->data;
641: cg->radius = radius;
642: return(0);
643: }
647: static PetscErrorCode KSPNASHGetNormD_NASH(KSP ksp, PetscReal *norm_d)
648: {
649: KSP_NASH *cg = (KSP_NASH*)ksp->data;
652: *norm_d = cg->norm_d;
653: return(0);
654: }
658: static PetscErrorCode KSPNASHGetObjFcn_NASH(KSP ksp, PetscReal *o_fcn)
659: {
660: KSP_NASH *cg = (KSP_NASH*)ksp->data;
663: *o_fcn = cg->o_fcn;
664: return(0);
665: }
669: PetscErrorCode KSPSetFromOptions_NASH(KSP ksp)
670: {
672: KSP_NASH *cg = (KSP_NASH*)ksp->data;
675: PetscOptionsHead("KSP NASH options");
677: PetscOptionsReal("-ksp_nash_radius", "Trust Region Radius", "KSPNASHSetRadius", cg->radius, &cg->radius, NULL);
679: PetscOptionsEList("-ksp_nash_dtype", "Norm used for direction", "", DType_Table, NASH_DIRECTION_TYPES, DType_Table[cg->dtype], &cg->dtype, NULL);
681: PetscOptionsTail();
682: return(0);
683: }
685: /*MC
686: KSPNASH - Code to run conjugate gradient method subject to a constraint
687: on the solution norm. This is used in Trust Region methods for
688: nonlinear equations, SNESNEWTONTR
690: Options Database Keys:
691: . -ksp_nash_radius <r> - Trust Region Radius
693: Notes: This is rarely used directly
695: Level: developer
697: Use preconditioned conjugate gradient to compute
698: an approximate minimizer of the quadratic function
700: q(s) = g^T * s + 0.5 * s^T * H * s
702: subject to the trust region constraint
704: || s || <= delta,
706: where
708: delta is the trust region radius,
709: g is the gradient vector,
710: H is the Hessian approximation, and
711: M is the positive definite preconditioner matrix.
713: KSPConvergedReason may be
714: $ KSP_CONVERGED_CG_NEG_CURVE if convergence is reached along a negative curvature direction,
715: $ KSP_CONVERGED_CG_CONSTRAINED if convergence is reached along a constrained step,
716: $ other KSP converged/diverged reasons
718: Notes:
719: The preconditioner supplied should be symmetric and positive definite.
721: .seealso: KSPCreate(), KSPSetType(), KSPType (for list of available types), KSP, KSPNASHSetRadius(), KSPNASHGetNormD(), KSPNASHGetObjFcn()
722: M*/
726: PETSC_EXTERN PetscErrorCode KSPCreate_NASH(KSP ksp)
727: {
729: KSP_NASH *cg;
732: PetscNewLog(ksp,&cg);
733: cg->radius = 0.0;
734: cg->dtype = NASH_UNPRECONDITIONED_DIRECTION;
736: ksp->data = (void*) cg;
737: KSPSetSupportedNorm(ksp,KSP_NORM_UNPRECONDITIONED,PC_LEFT,3);
738: KSPSetSupportedNorm(ksp,KSP_NORM_PRECONDITIONED,PC_LEFT,2);
739: KSPSetSupportedNorm(ksp,KSP_NORM_NATURAL,PC_LEFT,2);
741: /***************************************************************************/
742: /* Sets the functions that are associated with this data structure */
743: /* (in C++ this is the same as defining virtual functions). */
744: /***************************************************************************/
746: ksp->ops->setup = KSPSetUp_NASH;
747: ksp->ops->solve = KSPSolve_NASH;
748: ksp->ops->destroy = KSPDestroy_NASH;
749: ksp->ops->setfromoptions = KSPSetFromOptions_NASH;
750: ksp->ops->buildsolution = KSPBuildSolutionDefault;
751: ksp->ops->buildresidual = KSPBuildResidualDefault;
752: ksp->ops->view = 0;
754: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHSetRadius_C",KSPNASHSetRadius_NASH);
755: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHGetNormD_C",KSPNASHGetNormD_NASH);
756: PetscObjectComposeFunction((PetscObject)ksp,"KSPNASHGetObjFcn_C",KSPNASHGetObjFcn_NASH);
757: return(0);
758: }