1 2 #if !defined(__BJACOBI_H) 3 #define __BJACOBI_H 4 /* 5 Private data for block Jacobi and block Gauss-Seidel preconditioner. 6 */ 7 #include "petscksp.h" 8 #include "private/pcimpl.h" 9 10 /* 11 This data is general for all implementations 12 */ 13 typedef struct { 14 PetscInt n; /* number of global blocks */ 15 PetscInt n_local; /* number of blocks in this subcommunicator or in this process */ 16 PetscInt first_local; /* number of first block on processor */ 17 PetscBool use_true_local; /* use block from true matrix, not preconditioner matrix for local MatMult() */ 18 KSP *ksp; /* KSP contexts for blocks */ 19 void *data; /* implementation-specific data */ 20 PetscBool same_local_solves; /* flag indicating whether all local solvers are same (used for PCView()) */ 21 PetscInt *l_lens; /* lens of each block */ 22 PetscInt *g_lens; 23 Mat tp_mat,tp_pmat; /* diagonal block of matrix for this processor */ 24 PetscSubcomm psubcomm; /* for multiple processors per block */ 25 } PC_BJacobi; 26 27 /* 28 This data is specific for certain implementations 29 */ 30 31 /* This is for multiple blocks per processor */ 32 typedef struct { 33 Vec *x,*y; /* work vectors for solves on each block */ 34 PetscInt *starts; /* starting point of each block */ 35 Mat *mat,*pmat; /* submatrices for each block */ 36 IS *is; /* for gathering the submatrices */ 37 } PC_BJacobi_Multiblock; 38 39 /* This is for a single block per processor */ 40 typedef struct { 41 Vec x,y; 42 } PC_BJacobi_Singleblock; 43 44 /* This is for multiple processors per block */ 45 typedef struct { 46 KSP ksp; /* ksp used on each subcommunicator */ 47 PC pc; /* preconditioner used on each subcommunicator */ 48 Vec xsub,ysub; /* vectors of a subcommunicator to hold parallel vectors of ((PetscObject)pc)->comm */ 49 Mat submats; /* matrix and optional preconditioner matrix belong to a subcommunicator */ 50 PetscSubcomm psubcomm; 51 } PC_BJacobi_Multiproc; 52 #endif 53 54 55