1: #pragma once
3: #include <petscksp.h>
4: #include <petscbt.h>
6: /* special marks for interface graph: they cannot be enums
7: since PCBDDCGRAPH_SPECIAL_MARK ranges 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: /* coordinates (for corner detection) */
58: PetscBool active_coords;
59: PetscBool cloc;
60: PetscInt cdim, cnloc;
61: PetscReal *coords;
62: };
63: typedef struct _PCBDDCGraph *PCBDDCGraph;
65: struct _PCBDDCGraphCandidates {
66: PetscInt nfc, nec;
67: IS *Faces, *Edges, Vertices;
68: };
69: typedef struct _PCBDDCGraphCandidates *PCBDDCGraphCandidates;
71: /* Wrap to MatFactor solver in Schur complement mode. Provides
72: - standalone solver for interior variables
73: - forward and backward substitutions for correction solver
74: */
75: /* It assumes that interior variables are a contiguous set starting from 0 */
76: struct _PCBDDCReuseSolvers {
77: /* the factored matrix obtained from MatGetFactor(...,solver_package,...) */
78: Mat F;
79: /* placeholders for the solution and rhs on the whole set of dofs of A (size local_dofs - local_vertices)*/
80: Vec sol;
81: Vec rhs;
82: /* */
83: PetscBool has_vertices;
84: /* shell PCs to handle interior/correction solvers */
85: PC interior_solver;
86: PC correction_solver;
87: IS is_R;
88: /* objects to handle Schur complement solution */
89: Vec rhs_B;
90: Vec sol_B;
91: IS is_B;
92: VecScatter correction_scatter_B;
93: /* handle benign trick without change of basis on pressures */
94: PetscInt benign_n;
95: IS *benign_zerodiag_subs;
96: PetscScalar *benign_save_vals;
97: Mat benign_csAIB;
98: Mat benign_AIIm1ones;
99: Vec benign_corr_work;
100: Vec benign_dummy_schur_vec;
101: };
102: typedef struct _PCBDDCReuseSolvers *PCBDDCReuseSolvers;
104: /* structure to handle Schur complements on subsets */
105: struct _PCBDDCSubSchurs {
106: /* local Neumann matrix */
107: Mat A;
108: /* local Schur complement */
109: Mat S;
110: /* index sets */
111: IS is_I;
112: IS is_B;
113: /* whether Schur complements are explicitly computed with or not */
114: char mat_solver_type[64];
115: PetscBool schur_explicit;
116: /* BDDC or GDSW */
117: PetscBool gdsw;
118: /* matrices contained explicit schur complements cat together */
119: /* note that AIJ format is used but the values are inserted as in column major ordering */
120: Mat S_Ej_all;
121: Mat sum_S_Ej_all;
122: Mat sum_S_Ej_inv_all;
123: Mat sum_S_Ej_tilda_all;
124: IS is_Ej_all;
125: IS is_vertices;
126: IS is_dir;
127: /* l2g maps */
128: ISLocalToGlobalMapping l2gmap;
129: ISLocalToGlobalMapping BtoNmap;
130: /* number of local subproblems */
131: PetscInt n_subs;
132: /* connected components */
133: IS *is_subs;
134: PetscBT is_edge;
135: /* mat flags */
136: PetscBool is_symmetric;
137: PetscBool is_hermitian;
138: PetscBool is_posdef;
139: /* data structure to reuse MatFactor with Schur solver */
140: PCBDDCReuseSolvers reuse_solver;
141: /* change of variables */
142: KSP *change;
143: IS *change_primal_sub;
144: PetscBool change_with_qr;
145: /* prefix */
146: char *prefix;
147: /* */
148: PetscBool restrict_comm;
149: /* debug */
150: PetscBool debug;
151: };
152: typedef struct _PCBDDCSubSchurs *PCBDDCSubSchurs;
154: /* Structure for deluxe scaling */
155: struct _PCBDDCDeluxeScaling {
156: /* simple scaling on selected dofs (i.e. primal vertices and nodes on interface dirichlet boundaries) */
157: PetscInt n_simple;
158: PetscInt *idx_simple_B;
159: /* handle deluxe problems */
160: PetscInt seq_n;
161: PetscScalar *workspace;
162: VecScatter *seq_scctx;
163: Vec *seq_work1;
164: Vec *seq_work2;
165: Mat *seq_mat;
166: Mat *seq_mat_inv_sum;
167: KSP *change;
168: PetscBool change_with_qr;
169: };
170: typedef struct _PCBDDCDeluxeScaling *PCBDDCDeluxeScaling;
172: /* inexact solvers with nullspace correction */
173: struct _NullSpaceCorrection_ctx {
174: Mat basis_mat;
175: Mat inv_smat;
176: PC local_pc;
177: Vec *fw;
178: Vec *sw;
179: PetscScalar scale;
180: PetscLogEvent evapply;
181: PetscBool symm;
182: };
183: typedef struct _NullSpaceCorrection_ctx *NullSpaceCorrection_ctx;
185: /* MatShell context for benign mat mults */
186: struct _PCBDDCBenignMatMult_ctx {
187: Mat A;
188: PetscInt benign_n;
189: IS *benign_zerodiag_subs;
190: PetscScalar *work;
191: PetscBool apply_left;
192: PetscBool apply_right;
193: PetscBool apply_p0;
194: PetscBool free;
195: };
196: typedef struct _PCBDDCBenignMatMult_ctx *PCBDDCBenignMatMult_ctx;
198: /* feti-dp mat */
199: struct _FETIDPMat_ctx {
200: PetscInt n; /* local number of rows */
201: PetscInt N; /* global number of rows */
202: PetscInt n_lambda; /* global number of multipliers */
203: Vec lambda_local;
204: Vec temp_solution_B;
205: Vec temp_solution_D;
206: Mat B_delta;
207: Mat B_Ddelta;
208: PetscBool deluxe_nonred;
209: VecScatter l2g_lambda;
210: PC pc;
211: PetscBool fully_redundant;
212: /* saddle point */
213: VecScatter l2g_lambda_only;
214: Mat B_BB;
215: Mat B_BI;
216: Mat Bt_BB;
217: Mat Bt_BI;
218: Mat C;
219: VecScatter l2g_p;
220: VecScatter g2g_p;
221: Vec vP;
222: Vec xPg;
223: Vec yPg;
224: Vec rhs_flip;
225: IS pressure;
226: IS lagrange;
227: };
228: typedef struct _FETIDPMat_ctx *FETIDPMat_ctx;
230: /* feti-dp preconditioner */
231: struct _FETIDPPC_ctx {
232: Mat S_j;
233: Vec lambda_local;
234: Mat B_Ddelta;
235: VecScatter l2g_lambda;
236: PC pc;
237: /* saddle point */
238: Vec xPg;
239: Vec yPg;
240: };
241: typedef struct _FETIDPPC_ctx *FETIDPPC_ctx;
243: struct _BDdelta_DN {
244: Mat BD;
245: KSP kBD;
246: Vec work;
247: };
248: typedef struct _BDdelta_DN *BDdelta_DN;
250: /* Schur interface preconditioner */
251: struct _BDDCIPC_ctx {
252: VecScatter g2l;
253: PC bddc;
254: };
255: typedef struct _BDDCIPC_ctx *BDDCIPC_ctx;