Actual source code: pcbddcimpl.h

  1: #pragma once

  3: #include <petsc/private/pcisimpl.h>
  4: #include <petsc/private/pcbddcstructsimpl.h>

  6: #if !defined(PETSC_PCBDDC_MAXLEVELS)
  7:   #define PETSC_PCBDDC_MAXLEVELS 8
  8: #endif

 10: PETSC_EXTERN PetscLogEvent PC_BDDC_Topology[PETSC_PCBDDC_MAXLEVELS];
 11: PETSC_EXTERN PetscLogEvent PC_BDDC_LocalSolvers[PETSC_PCBDDC_MAXLEVELS];
 12: PETSC_EXTERN PetscLogEvent PC_BDDC_LocalWork[PETSC_PCBDDC_MAXLEVELS];
 13: PETSC_EXTERN PetscLogEvent PC_BDDC_CorrectionSetUp[PETSC_PCBDDC_MAXLEVELS];
 14: PETSC_EXTERN PetscLogEvent PC_BDDC_CoarseSetUp[PETSC_PCBDDC_MAXLEVELS];
 15: PETSC_EXTERN PetscLogEvent PC_BDDC_ApproxSetUp[PETSC_PCBDDC_MAXLEVELS];
 16: PETSC_EXTERN PetscLogEvent PC_BDDC_ApproxApply[PETSC_PCBDDC_MAXLEVELS];
 17: PETSC_EXTERN PetscLogEvent PC_BDDC_CoarseSolver[PETSC_PCBDDC_MAXLEVELS];
 18: PETSC_EXTERN PetscLogEvent PC_BDDC_AdaptiveSetUp[PETSC_PCBDDC_MAXLEVELS];
 19: PETSC_EXTERN PetscLogEvent PC_BDDC_Scaling[PETSC_PCBDDC_MAXLEVELS];
 20: PETSC_EXTERN PetscLogEvent PC_BDDC_Schurs[PETSC_PCBDDC_MAXLEVELS];
 21: PETSC_EXTERN PetscLogEvent PC_BDDC_Solves[PETSC_PCBDDC_MAXLEVELS][3];

 23: /* Private context (data structure) for the BDDC preconditioner.  */
 24: typedef struct {
 25:   /* First MUST come the following line, for the stuff that is common to FETI and Neumann-Neumann. */
 26:   PC_IS pcis;
 27:   /* Coarse stuffs needed by BDDC application in KSP */
 28:   Vec        coarse_vec;
 29:   KSP        coarse_ksp;
 30:   Mat        coarse_phi_B;
 31:   Mat        coarse_phi_D;
 32:   Mat        coarse_psi_B;
 33:   Mat        coarse_psi_D;
 34:   PetscInt   local_primal_size;
 35:   PetscInt   coarse_size;
 36:   PetscInt  *global_primal_indices;
 37:   VecScatter coarse_loc_to_glob;
 38:   /* Local stuffs needed by BDDC application in KSP */
 39:   Vec        vec1_P;
 40:   Vec        vec1_C;
 41:   Mat        local_auxmat1;
 42:   Mat        local_auxmat2;
 43:   Vec        vec1_R;
 44:   Vec        vec2_R;
 45:   IS         is_R_local;
 46:   VecScatter R_to_B;
 47:   VecScatter R_to_D;
 48:   KSP        ksp_R;
 49:   KSP        ksp_D;
 50:   /* Quantities defining constraining details (local) of the preconditioner */
 51:   /* These quantities define the preconditioner itself */
 52:   PetscInt  n_vertices;
 53:   Mat       ConstraintMatrix;
 54:   PetscBool new_primal_space;
 55:   PetscBool new_primal_space_local;
 56:   PetscInt *primal_indices_local_idxs;
 57:   PetscInt  local_primal_size_cc;
 58:   PetscInt *local_primal_ref_node;
 59:   PetscInt *local_primal_ref_mult;
 60:   PetscBool use_change_of_basis;
 61:   PetscBool use_change_on_faces;
 62:   PetscBool fake_change;
 63:   Mat       ChangeOfBasisMatrix;
 64:   Mat       user_ChangeOfBasisMatrix;
 65:   PetscBool change_interior;
 66:   Mat       switch_static_change;
 67:   Vec       work_change;
 68:   Vec       original_rhs;
 69:   Vec       temp_solution;
 70:   Mat       local_mat;
 71:   PetscBool use_exact_dirichlet_trick;
 72:   PetscBool exact_dirichlet_trick_app;
 73:   PetscBool ksp_guess_nonzero;
 74:   PetscBool rhs_change;
 75:   PetscBool temp_solution_used;
 76:   /* benign subspace trick */
 77:   PetscBool    benign_saddle_point;
 78:   PetscBool    benign_have_null;
 79:   PetscBool    benign_skip_correction;
 80:   PetscBool    benign_compute_correction;
 81:   Mat          benign_change;
 82:   Mat          benign_original_mat;
 83:   IS          *benign_zerodiag_subs;
 84:   Vec          benign_vec;
 85:   Mat          benign_B0;
 86:   PetscSF      benign_sf;
 87:   PetscScalar *benign_p0;
 88:   PetscInt     benign_n;
 89:   PetscInt    *benign_p0_lidx;
 90:   PetscInt    *benign_p0_gidx;
 91:   PetscBool    benign_null;
 92:   PetscBool    benign_change_explicit;
 93:   PetscBool    benign_apply_coarse_only;

 95:   /* Some defaults on selecting vertices and constraints*/
 96:   PetscBool use_local_adj;
 97:   PetscInt  local_adj_square;
 98:   PetscBool use_vertices;
 99:   PetscBool use_faces;
100:   PetscBool use_edges;

102:   /* Some customization is possible */
103:   PetscBool              corner_selection;
104:   PetscBool              corner_selected;
105:   PetscBool              recompute_topography;
106:   PetscBool              graphanalyzed;
107:   PCBDDCGraph            mat_graph;
108:   PetscInt               graphmaxcount;
109:   MatNullSpace           onearnullspace;
110:   PetscObjectState      *onearnullvecs_state;
111:   PetscBool              NullSpace_corr[4];
112:   IS                     user_primal_vertices;
113:   IS                     user_primal_vertices_local;
114:   PetscBool              use_nnsp;
115:   PetscBool              use_nnsp_true;
116:   PetscBool              use_qr_single;
117:   PetscBool              user_provided_isfordofs;
118:   PetscInt               n_ISForDofs;
119:   PetscInt               n_ISForDofsLocal;
120:   IS                    *ISForDofs;
121:   IS                    *ISForDofsLocal;
122:   IS                     NeumannBoundaries;
123:   IS                     NeumannBoundariesLocal;
124:   IS                     DirichletBoundaries;
125:   IS                     DirichletBoundariesLocal;
126:   PetscBool              eliminate_dirdofs;
127:   PetscBool              switch_static;
128:   PetscInt               coarsening_ratio;
129:   PetscInt               coarse_adj_red;
130:   PetscInt               current_level;
131:   PetscInt               max_levels;
132:   PetscInt               coarse_eqs_per_proc;
133:   PetscInt               coarse_eqs_limit;
134:   IS                     coarse_subassembling;
135:   PetscBool              use_coarse_estimates;
136:   PetscBool              symmetric_primal;
137:   PetscInt               vertex_size;
138:   PCBDDCInterfaceExtType interface_extension;

140:   /* no-net-flux */
141:   PetscBool compute_nonetflux;
142:   Mat       divudotp;
143:   PetscBool divudotp_trans;
144:   IS        divudotp_vl2l;

146:   /* nedelec */
147:   Mat       discretegradient;
148:   PetscInt  nedorder;
149:   PetscBool conforming;
150:   PetscInt  nedfield;
151:   PetscBool nedglobal;
152:   Mat       nedcG;
153:   IS        nedclocal;

155:   /* local disconnected subdomains */
156:   PetscBool detect_disconnected;
157:   PetscBool detect_disconnected_filter;
158:   PetscInt  n_local_subs;
159:   IS       *local_subs;

161:   /* scaling */
162:   Vec                 work_scaling;
163:   PetscBool           use_deluxe_scaling;
164:   PCBDDCDeluxeScaling deluxe_ctx;
165:   PetscBool           deluxe_zerorows;
166:   PetscBool           deluxe_singlemat;

168:   /* schur complements on interface's subsets */
169:   PCBDDCSubSchurs sub_schurs;
170:   PetscBool       sub_schurs_rebuild;
171:   PetscBool       sub_schurs_exact_schur;
172:   PetscInt        sub_schurs_layers;
173:   PetscBool       sub_schurs_use_useradj;
174:   PetscBool       computed_rowadj;

176:   /* adaptive selection of constraints */
177:   PetscBool    adaptive_selection;
178:   PetscBool    adaptive_userdefined;
179:   PetscReal    adaptive_threshold[2];
180:   PetscInt     adaptive_nmin;
181:   PetscInt     adaptive_nmax;
182:   PetscInt    *adaptive_constraints_n;
183:   PetscInt    *adaptive_constraints_idxs;
184:   PetscInt    *adaptive_constraints_idxs_ptr;
185:   PetscScalar *adaptive_constraints_data;
186:   PetscInt    *adaptive_constraints_data_ptr;

188:   /* For verbose output of some bddc data structures */
189:   PetscInt    dbg_flag;
190:   PetscViewer dbg_viewer;
191: } PC_BDDC;