Actual source code: bddc.h

petsc-3.11.4 2019-09-28
Report Typos and Errors

  4:  #include <../src/ksp/pc/impls/is/pcis.h>
  5:  #include <../src/ksp/pc/impls/bddc/bddcstructs.h>

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

 11: PETSC_EXTERN PetscLogEvent PC_BDDC_Topology[PETSC_PCBDDC_MAXLEVELS];
 12: PETSC_EXTERN PetscLogEvent PC_BDDC_LocalSolvers[PETSC_PCBDDC_MAXLEVELS];
 13: PETSC_EXTERN PetscLogEvent PC_BDDC_LocalWork[PETSC_PCBDDC_MAXLEVELS];
 14: PETSC_EXTERN PetscLogEvent PC_BDDC_CorrectionSetUp[PETSC_PCBDDC_MAXLEVELS];
 15: PETSC_EXTERN PetscLogEvent PC_BDDC_CoarseSetUp[PETSC_PCBDDC_MAXLEVELS];
 16: PETSC_EXTERN PetscLogEvent PC_BDDC_CoarseSolver[PETSC_PCBDDC_MAXLEVELS];
 17: PETSC_EXTERN PetscLogEvent PC_BDDC_AdaptiveSetUp[PETSC_PCBDDC_MAXLEVELS];
 18: PETSC_EXTERN PetscLogEvent PC_BDDC_Scaling[PETSC_PCBDDC_MAXLEVELS];
 19: PETSC_EXTERN PetscLogEvent PC_BDDC_Schurs[PETSC_PCBDDC_MAXLEVELS];

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

 93:   /* Some defaults on selecting vertices and constraints*/
 94:   PetscBool     use_local_adj;
 95:   PetscBool     use_vertices;
 96:   PetscBool     use_faces;
 97:   PetscBool     use_edges;

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

135:   /* no-net-flux */
136:   PetscBool compute_nonetflux;
137:   Mat       divudotp;
138:   PetscBool divudotp_trans;
139:   IS        divudotp_vl2l;

141:   /* nedelec */
142:   Mat       discretegradient;
143:   PetscInt  nedorder;
144:   PetscBool conforming;
145:   PetscInt  nedfield;
146:   PetscBool nedglobal;
147:   Mat       nedcG;
148:   IS        nedclocal;

150:   /* local disconnected subdomains */
151:   PetscBool detect_disconnected;
152:   PetscBool detect_disconnected_filter;
153:   PetscInt  n_local_subs;
154:   IS        *local_subs;

156:   /* scaling */
157:   Vec                 work_scaling;
158:   PetscBool           use_deluxe_scaling;
159:   PCBDDCDeluxeScaling deluxe_ctx;
160:   PetscBool           deluxe_zerorows;
161:   PetscBool           deluxe_singlemat;

163:   /* schur complements on interface's subsets */
164:   PCBDDCSubSchurs sub_schurs;
165:   PetscBool       sub_schurs_rebuild;
166:   PetscBool       sub_schurs_exact_schur;
167:   PetscInt        sub_schurs_layers;
168:   PetscBool       sub_schurs_use_useradj;
169:   PetscBool       computed_rowadj;

171:   /* adaptive selection of constraints */
172:   PetscBool    adaptive_selection;
173:   PetscBool    adaptive_userdefined;
174:   PetscReal    adaptive_threshold[2];
175:   PetscInt     adaptive_nmin;
176:   PetscInt     adaptive_nmax;
177:   PetscInt*    adaptive_constraints_n;
178:   PetscInt*    adaptive_constraints_idxs;
179:   PetscInt*    adaptive_constraints_idxs_ptr;
180:   PetscScalar* adaptive_constraints_data;
181:   PetscInt*    adaptive_constraints_data_ptr;

183:   /* For verbose output of some bddc data structures */
184:   PetscInt    dbg_flag;
185:   PetscViewer dbg_viewer;
186: } PC_BDDC;

188: #endif /* __pcbddc_h */