2: #include <../src/ksp/ksp/impls/cg/cgimpl.h> /*I "petscksp.h" I*/
6: /*@
7: KSPCGSetType - Sets the variant of the conjugate gradient method to
8: use for solving a linear system with a complex coefficient matrix.
9: This option is irrelevant when solving a real system.
11: Logically Collective on KSP 13: Input Parameters:
14: + ksp - the iterative context
15: - type - the variant of CG to use, one of
16: .vb
17: KSP_CG_HERMITIAN - complex, Hermitian matrix (default)
18: KSP_CG_SYMMETRIC - complex, symmetric matrix
19: .ve
21: Level: intermediate
22: 23: Options Database Keys:
24: + -ksp_cg_Hermitian - Indicates Hermitian matrix
25: - -ksp_cg_symmetric - Indicates symmetric matrix
27: Note:
28: By default, the matrix is assumed to be complex, Hermitian.
30: .keywords: CG, conjugate gradient, Hermitian, symmetric, set, type
31: @*/
32: PetscErrorCodeKSPCGSetType(KSP ksp,KSPCGType type) 33: {
38: PetscTryMethod(ksp,"KSPCGSetType_C",(KSP,KSPCGType),(ksp,type));
39: return(0);
40: }
44: /*@
45: KSPCGUseSingleReduction - Merge the two inner products needed in CG into a single MPI_Allreduce() call.
47: Logically Collective on KSP 49: Input Parameters:
50: + ksp - the iterative context
51: - flg - turn on or off the single reduction
53: Options Database:
54: . -ksp_cg_single_reduction
56: Level: intermediate
58: The algorithm used in this case is described as Method 1 in Lapack Working Note 56, "Conjugate Gradient Algorithms with Reduced Synchronization Overhead
59: Distributed Memory Multiprocessors", by E. F. D'Azevedo, V. L. Eijkhout, and C. H. Romine, December 3, 1999. V. Eijkhout creates the algorithm
60: initially to Chronopoulos and Gear.
62: It requires two extra work vectors than the conventional implementation in PETSc.
63: 64: .keywords: CG, conjugate gradient, Hermitian, symmetric, set, type
65: @*/
66: PetscErrorCodeKSPCGUseSingleReduction(KSP ksp,PetscBool flg) 67: {
73: PetscTryMethod(ksp,"KSPCGUseSingleReduction_C",(KSP,PetscBool),(ksp,flg));
74: return(0);
75: }