Actual source code: bddc.h
petsc-3.4.5 2014-06-29
5: #include <../src/ksp/pc/impls/is/pcis.h>
7: /* BDDC requires metis 5.0.1 for multilevel */
8: #include "metis.h"
9: #define MetisInt idx_t
10: #define MetisScalar real_t
13: /* Structure for graph partitioning (adapted from Metis) */
14: struct _PCBDDCGraph {
15: PetscInt nvtxs;
16: PetscInt ncmps;
17: PetscInt *xadj;
18: PetscInt *adjncy;
19: PetscInt **neighbours_set;
20: PetscInt *where;
21: PetscInt *which_dof;
22: PetscInt *cptr;
23: PetscInt *queue;
24: PetscInt *count;
25: PetscMPIInt *where_ncmps;
26: PetscBool *touched;
27: };
29: typedef enum {SCATTERS_BDDC,GATHERS_BDDC} CoarseCommunicationsType;
30: typedef struct _PCBDDCGraph *PCBDDCGraph;
32: /* Private context (data structure) for the BDDC preconditioner. */
33: typedef struct {
34: /* First MUST come the folowing line, for the stuff that is common to FETI and Neumann-Neumann. */
35: PC_IS pcis;
37: /* Coarse stuffs needed by BDDC application in KSP */
38: PetscInt coarse_size;
39: Mat coarse_mat;
40: Vec coarse_vec;
41: Vec coarse_rhs;
42: KSP coarse_ksp;
43: Mat coarse_phi_B;
44: Mat coarse_phi_D;
45: PetscInt local_primal_size;
46: PetscInt *local_primal_indices;
47: PetscMPIInt *local_primal_displacements;
48: PetscMPIInt *local_primal_sizes;
49: PetscMPIInt replicated_primal_size;
50: PetscMPIInt *replicated_local_primal_indices;
51: PetscScalar *replicated_local_primal_values;
52: VecScatter coarse_loc_to_glob;
54: /* Local stuffs needed by BDDC application in KSP */
55: Vec vec1_P;
56: Vec vec1_C;
57: Mat local_auxmat1;
58: Mat local_auxmat2;
59: Vec vec1_R;
60: Vec vec2_R;
61: VecScatter R_to_B;
62: VecScatter R_to_D;
63: KSP ksp_R;
64: KSP ksp_D;
65: Vec vec4_D;
67: /* Quantities defining constraining details (local) of the preconditioner */
68: /* These quantities define the preconditioner itself */
69: IS *ISForFaces;
70: IS *ISForEdges;
71: IS ISForVertices;
72: PetscInt n_ISForFaces;
73: PetscInt n_ISForEdges;
74: PetscInt n_constraints;
75: PetscInt n_vertices;
76: Mat ConstraintMatrix;
77: PetscBool usechangeofbasis;
78: PetscBool usechangeonfaces;
79: Mat ChangeOfBasisMatrix;
80: Vec original_rhs;
81: Vec temp_solution;
82: Mat local_mat;
83: PetscBool use_exact_dirichlet;
85: /* Some defaults on selecting vertices and constraints*/
86: PetscBool vertices_flag;
87: PetscBool constraints_flag;
88: PetscBool faces_flag;
89: PetscBool edges_flag;
91: /* Some customization is possible */
92: PCBDDCGraph mat_graph;
93: MatNullSpace NullSpace;
94: MatNullSpace CoarseNullSpace;
95: PetscBool use_nnsp_true;
96: PetscInt n_ISForDofs;
97: IS *ISForDofs;
98: IS NeumannBoundaries;
99: IS DirichletBoundaries;
100: PetscBool inexact_prec_type;
101: CoarseProblemType coarse_problem_type;
102: CoarseCommunicationsType coarse_communications_type;
103: PetscInt coarsening_ratio;
104: PetscInt current_level;
105: PetscInt max_levels;
107: /* For verbose output of some bddc data structures */
108: PetscBool dbg_flag;
109: PetscViewer dbg_viewer;
110: } PC_BDDC;
112: /* prototypes for functions contained in bddc.c */
113: static PetscErrorCode PCBDDCSetUseExactDirichlet(PC,PetscBool);
114: static PetscErrorCode PCBDDCSetLevel(PC,PetscInt);
115: static PetscErrorCode PCBDDCAdaptNullSpace(PC);
116: static PetscErrorCode PCBDDCCoarseSetUp(PC);
117: static PetscErrorCode PCBDDCFindConnectedComponents(PCBDDCGraph,PetscInt);
118: static PetscErrorCode PCBDDCSetupCoarseEnvironment(PC,PetscScalar*);
119: static PetscErrorCode PCBDDCManageLocalBoundaries(PC);
120: static PetscErrorCode PCBDDCApplyInterfacePreconditioner(PC);
121: static PetscErrorCode PCBDDCSolveSaddlePoint(PC);
122: static PetscErrorCode PCBDDCScatterCoarseDataBegin(PC,Vec,Vec,InsertMode,ScatterMode);
123: static PetscErrorCode PCBDDCScatterCoarseDataEnd(PC,Vec,Vec,InsertMode,ScatterMode);
124: static PetscErrorCode PCBDDCCreateConstraintMatrix(PC);
126: /* feti-dp */
127: typedef struct {
128: PetscInt n_lambda;
129: Vec lambda_local;
130: Vec temp_solution_B;
131: Vec temp_solution_D;
132: Mat B_delta;
133: Mat B_Ddelta;
134: VecScatter l2g_lambda;
135: PC pc;
136: } FETIDPMat_ctx;
138: typedef struct {
139: Vec lambda_local;
140: Mat B_Ddelta;
141: VecScatter l2g_lambda;
142: PC pc;
143: } FETIDPPC_ctx;
145: /* inexact solvers with nullspace correction */
146: typedef struct {
147: Mat basis_mat;
148: Mat Kbasis_mat;
149: Mat Lbasis_mat;
150: PC local_pc;
151: Vec work_small_1;
152: Vec work_small_2;
153: Vec work_full_1;
154: Vec work_full_2;
155: } NullSpaceCorrection_ctx;
157: static PetscErrorCode PCBDDCAdaptLocalProblem(PC,IS);
158: static PetscErrorCode PCBDDCDestroyNullSpaceCorrectionPC(PC);
159: static PetscErrorCode PCBDDCApplyNullSpaceCorrectionPC(PC,Vec,Vec);
160: static PetscErrorCode PCBDDCCreateFETIDPMatContext(PC,FETIDPMat_ctx**);
161: static PetscErrorCode PCBDDCDestroyFETIDPMat(Mat);
162: static PetscErrorCode PCBDDCSetupFETIDPMatContext(FETIDPMat_ctx*);
163: static PetscErrorCode PCBDDCCreateFETIDPPCContext(PC,FETIDPPC_ctx**);
164: static PetscErrorCode PCBDDCDestroyFETIDPPC(PC);
165: static PetscErrorCode PCBDDCSetupFETIDPPCContext(Mat,FETIDPPC_ctx*);
166: static PetscErrorCode FETIDPPCApply(PC,Vec,Vec);
167: static PetscErrorCode FETIDPMatMult(Mat,Vec,Vec);
169: #endif /* __pcbddc_h */