Actual source code: cgtype.c
  1: #include <../src/ksp/ksp/impls/cg/cgimpl.h>
  3: /*@
  4:   KSPCGSetType - Sets the variant of the conjugate gradient method to
  5:   use for solving a linear system with a complex coefficient matrix.
  6:   This option is irrelevant when solving a real system.
  8:   Logically Collective
 10:   Input Parameters:
 11: + ksp  - the iterative context
 12: - type - the variant of CG to use, one of
 13: .vb
 14:       KSP_CG_HERMITIAN - complex, Hermitian matrix (default)
 15:       KSP_CG_SYMMETRIC - complex, symmetric matrix
 16: .ve
 18:   Options Database Keys:
 19: + -ksp_cg_type hermitian - Indicates Hermitian matrix
 20: - -ksp_cg_type symmetric - Indicates symmetric matrix
 22:   Level: intermediate
 24:   Note:
 25:   By default, the matrix is assumed to be complex, Hermitian.
 27: .seealso: [](ch_ksp), `KSP`, `KSPCG`
 28: @*/
 29: PetscErrorCode KSPCGSetType(KSP ksp, KSPCGType type)
 30: {
 31:   PetscFunctionBegin;
 33:   PetscTryMethod(ksp, "KSPCGSetType_C", (KSP, KSPCGType), (ksp, type));
 34:   PetscFunctionReturn(PETSC_SUCCESS);
 35: }
 37: /*@
 38:   KSPCGUseSingleReduction - Merge the two inner products needed in `KSPCG` into a single `MPI_Allreduce()` call.
 40:   Logically Collective
 42:   Input Parameters:
 43: + ksp - the iterative context
 44: - flg - turn on or off the single reduction
 46:   Options Database Key:
 47: . -ksp_cg_single_reduction <bool> - Merge inner products into single `MPI_Allreduce()`
 49:   Level: intermediate
 51:   Notes:
 52:   The algorithm used in this case is described as Method 1 in {cite}`d1993conjugate`. V. Eijkhout credits the algorithm initially to Chronopoulos and Gear.
 54:   It requires two extra work vectors than the conventional implementation in PETSc.
 56:   See also `KSPPIPECG`, `KSPPIPECR`, and `KSPGROPPCG` that use non-blocking reductions. [](sec_pipelineksp),
 58: .seealso: [](ch_ksp), [](sec_pipelineksp), `KSP`, `KSPCG`, `KSPGMRES`, `KSPPIPECG`, `KSPPIPECR`, `and KSPGROPPCG`
 59: @*/
 60: PetscErrorCode KSPCGUseSingleReduction(KSP ksp, PetscBool flg)
 61: {
 62:   PetscFunctionBegin;
 65:   PetscTryMethod(ksp, "KSPCGUseSingleReduction_C", (KSP, PetscBool), (ksp, flg));
 66:   PetscFunctionReturn(PETSC_SUCCESS);
 67: }
 69: /*@
 70:   KSPCGSetRadius - Sets the radius of the trust region used by the `KSPCG` when the solver is used inside `SNESNEWTONTR`
 72:   Logically Collective
 74:   Input Parameters:
 75: + ksp    - the iterative context
 76: - radius - the trust region radius (0 is the default that disable the use of the radius)
 78:   Level: advanced
 80:   Note:
 81:   When radius is greater then 0, the Steihaugh-Toint trick is used
 83: .seealso: [](ch_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`, `SNESNEWTONTR`
 84: @*/
 85: PetscErrorCode KSPCGSetRadius(KSP ksp, PetscReal radius)
 86: {
 87:   PetscFunctionBegin;
 90:   PetscTryMethod(ksp, "KSPCGSetRadius_C", (KSP, PetscReal), (ksp, radius));
 91:   PetscFunctionReturn(PETSC_SUCCESS);
 92: }
 94: /*@
 95:   KSPCGSetObjectiveTarget - Sets the target value for the CG quadratic model
 97:   Logically Collective
 99:   Input Parameters:
100: + ksp - the iterative context
101: - obj - the objective value (0 is the default)
103:   Level: advanced
105:   Notes:
106:   The `KSPSolve()` will stop when the current value of the objective function $\frac{1}{2} x^H_k A x_k - b^H x_k$
107:   is smaller than `obj` if `obj` is negative. Otherwise the test is ignored.
108:   This is used for example inside `SNESNEWTONTR`.
110: .seealso: [](ch_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`, `SNESNEWTONTR`
111: @*/
112: PetscErrorCode KSPCGSetObjectiveTarget(KSP ksp, PetscReal obj)
113: {
114:   PetscFunctionBegin;
117:   PetscTryMethod(ksp, "KSPCGSetObjectiveTarget_C", (KSP, PetscReal), (ksp, obj));
118:   PetscFunctionReturn(PETSC_SUCCESS);
119: }
121: /*@
122:   KSPCGGetNormD - Get norm of the direction when the solver is used inside `SNESNEWTONTR`
124:   Not collective
126:   Input Parameters:
127: + ksp    - the iterative context
128: - norm_d - the norm of the direction
130:   Level: advanced
132: .seealso: [](ch_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`, `SNESNEWTONTR`
133: @*/
134: PetscErrorCode KSPCGGetNormD(KSP ksp, PetscReal *norm_d)
135: {
136:   PetscFunctionBegin;
138:   PetscAssertPointer(norm_d, 2);
139:   PetscUseMethod(ksp, "KSPCGGetNormD_C", (KSP, PetscReal *), (ksp, norm_d));
140:   PetscFunctionReturn(PETSC_SUCCESS);
141: }
143: /*@
144:   KSPCGGetObjFcn - Get the conjugate gradient objective function value
146:   Not collective
148:   Input Parameters:
149: + ksp   - the iterative context
150: - o_fcn - the objective function value
152:   Level: advanced
154:   Note:
155:   This function will return the current objective function value $\frac{1}{2} x^H_k A x_k - b^H x_k$.
156:   if called during `KSPSolve()` (e.g. during a monitor call).
157:   When called outside of a `KSPSolve()`, it will return the last computed value inside the solver.
159: .seealso: [](ch_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`, `KSPMonitorSet`
160: @*/
161: PetscErrorCode KSPCGGetObjFcn(KSP ksp, PetscReal *o_fcn)
162: {
163:   PetscFunctionBegin;
165:   PetscAssertPointer(o_fcn, 2);
166:   PetscUseMethod(ksp, "KSPCGGetObjFcn_C", (KSP, PetscReal *), (ksp, o_fcn));
167:   PetscFunctionReturn(PETSC_SUCCESS);
168: }