14b9ad928SBarry Smith 24b9ad928SBarry Smith #if !defined(__BJACOBI_H) 34b9ad928SBarry Smith #define __BJACOBI_H 44b9ad928SBarry Smith /* 54b9ad928SBarry Smith Private data for block Jacobi and block Gauss-Seidel preconditioner. 64b9ad928SBarry Smith */ 7c6db04a5SJed Brown #include <petscksp.h> 8b45d2f2cSJed Brown #include <petsc-private/pcimpl.h> 94b9ad928SBarry Smith 104b9ad928SBarry Smith /* 114b9ad928SBarry Smith This data is general for all implementations 124b9ad928SBarry Smith */ 134b9ad928SBarry Smith typedef struct { 145a7e1895SHong Zhang PetscInt n; /* number of global blocks */ 155a7e1895SHong Zhang PetscInt n_local; /* number of blocks in this subcommunicator or in this process */ 1613f74950SBarry Smith PetscInt first_local; /* number of first block on processor */ 17ace3abfcSBarry Smith PetscBool use_true_local; /* use block from true matrix, not preconditioner matrix for local MatMult() */ 182bbb50dfSHong Zhang KSP *ksp; /* KSP contexts for blocks or for subcommunicator */ 194b9ad928SBarry Smith void *data; /* implementation-specific data */ 20ace3abfcSBarry Smith PetscBool same_local_solves; /* flag indicating whether all local solvers are same (used for PCView()) */ 2113f74950SBarry Smith PetscInt *l_lens; /* lens of each block */ 2213f74950SBarry Smith PetscInt *g_lens; 235a7e1895SHong Zhang PetscSubcomm psubcomm; /* for multiple processors per block */ 244b9ad928SBarry Smith } PC_BJacobi; 254b9ad928SBarry Smith 264b9ad928SBarry Smith /* 274b9ad928SBarry Smith This data is specific for certain implementations 284b9ad928SBarry Smith */ 294b9ad928SBarry Smith 304b9ad928SBarry Smith /* This is for multiple blocks per processor */ 314b9ad928SBarry Smith typedef struct { 324b9ad928SBarry Smith Vec *x,*y; /* work vectors for solves on each block */ 3313f74950SBarry Smith PetscInt *starts; /* starting point of each block */ 344b9ad928SBarry Smith Mat *mat,*pmat; /* submatrices for each block */ 354b9ad928SBarry Smith IS *is; /* for gathering the submatrices */ 364b9ad928SBarry Smith } PC_BJacobi_Multiblock; 374b9ad928SBarry Smith 384b9ad928SBarry Smith /* This is for a single block per processor */ 394b9ad928SBarry Smith typedef struct { 404b9ad928SBarry Smith Vec x,y; 414b9ad928SBarry Smith } PC_BJacobi_Singleblock; 424b9ad928SBarry Smith 435a7e1895SHong Zhang /* This is for multiple processors per block */ 445a7e1895SHong Zhang typedef struct { 455a7e1895SHong Zhang PC pc; /* preconditioner used on each subcommunicator */ 46*ce94432eSBarry Smith Vec xsub,ysub; /* vectors of a subcommunicator to hold parallel vectors of PetscObjectComm((PetscObject)pc) */ 475a7e1895SHong Zhang Mat submats; /* matrix and optional preconditioner matrix belong to a subcommunicator */ 48d6037b41SHong Zhang PetscSubcomm psubcomm; 495a7e1895SHong Zhang } PC_BJacobi_Multiproc; 504b9ad928SBarry Smith #endif 514b9ad928SBarry Smith 524b9ad928SBarry Smith 53