Actual source code: cgtype.c
2: #include <../src/ksp/ksp/impls/cg/cgimpl.h>
4: /*@
5: KSPCGSetType - Sets the variant of the conjugate gradient method to
6: use for solving a linear system with a complex coefficient matrix.
7: This option is irrelevant when solving a real system.
9: Logically Collective
11: Input Parameters:
12: + ksp - the iterative context
13: - type - the variant of CG to use, one of
14: .vb
15: KSP_CG_HERMITIAN - complex, Hermitian matrix (default)
16: KSP_CG_SYMMETRIC - complex, symmetric matrix
17: .ve
19: Level: intermediate
21: Options Database Keys:
22: + -ksp_cg_hermitian - Indicates Hermitian matrix
23: - -ksp_cg_symmetric - Indicates symmetric matrix
25: Note:
26: By default, the matrix is assumed to be complex, Hermitian.
28: .seealso: [](chapter_ksp), `KSP`, `KSPCG`
29: @*/
30: PetscErrorCode KSPCGSetType(KSP ksp, KSPCGType type)
31: {
33: PetscTryMethod(ksp, "KSPCGSetType_C", (KSP, KSPCGType), (ksp, type));
34: return 0;
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 [1]. 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: References:
59: . [1] - Lapack Working Note 56, "Conjugate Gradient Algorithms with Reduced Synchronization Overhead
60: Distributed Memory Multiprocessors", by E. F. D'Azevedo, V. L. Eijkhout, and C. H. Romine, December 3, 1999.
62: .seealso: [](chapter_ksp), [](sec_pipelineksp), `KSP`, `KSPCG`, `KSPGMRES`, `KSPPIPECG`, `KSPPIPECR`, and `KSPGROPPCG`
63: @*/
64: PetscErrorCode KSPCGUseSingleReduction(KSP ksp, PetscBool flg)
65: {
68: PetscTryMethod(ksp, "KSPCGUseSingleReduction_C", (KSP, PetscBool), (ksp, flg));
69: return 0;
70: }
72: /*@
73: KSPCGSetRadius - Sets the radius of the trust region when the solver is used inside `SNESNEWTONTR`
75: Logically Collective
77: Input Parameters:
78: + ksp - the iterative context
79: - radius - the trust region radius (Infinity is the default)
81: Level: advanced
83: .seealso: [](chapter_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`
84: @*/
85: PetscErrorCode KSPCGSetRadius(KSP ksp, PetscReal radius)
86: {
89: PetscTryMethod(ksp, "KSPCGSetRadius_C", (KSP, PetscReal), (ksp, radius));
90: return 0;
91: }
93: /*@
94: KSPCGGetNormD - Got norm of the direction when the solver is used inside `SNESNEWTONTR`
96: Collective
98: Input Parameters:
99: + ksp - the iterative context
100: - norm_d - the norm of the direction
102: Level: advanced
104: .seealso: [](chapter_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`
105: @*/
106: PetscErrorCode KSPCGGetNormD(KSP ksp, PetscReal *norm_d)
107: {
109: PetscUseMethod(ksp, "KSPCGGetNormD_C", (KSP, PetscReal *), (ksp, norm_d));
110: return 0;
111: }
113: /*@
114: KSPCGGetObjFcn - Get objective function value when the solver is used inside `SNESNEWTONTR`
116: Collective
118: Input Parameters:
119: + ksp - the iterative context
120: - o_fcn - the objective function value
122: Level: advanced
124: .seealso: [](chapter_ksp), `KSP`, `KSPCG`, `KSPNASH`, `KSPSTCG`, `KSPGLTR`
125: @*/
126: PetscErrorCode KSPCGGetObjFcn(KSP ksp, PetscReal *o_fcn)
127: {
129: PetscUseMethod(ksp, "KSPCGGetObjFcn_C", (KSP, PetscReal *), (ksp, o_fcn));
130: return 0;
131: }