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