Actual source code: cgtype.c

  1: #define PETSCKSP_DLL

 3:  #include ../src/ksp/ksp/impls/cg/cgimpl.h

  7: /*@
  8:     KSPCGSetType - Sets the variant of the conjugate gradient method to
  9:     use for solving a linear system with a complex coefficient matrix.
 10:     This option is irrelevant when solving a real system.

 12:     Collective on KSP

 14:     Input Parameters:
 15: +   ksp - the iterative context
 16: -   type - the variant of CG to use, one of
 17: .vb
 18:       KSP_CG_HERMITIAN - complex, Hermitian matrix (default)
 19:       KSP_CG_SYMMETRIC - complex, symmetric matrix
 20: .ve

 22:     Level: intermediate
 23:     
 24:     Options Database Keys:
 25: +   -ksp_cg_Hermitian - Indicates Hermitian matrix
 26: -   -ksp_cg_symmetric - Indicates symmetric matrix

 28:     Note:
 29:     By default, the matrix is assumed to be complex, Hermitian.

 31: .keywords: CG, conjugate gradient, Hermitian, symmetric, set, type
 32: @*/
 33: PetscErrorCode  KSPCGSetType(KSP ksp,KSPCGType type)
 34: {
 35:   PetscErrorCode ierr,(*f)(KSP,KSPCGType);

 39:   PetscObjectQueryFunction((PetscObject)ksp,"KSPCGSetType_C",(void (**)(void))&f);
 40:   if (f) {
 41:     (*f)(ksp,type);
 42:   }
 43:   return(0);
 44: }

 48: /*@
 49:     KSPCGUseSingleReduction - Merge the two inner products needed in CG into a single 
 50:      MPI_Allreduce() call.

 52:     Collective on KSP

 54:     Input Parameters:
 55: +   ksp - the iterative context
 56: -   flg - turn on or off the single reduction

 58:     Options Database:
 59: .   -ksp_cg_single_reduction 

 61:     Level: intermediate

 63:      The algorithm used in this case is described as Method 1 in Lapack Working Note 56, "Conjugate Gradient Algorithms with Reduced Synchronization Overhead 
 64:      Distributed Memory Multiprocessors", by E. F. D'Azevedo, V. L. Eijkhout, and C. H. Romine, December 3, 1999. V. Eijkhout creates the algorithm 
 65:      initially to Chronopoulos and Gear.

 67:      It requires two extra work vectors than the conventional implementation in PETSc. 
 68:     
 69: .keywords: CG, conjugate gradient, Hermitian, symmetric, set, type
 70: @*/
 71: PetscErrorCode  KSPCGUseSingleReduction(KSP ksp,PetscTruth flg)
 72: {
 73:   PetscErrorCode ierr,(*f)(KSP,PetscTruth);

 77:   PetscObjectQueryFunction((PetscObject)ksp,"KSPCGUseSingleReduction_C",(void (**)(void))&f);
 78:   if (f) {
 79:     (*f)(ksp,flg);
 80:   }
 81:   return(0);
 82: }