4: #include <petscksp.h> 5: #include <petscbt.h> 7: /* special marks for interface graph: they cannot be enums, since special marks should in principle range from -4 to -max_int */ 8: #define PCBDDCGRAPH_NEUMANN_MARK -1 9: #define PCBDDCGRAPH_DIRICHLET_MARK -2 10: #define PCBDDCGRAPH_LOCAL_PERIODIC_MARK -3 11: #define PCBDDCGRAPH_SPECIAL_MARK -4 13: /* Structure for local graph partitioning */ 14: struct _PCBDDCGraph { 15: PetscBool setupcalled; 16: /* graph information */ 17: ISLocalToGlobalMapping l2gmap; 18: PetscInt nvtxs; 19: PetscInt nvtxs_global; 20: PetscBT touched; 21: PetscInt *count; 22: PetscInt **neighbours_set; 23: PetscInt *subset; 24: PetscInt *which_dof; 25: PetscInt *special_dof; 26: PetscInt custom_minimal_size; 27: PetscBool twodim; 28: PetscBool twodimset; 29: PetscBool has_dirichlet; 30: IS dirdofs; 31: IS dirdofsB; 32: PetscInt commsizelimit; 33: PetscInt maxcount; 34: /* data for connected components */ 35: PetscInt ncc; 36: PetscInt *cptr; 37: PetscInt *queue; 38: PetscBool queue_sorted; 39: /* data for interface subsets */ 40: PetscInt n_subsets; 41: PetscInt *subset_size; 42: PetscInt **subset_idxs; 43: PetscInt *subset_ncc; 44: PetscInt *subset_ref_node; 45: /* data for periodic dofs */ 46: PetscInt *mirrors; 47: PetscInt **mirrors_set; 48: /* placeholders for connectivity relation between dofs */ 49: PetscInt nvtxs_csr; 50: PetscInt *xadj; 51: PetscInt *adjncy; 52: PetscBool freecsr; 53: /* data for local subdomains (if any have been detected) 54: these are not intended to be exposed */ 55: PetscInt n_local_subs; 56: PetscInt *local_subs; 57: }; 58: typedef struct _PCBDDCGraph *PCBDDCGraph; 60: /* Wrap to MatFactor solver in Schur complement mode. Provides 61: - standalone solver for interior variables 62: - forward and backward substitutions for correction solver 63: */ 64: /* It assumes that interior variables are a contiguous set starting from 0 */ 65: struct _PCBDDCReuseSolvers { 66: /* the factored matrix obtained from MatGetFactor(...,solver_package,...) */ 67: Mat F; 68: /* placeholders for the solution and rhs on the whole set of dofs of A (size local_dofs - local_vertices)*/ 69: Vec sol; 70: Vec rhs; 71: /* */ 72: PetscBool has_vertices; 73: /* shell PCs to handle interior/correction solvers */ 74: PC interior_solver; 75: PC correction_solver; 76: IS is_R; 77: /* objects to handle Schur complement solution */ 78: Vec rhs_B; 79: Vec sol_B; 80: IS is_B; 81: VecScatter correction_scatter_B; 82: /* handle benign trick without change of basis on pressures */ 83: PetscInt benign_n; 84: IS *benign_zerodiag_subs; 85: PetscScalar *benign_save_vals; 86: Mat benign_csAIB; 87: Mat benign_AIIm1ones; 88: Vec benign_corr_work; 89: Vec benign_dummy_schur_vec; 90: }; 91: typedef struct _PCBDDCReuseSolvers *PCBDDCReuseSolvers; 93: /* structure to handle Schur complements on subsets */ 94: struct _PCBDDCSubSchurs { 95: /* local Neumann matrix */ 96: Mat A; 97: /* local Schur complement */ 98: Mat S; 99: /* index sets */ 100: IS is_I; 101: IS is_B; 102: /* whether Schur complements are explicitly computed with or not */ 103: PetscBool schur_explicit; 104: /* matrices cointained explicit schur complements cat together */ 105: /* note that AIJ format is used but the values are inserted as in column major ordering */ 106: Mat S_Ej_all; 107: Mat sum_S_Ej_all; 108: Mat sum_S_Ej_inv_all; 109: Mat sum_S_Ej_tilda_all; 110: IS is_Ej_all; 111: IS is_vertices; 112: IS is_dir; 113: /* l2g maps */ 114: ISLocalToGlobalMapping l2gmap; 115: ISLocalToGlobalMapping BtoNmap; 116: /* number of local subproblems */ 117: PetscInt n_subs; 118: /* connected components */ 119: IS* is_subs; 120: PetscBT is_edge; 121: /* mat flags */ 122: PetscBool is_hermitian; 123: PetscBool is_posdef; 124: /* data structure to reuse MatFactor with Schur solver */ 125: PCBDDCReuseSolvers reuse_solver; 126: /* change of variables */ 127: KSP *change; 128: IS *change_primal_sub; 129: PetscBool change_with_qr; 130: /* prefix */ 131: char *prefix; 132: }; 133: typedef struct _PCBDDCSubSchurs *PCBDDCSubSchurs; 135: /* Structure for deluxe scaling */ 136: struct _PCBDDCDeluxeScaling { 137: /* simple scaling on selected dofs (i.e. primal vertices and nodes on interface dirichlet boundaries) */ 138: PetscInt n_simple; 139: PetscInt* idx_simple_B; 140: /* handle deluxe problems */ 141: PetscInt seq_n; 142: PetscScalar *workspace; 143: VecScatter *seq_scctx; 144: Vec *seq_work1; 145: Vec *seq_work2; 146: Mat *seq_mat; 147: Mat *seq_mat_inv_sum; 148: KSP *change; 149: PetscBool change_with_qr; 150: }; 151: typedef struct _PCBDDCDeluxeScaling *PCBDDCDeluxeScaling; 153: /* inexact solvers with nullspace correction */ 154: struct _NullSpaceCorrection_ctx { 155: Mat basis_mat; 156: Mat Kbasis_mat; 157: Mat Lbasis_mat; 158: PC local_pc; 159: Vec work_small_1; 160: Vec work_small_2; 161: Vec work_full_1; 162: Vec work_full_2; 163: PetscBool apply_scaling; 164: PetscScalar scale; 165: }; 166: typedef struct _NullSpaceCorrection_ctx *NullSpaceCorrection_ctx; 168: /* MatShell context for benign mat mults */ 169: struct _PCBDDCBenignMatMult_ctx { 170: Mat A; 171: PetscInt benign_n; 172: IS *benign_zerodiag_subs; 173: PetscScalar *work; 174: PetscBool apply_left; 175: PetscBool apply_right; 176: PetscBool apply_p0; 177: PetscBool free; 178: }; 179: typedef struct _PCBDDCBenignMatMult_ctx *PCBDDCBenignMatMult_ctx; 181: /* feti-dp mat */ 182: struct _FETIDPMat_ctx { 183: PetscInt n; /* local number of rows */ 184: PetscInt N; /* global number of rows */ 185: PetscInt n_lambda; /* global number of multipliers */ 186: Vec lambda_local; 187: Vec temp_solution_B; 188: Vec temp_solution_D; 189: Mat B_delta; 190: Mat B_Ddelta; 191: PetscBool deluxe_nonred; 192: VecScatter l2g_lambda; 193: PC pc; 194: PetscBool fully_redundant; 195: /* saddle point */ 196: VecScatter l2g_lambda_only; 197: Mat B_BB; 198: Mat B_BI; 199: Mat Bt_BB; 200: Mat Bt_BI; 201: Mat C; 202: VecScatter l2g_p; 203: VecScatter g2g_p; 204: Vec vP; 205: Vec xPg; 206: Vec yPg; 207: Vec rhs_flip; 208: IS pressure; 209: IS lagrange; 210: }; 211: typedef struct _FETIDPMat_ctx *FETIDPMat_ctx; 213: /* feti-dp preconditioner */ 214: struct _FETIDPPC_ctx { 215: Mat S_j; 216: Vec lambda_local; 217: Mat B_Ddelta; 218: VecScatter l2g_lambda; 219: PC pc; 220: /* saddle point */ 221: Vec xPg; 222: Vec yPg; 223: }; 224: typedef struct _FETIDPPC_ctx *FETIDPPC_ctx; 226: struct _BDdelta_DN { 227: Mat BD; 228: KSP kBD; 229: Vec work; 230: }; 231: typedef struct _BDdelta_DN *BDdelta_DN; 233: #endif