36#include "PCTwoLevelsBlockDiagonal.hpp"
75 Mat system_matrix, dummy;
76#if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>=5)
77 KSPGetOperators(rKspObject, &system_matrix, &dummy);
80 KSPGetOperators(rKspObject, &system_matrix, &dummy, &flag);
84 MatGetSize(system_matrix, &num_rows, &num_columns);
85 assert(num_rows==num_columns);
87 PetscInt num_local_rows, num_local_columns;
88 MatGetLocalSize(system_matrix, &num_local_rows, &num_local_columns);
92 if ((num_rows%2 != 0) || (num_local_rows%2 != 0))
94 TERMINATE(
"Wrong matrix parallel layout detected in PCLDUFactorisation.");
98 unsigned subvector_num_rows = num_rows/2;
99 unsigned subvector_local_rows = num_local_rows/2;
101 unsigned subvector_num_rows_tissue = subvector_num_rows - rBathNodes.size();
102 unsigned subvector_local_rows_tissue = subvector_num_rows_tissue;
104 unsigned subvector_num_rows_bath = rBathNodes.size();
105 unsigned subvector_local_rows_bath = subvector_num_rows_bath;
118 ISCreateStride(PETSC_COMM_WORLD, num_rows/2, 0, 2, &A11_all_rows);
124 ISCreateStride(PETSC_COMM_WORLD, A22_size, 1, 2, &A22_all_rows);
128 for (
unsigned index=0; index<rBathNodes.size(); index++)
130 phi_e_bath_rows[index] = 2*rBathNodes[index] + 1;
132#if (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR >= 2)
133 ISCreateGeneral(PETSC_COMM_WORLD, rBathNodes.size(), phi_e_bath_rows, PETSC_USE_POINTER, &A22_bath_rows);
135 ISCreateGeneralWithArray(PETSC_COMM_WORLD, rBathNodes.size(), phi_e_bath_rows, &A22_bath_rows);
139 ISDifference(A22_all_rows, A22_bath_rows, &A22_tissue_rows);
147 IS& A11_rows=A11_all_rows;
148 IS& A22_B1_rows=A22_tissue_rows;
149 IS& A22_B2_rows=A22_bath_rows;
152 ISCreateStride(PETSC_COMM_WORLD, num_rows/2, 0, 1, &all_vector);
155 ISCreateStride(PETSC_COMM_WORLD, (num_rows/2)-rBathNodes.size(), 0, 1, &tissue_vector);
158 ISCreateStride(PETSC_COMM_WORLD, rBathNodes.size(), 0, 1, &bath_vector);
177 assert(global_size == num_rows/2);
180 IS& A11_columns=A11_all_rows;
181 ISCreateStride(PETSC_COMM_WORLD, high-low, 2*low, 2, &A11_local_rows);
183#if (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR >= 8)
184 MatCreateSubMatrix(system_matrix, A11_local_rows, A11_columns,
187 MatGetSubMatrix(system_matrix, A11_local_rows, A11_columns,
203 IS& A22_B1_local_rows = A22_tissue_rows;
204 IS& A22_B1_columns = A22_tissue_rows;
206#if (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR >= 8)
207 MatCreateSubMatrix(system_matrix, A22_B1_local_rows, A22_B1_columns,
210 MatGetSubMatrix(system_matrix, A22_B1_local_rows, A22_B1_columns,
225 IS& A22_B2_local_rows = A22_bath_rows;
226 IS& A22_B2_columns = A22_bath_rows;
228#if (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR >= 8)
229 MatCreateSubMatrix(system_matrix, A22_B2_local_rows, A22_B2_columns,
232 MatGetSubMatrix(system_matrix, A22_B2_local_rows, A22_B2_columns,
240 delete[] phi_e_bath_rows;
245#if (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 2)
266#if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>=5)
277#if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>=5)
295#if (PETSC_VERSION_MAJOR==3 && PETSC_VERSION_MINOR>=5)
304#if (PETSC_VERSION_MAJOR == 3 && PETSC_VERSION_MINOR >= 1)
305PetscErrorCode PCTwoLevelsBlockDiagonalApply(PC pc_object,
Vec x,
Vec y)
309 PCShellGetContext(pc_object, &pc_context);
311PetscErrorCode PCTwoLevelsBlockDiagonalApply(
void* pc_context,
Vec x,
Vec y)
317 assert(block_diag_context!=
nullptr);
323#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
332#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
341#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
362#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
371#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
380#if ((PETSC_VERSION_MAJOR == 3) || (PETSC_VERSION_MAJOR == 2 && PETSC_VERSION_MINOR == 3 && PETSC_VERSION_SUBMINOR == 3))
#define TERMINATE(message)
void PCTwoLevelsBlockDiagonalCreate(KSP &rKspObject, std::vector< PetscInt > &rBathNodes)
void PCTwoLevelsBlockDiagonalSetUp()
PCTwoLevelsBlockDiagonalContext mPCContext
~PCTwoLevelsBlockDiagonal()
PCTwoLevelsBlockDiagonal(KSP &rKspObject, std::vector< PetscInt > &rBathNodes)
VecScatter A22_B2_scatter_ctx
Mat A22_B1_matrix_subblock
VecScatter A22_B1_scatter_ctx
VecScatter A11_scatter_ctx
Mat A22_B2_matrix_subblock