Actual source code: nash.c
1: #include <../src/ksp/ksp/impls/cg/nash/nashimpl.h>
3: #define NASH_PRECONDITIONED_DIRECTION 0
4: #define NASH_UNPRECONDITIONED_DIRECTION 1
5: #define NASH_DIRECTION_TYPES 2
7: static const char *DType_Table[64] = {"preconditioned", "unpreconditioned"};
9: static PetscErrorCode KSPCGSolve_NASH(KSP ksp)
10: {
11: #if defined(PETSC_USE_COMPLEX)
12: SETERRQ(PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "NASH is not available for complex systems");
13: #else
14: KSPCG_NASH *cg = (KSPCG_NASH *)ksp->data;
15: Mat Qmat, Mmat;
16: Vec r, z, p, d;
17: PC pc;
19: PetscReal norm_r, norm_d, norm_dp1, norm_p, dMp;
20: PetscReal alpha, beta, kappa, rz, rzm1;
21: PetscReal rr, r2, step;
23: PetscInt max_cg_its;
25: PetscBool diagonalscale;
27: PetscFunctionBegin;
28: /***************************************************************************/
29: /* Check the arguments and parameters. */
30: /***************************************************************************/
32: PetscCall(PCGetDiagonalScale(ksp->pc, &diagonalscale));
33: PetscCheck(!diagonalscale, PetscObjectComm((PetscObject)ksp), PETSC_ERR_SUP, "Krylov method %s does not support diagonal scaling", ((PetscObject)ksp)->type_name);
34: PetscCheck(cg->radius >= 0.0, PetscObjectComm((PetscObject)ksp), PETSC_ERR_ARG_OUTOFRANGE, "Input error: radius < 0");
36: /***************************************************************************/
37: /* Get the workspace vectors and initialize variables */
38: /***************************************************************************/
40: r2 = cg->radius * cg->radius;
41: r = ksp->work[0];
42: z = ksp->work[1];
43: p = ksp->work[2];
44: d = ksp->vec_sol;
45: pc = ksp->pc;
47: PetscCall(PCGetOperators(pc, &Qmat, &Mmat));
49: PetscCall(VecGetSize(d, &max_cg_its));
50: max_cg_its = PetscMin(max_cg_its, ksp->max_it);
51: ksp->its = 0;
53: /***************************************************************************/
54: /* Initialize objective function and direction. */
55: /***************************************************************************/
57: cg->o_fcn = 0.0;
59: PetscCall(VecSet(d, 0.0)); /* d = 0 */
60: cg->norm_d = 0.0;
62: /***************************************************************************/
63: /* Begin the conjugate gradient method. Check the right-hand side for */
64: /* numerical problems. The check for not-a-number and infinite values */
65: /* need be performed only once. */
66: /***************************************************************************/
68: PetscCall(VecCopy(ksp->vec_rhs, r)); /* r = -grad */
69: PetscCall(VecDot(r, r, &rr)); /* rr = r^T r */
70: KSPCheckDot(ksp, rr);
72: /***************************************************************************/
73: /* Check the preconditioner for numerical problems and for positive */
74: /* definiteness. The check for not-a-number and infinite values need be */
75: /* performed only once. */
76: /***************************************************************************/
78: PetscCall(KSP_PCApply(ksp, r, z)); /* z = inv(M) r */
79: PetscCall(VecDot(r, z, &rz)); /* rz = r^T inv(M) r */
80: if (PetscIsInfOrNanScalar(rz)) {
81: /*************************************************************************/
82: /* The preconditioner contains not-a-number or an infinite value. */
83: /* Return the gradient direction intersected with the trust region. */
84: /*************************************************************************/
86: ksp->reason = KSP_DIVERGED_NANORINF;
87: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: bad preconditioner: rz=%g\n", (double)rz));
89: if (cg->radius) {
90: if (r2 >= rr) {
91: alpha = 1.0;
92: cg->norm_d = PetscSqrtReal(rr);
93: } else {
94: alpha = PetscSqrtReal(r2 / rr);
95: cg->norm_d = cg->radius;
96: }
98: PetscCall(VecAXPY(d, alpha, r)); /* d = d + alpha r */
100: /***********************************************************************/
101: /* Compute objective function. */
102: /***********************************************************************/
104: PetscCall(KSP_MatMult(ksp, Qmat, d, z));
105: PetscCall(VecAYPX(z, -0.5, ksp->vec_rhs));
106: PetscCall(VecDot(d, z, &cg->o_fcn));
107: cg->o_fcn = -cg->o_fcn;
108: ++ksp->its;
109: }
110: PetscFunctionReturn(PETSC_SUCCESS);
111: }
113: if (rz < 0.0) {
114: /*************************************************************************/
115: /* The preconditioner is indefinite. Because this is the first */
116: /* and we do not have a direction yet, we use the gradient step. Note */
117: /* that we cannot use the preconditioned norm when computing the step */
118: /* because the matrix is indefinite. */
119: /*************************************************************************/
121: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
122: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: indefinite preconditioner: rz=%g\n", (double)rz));
124: if (cg->radius) {
125: if (r2 >= rr) {
126: alpha = 1.0;
127: cg->norm_d = PetscSqrtReal(rr);
128: } else {
129: alpha = PetscSqrtReal(r2 / rr);
130: cg->norm_d = cg->radius;
131: }
133: PetscCall(VecAXPY(d, alpha, r)); /* d = d + alpha r */
135: /***********************************************************************/
136: /* Compute objective function. */
137: /***********************************************************************/
139: PetscCall(KSP_MatMult(ksp, Qmat, d, z));
140: PetscCall(VecAYPX(z, -0.5, ksp->vec_rhs));
141: PetscCall(VecDot(d, z, &cg->o_fcn));
142: cg->o_fcn = -cg->o_fcn;
143: ++ksp->its;
144: }
145: PetscFunctionReturn(PETSC_SUCCESS);
146: }
148: /***************************************************************************/
149: /* As far as we know, the preconditioner is positive semidefinite. */
150: /* Compute and log the residual. Check convergence because this */
151: /* initializes things, but do not terminate until at least one conjugate */
152: /* gradient iteration has been performed. */
153: /***************************************************************************/
155: switch (ksp->normtype) {
156: case KSP_NORM_PRECONDITIONED:
157: PetscCall(VecNorm(z, NORM_2, &norm_r)); /* norm_r = |z| */
158: break;
160: case KSP_NORM_UNPRECONDITIONED:
161: norm_r = PetscSqrtReal(rr); /* norm_r = |r| */
162: break;
164: case KSP_NORM_NATURAL:
165: norm_r = PetscSqrtReal(rz); /* norm_r = |r|_M */
166: break;
168: default:
169: norm_r = 0.0;
170: break;
171: }
173: PetscCall(KSPLogResidualHistory(ksp, norm_r));
174: PetscCall(KSPMonitor(ksp, ksp->its, norm_r));
175: ksp->rnorm = norm_r;
177: PetscCall((*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP));
179: /***************************************************************************/
180: /* Compute the first direction and update the iteration. */
181: /***************************************************************************/
183: PetscCall(VecCopy(z, p)); /* p = z */
184: PetscCall(KSP_MatMult(ksp, Qmat, p, z)); /* z = Q * p */
185: ++ksp->its;
187: /***************************************************************************/
188: /* Check the matrix for numerical problems. */
189: /***************************************************************************/
191: PetscCall(VecDot(p, z, &kappa)); /* kappa = p^T Q p */
192: if (PetscIsInfOrNanScalar(kappa)) {
193: /*************************************************************************/
194: /* The matrix produced not-a-number or an infinite value. In this case, */
195: /* we must stop and use the gradient direction. This condition need */
196: /* only be checked once. */
197: /*************************************************************************/
199: ksp->reason = KSP_DIVERGED_NANORINF;
200: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: bad matrix: kappa=%g\n", (double)kappa));
202: if (cg->radius) {
203: if (r2 >= rr) {
204: alpha = 1.0;
205: cg->norm_d = PetscSqrtReal(rr);
206: } else {
207: alpha = PetscSqrtReal(r2 / rr);
208: cg->norm_d = cg->radius;
209: }
211: PetscCall(VecAXPY(d, alpha, r)); /* d = d + alpha r */
213: /***********************************************************************/
214: /* Compute objective function. */
215: /***********************************************************************/
217: PetscCall(KSP_MatMult(ksp, Qmat, d, z));
218: PetscCall(VecAYPX(z, -0.5, ksp->vec_rhs));
219: PetscCall(VecDot(d, z, &cg->o_fcn));
220: cg->o_fcn = -cg->o_fcn;
221: ++ksp->its;
222: }
223: PetscFunctionReturn(PETSC_SUCCESS);
224: }
226: /***************************************************************************/
227: /* Initialize variables for calculating the norm of the direction. */
228: /***************************************************************************/
230: dMp = 0.0;
231: norm_d = 0.0;
232: switch (cg->dtype) {
233: case NASH_PRECONDITIONED_DIRECTION:
234: norm_p = rz;
235: break;
237: default:
238: PetscCall(VecDot(p, p, &norm_p));
239: break;
240: }
242: /***************************************************************************/
243: /* Check for negative curvature. */
244: /***************************************************************************/
246: if (kappa <= 0.0) {
247: /*************************************************************************/
248: /* In this case, the matrix is indefinite and we have encountered a */
249: /* direction of negative curvature. Because negative curvature occurs */
250: /* during the first step, we must follow a direction. */
251: /*************************************************************************/
253: ksp->reason = ksp->converged_neg_curve ? KSP_CONVERGED_NEG_CURVE : KSP_DIVERGED_INDEFINITE_MAT;
254: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: negative curvature: kappa=%g\n", (double)kappa));
256: if (cg->radius && norm_p > 0.0) {
257: /***********************************************************************/
258: /* Follow direction of negative curvature to the boundary of the */
259: /* trust region. */
260: /***********************************************************************/
262: step = PetscSqrtReal(r2 / norm_p);
263: cg->norm_d = cg->radius;
265: PetscCall(VecAXPY(d, step, p)); /* d = d + step p */
267: /***********************************************************************/
268: /* Update objective function. */
269: /***********************************************************************/
271: cg->o_fcn += step * (0.5 * step * kappa - rz);
272: } else if (cg->radius) {
273: /***********************************************************************/
274: /* The norm of the preconditioned direction is zero; use the gradient */
275: /* step. */
276: /***********************************************************************/
278: if (r2 >= rr) {
279: alpha = 1.0;
280: cg->norm_d = PetscSqrtReal(rr);
281: } else {
282: alpha = PetscSqrtReal(r2 / rr);
283: cg->norm_d = cg->radius;
284: }
286: PetscCall(VecAXPY(d, alpha, r)); /* d = d + alpha r */
288: /***********************************************************************/
289: /* Compute objective function. */
290: /***********************************************************************/
292: PetscCall(KSP_MatMult(ksp, Qmat, d, z));
293: PetscCall(VecAYPX(z, -0.5, ksp->vec_rhs));
294: PetscCall(VecDot(d, z, &cg->o_fcn));
295: cg->o_fcn = -cg->o_fcn;
296: ++ksp->its;
297: }
298: PetscFunctionReturn(PETSC_SUCCESS);
299: }
301: /***************************************************************************/
302: /* Run the conjugate gradient method until either the problem is solved, */
303: /* we encounter the boundary of the trust region, or the conjugate */
304: /* gradient method breaks down. */
305: /***************************************************************************/
307: while (1) {
308: /*************************************************************************/
309: /* Know that kappa is nonzero, because we have not broken down, so we */
310: /* can compute the steplength. */
311: /*************************************************************************/
313: alpha = rz / kappa;
315: /*************************************************************************/
316: /* Compute the steplength and check for intersection with the trust */
317: /* region. */
318: /*************************************************************************/
320: norm_dp1 = norm_d + alpha * (2.0 * dMp + alpha * norm_p);
321: if (cg->radius && norm_dp1 >= r2) {
322: /***********************************************************************/
323: /* In this case, the matrix is positive definite as far as we know. */
324: /* However, the full step goes beyond the trust region. */
325: /***********************************************************************/
327: ksp->reason = KSP_CONVERGED_STEP_LENGTH;
328: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: constrained step: radius=%g\n", (double)cg->radius));
330: if (norm_p > 0.0) {
331: /*********************************************************************/
332: /* Follow the direction to the boundary of the trust region. */
333: /*********************************************************************/
335: step = (PetscSqrtReal(dMp * dMp + norm_p * (r2 - norm_d)) - dMp) / norm_p;
336: cg->norm_d = cg->radius;
338: PetscCall(VecAXPY(d, step, p)); /* d = d + step p */
340: /*********************************************************************/
341: /* Update objective function. */
342: /*********************************************************************/
344: cg->o_fcn += step * (0.5 * step * kappa - rz);
345: } else {
346: /*********************************************************************/
347: /* The norm of the direction is zero; there is nothing to follow. */
348: /*********************************************************************/
349: }
350: break;
351: }
353: /*************************************************************************/
354: /* Now we can update the direction and residual. */
355: /*************************************************************************/
357: PetscCall(VecAXPY(d, alpha, p)); /* d = d + alpha p */
358: PetscCall(VecAXPY(r, -alpha, z)); /* r = r - alpha Q p */
359: PetscCall(KSP_PCApply(ksp, r, z)); /* z = inv(M) r */
361: switch (cg->dtype) {
362: case NASH_PRECONDITIONED_DIRECTION:
363: norm_d = norm_dp1;
364: break;
366: default:
367: PetscCall(VecDot(d, d, &norm_d));
368: break;
369: }
370: cg->norm_d = PetscSqrtReal(norm_d);
372: /*************************************************************************/
373: /* Update objective function. */
374: /*************************************************************************/
376: cg->o_fcn -= 0.5 * alpha * rz;
378: /*************************************************************************/
379: /* Check that the preconditioner appears positive semidefinite. */
380: /*************************************************************************/
382: rzm1 = rz;
383: PetscCall(VecDot(r, z, &rz)); /* rz = r^T z */
384: if (rz < 0.0) {
385: /***********************************************************************/
386: /* The preconditioner is indefinite. */
387: /***********************************************************************/
389: ksp->reason = KSP_DIVERGED_INDEFINITE_PC;
390: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: cg indefinite preconditioner: rz=%g\n", (double)rz));
391: break;
392: }
394: /*************************************************************************/
395: /* As far as we know, the preconditioner is positive semidefinite. */
396: /* Compute the residual and check for convergence. */
397: /*************************************************************************/
399: switch (ksp->normtype) {
400: case KSP_NORM_PRECONDITIONED:
401: PetscCall(VecNorm(z, NORM_2, &norm_r)); /* norm_r = |z| */
402: break;
404: case KSP_NORM_UNPRECONDITIONED:
405: PetscCall(VecNorm(r, NORM_2, &norm_r)); /* norm_r = |r| */
406: break;
408: case KSP_NORM_NATURAL:
409: norm_r = PetscSqrtReal(rz); /* norm_r = |r|_M */
410: break;
412: default:
413: norm_r = 0.;
414: break;
415: }
417: PetscCall(KSPLogResidualHistory(ksp, norm_r));
418: PetscCall(KSPMonitor(ksp, ksp->its, norm_r));
419: ksp->rnorm = norm_r;
421: PetscCall((*ksp->converged)(ksp, ksp->its, norm_r, &ksp->reason, ksp->cnvP));
422: if (ksp->reason) {
423: /***********************************************************************/
424: /* The method has converged. */
425: /***********************************************************************/
427: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: truncated step: rnorm=%g, radius=%g\n", (double)norm_r, (double)cg->radius));
428: break;
429: }
431: /*************************************************************************/
432: /* We have not converged yet. Check for breakdown. */
433: /*************************************************************************/
435: beta = rz / rzm1;
436: if (PetscAbsReal(beta) <= 0.0) {
437: /***********************************************************************/
438: /* Conjugate gradients has broken down. */
439: /***********************************************************************/
441: ksp->reason = KSP_DIVERGED_BREAKDOWN;
442: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: breakdown: beta=%g\n", (double)beta));
443: break;
444: }
446: /*************************************************************************/
447: /* Check iteration limit. */
448: /*************************************************************************/
450: if (ksp->its >= max_cg_its) {
451: ksp->reason = KSP_DIVERGED_ITS;
452: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: iterlim: its=%" PetscInt_FMT "\n", ksp->its));
453: break;
454: }
456: /*************************************************************************/
457: /* Update p and the norms. */
458: /*************************************************************************/
460: PetscCall(VecAYPX(p, beta, z)); /* p = z + beta p */
462: switch (cg->dtype) {
463: case NASH_PRECONDITIONED_DIRECTION:
464: dMp = beta * (dMp + alpha * norm_p);
465: norm_p = beta * (rzm1 + beta * norm_p);
466: break;
468: default:
469: PetscCall(VecDot(d, p, &dMp));
470: PetscCall(VecDot(p, p, &norm_p));
471: break;
472: }
474: /*************************************************************************/
475: /* Compute the new direction and update the iteration. */
476: /*************************************************************************/
478: PetscCall(KSP_MatMult(ksp, Qmat, p, z)); /* z = Q * p */
479: PetscCall(VecDot(p, z, &kappa)); /* kappa = p^T Q p */
480: ++ksp->its;
482: /*************************************************************************/
483: /* Check for negative curvature. */
484: /*************************************************************************/
486: if (kappa <= 0.0) {
487: /***********************************************************************/
488: /* In this case, the matrix is indefinite and we have encountered */
489: /* a direction of negative curvature. Stop at the base. */
490: /***********************************************************************/
492: ksp->reason = ksp->converged_neg_curve ? KSP_CONVERGED_NEG_CURVE : KSP_DIVERGED_INDEFINITE_MAT;
493: PetscCall(PetscInfo(ksp, "KSPCGSolve_NASH: negative curvature: kappa=%g\n", (double)kappa));
494: break;
495: }
496: }
497: PetscFunctionReturn(PETSC_SUCCESS);
498: #endif
499: }
501: static PetscErrorCode KSPCGSetUp_NASH(KSP ksp)
502: {
503: /***************************************************************************/
504: /* Set work vectors needed by conjugate gradient method and allocate */
505: /***************************************************************************/
507: PetscFunctionBegin;
508: PetscCall(KSPSetWorkVecs(ksp, 3));
509: PetscFunctionReturn(PETSC_SUCCESS);
510: }
512: static PetscErrorCode KSPCGDestroy_NASH(KSP ksp)
513: {
514: PetscFunctionBegin;
515: /***************************************************************************/
516: /* Clear composed functions */
517: /***************************************************************************/
519: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGSetRadius_C", NULL));
520: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGGetNormD_C", NULL));
521: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGGetObjFcn_C", NULL));
523: /***************************************************************************/
524: /* Destroy KSP object. */
525: /***************************************************************************/
527: PetscCall(KSPDestroyDefault(ksp));
528: PetscFunctionReturn(PETSC_SUCCESS);
529: }
531: static PetscErrorCode KSPCGSetRadius_NASH(KSP ksp, PetscReal radius)
532: {
533: KSPCG_NASH *cg = (KSPCG_NASH *)ksp->data;
535: PetscFunctionBegin;
536: cg->radius = radius;
537: PetscFunctionReturn(PETSC_SUCCESS);
538: }
540: static PetscErrorCode KSPCGGetNormD_NASH(KSP ksp, PetscReal *norm_d)
541: {
542: KSPCG_NASH *cg = (KSPCG_NASH *)ksp->data;
544: PetscFunctionBegin;
545: *norm_d = cg->norm_d;
546: PetscFunctionReturn(PETSC_SUCCESS);
547: }
549: static PetscErrorCode KSPCGGetObjFcn_NASH(KSP ksp, PetscReal *o_fcn)
550: {
551: KSPCG_NASH *cg = (KSPCG_NASH *)ksp->data;
553: PetscFunctionBegin;
554: *o_fcn = cg->o_fcn;
555: PetscFunctionReturn(PETSC_SUCCESS);
556: }
558: static PetscErrorCode KSPCGSetFromOptions_NASH(KSP ksp, PetscOptionItems *PetscOptionsObject)
559: {
560: KSPCG_NASH *cg = (KSPCG_NASH *)ksp->data;
562: PetscFunctionBegin;
563: PetscOptionsHeadBegin(PetscOptionsObject, "KSPCG NASH options");
565: PetscCall(PetscOptionsReal("-ksp_cg_radius", "Trust Region Radius", "KSPCGSetRadius", cg->radius, &cg->radius, NULL));
567: PetscCall(PetscOptionsEList("-ksp_cg_dtype", "Norm used for direction", "", DType_Table, NASH_DIRECTION_TYPES, DType_Table[cg->dtype], &cg->dtype, NULL));
569: PetscOptionsHeadEnd();
570: PetscFunctionReturn(PETSC_SUCCESS);
571: }
573: /*MC
574: KSPNASH - Code to run conjugate gradient method subject to a constraint on the solution norm in a trust region method {cite}`nash1984newton`
576: Options Database Keys:
577: . -ksp_cg_radius <r> - Trust Region Radius
579: Level: developer
581: Notes:
582: This is rarely used directly, it is used in Trust Region methods for nonlinear equations, `SNESNEWTONTR`
584: Uses preconditioned conjugate gradient to compute
585: an approximate minimizer of the quadratic function
587: $$
588: q(s) = g^T * s + 0.5 * s^T * H * s
589: $$
591: subject to the trust region constraint
593: $$
594: || s || \le delta,
595: $$
597: where
598: .vb
599: delta is the trust region radius,
600: g is the gradient vector,
601: H is the Hessian approximation, and
602: M is the positive definite preconditioner matrix.
603: .ve
605: `KSPConvergedReason` may include
606: + `KSP_CONVERGED_NEG_CURVE` - if convergence is reached along a negative curvature direction,
607: - `KSP_CONVERGED_STEP_LENGTH` - if convergence is reached along a constrained step,
609: The preconditioner supplied must be symmetric and positive definite.
611: .seealso: [](ch_ksp), `KSPQCG`, `KSPGLTR`, `KSPSTCG`, `KSPCreate()`, `KSPSetType()`, `KSPType`, `KSP`, `KSPCGSetRadius()`, `KSPCGGetNormD()`, `KSPCGGetObjFcn()`
612: M*/
614: PETSC_EXTERN PetscErrorCode KSPCreate_NASH(KSP ksp)
615: {
616: KSPCG_NASH *cg;
618: PetscFunctionBegin;
619: PetscCall(PetscNew(&cg));
620: cg->radius = 0.0;
621: cg->dtype = NASH_UNPRECONDITIONED_DIRECTION;
623: ksp->data = (void *)cg;
624: PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_UNPRECONDITIONED, PC_LEFT, 3));
625: PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_PRECONDITIONED, PC_LEFT, 2));
626: PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_NATURAL, PC_LEFT, 2));
627: PetscCall(KSPSetSupportedNorm(ksp, KSP_NORM_NONE, PC_LEFT, 1));
628: PetscCall(KSPSetConvergedNegativeCurvature(ksp, PETSC_TRUE));
630: /***************************************************************************/
631: /* Sets the functions that are associated with this data structure */
632: /* (in C++ this is the same as defining virtual functions). */
633: /***************************************************************************/
635: ksp->ops->setup = KSPCGSetUp_NASH;
636: ksp->ops->solve = KSPCGSolve_NASH;
637: ksp->ops->destroy = KSPCGDestroy_NASH;
638: ksp->ops->setfromoptions = KSPCGSetFromOptions_NASH;
639: ksp->ops->buildsolution = KSPBuildSolutionDefault;
640: ksp->ops->buildresidual = KSPBuildResidualDefault;
641: ksp->ops->view = NULL;
643: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGSetRadius_C", KSPCGSetRadius_NASH));
644: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGGetNormD_C", KSPCGGetNormD_NASH));
645: PetscCall(PetscObjectComposeFunction((PetscObject)ksp, "KSPCGGetObjFcn_C", KSPCGGetObjFcn_NASH));
646: PetscFunctionReturn(PETSC_SUCCESS);
647: }