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