5 #include <amdis/linearalgebra/Constraints.hpp> 6 #include <amdis/linearalgebra/petsc/MatrixBackend.hpp> 7 #include <amdis/linearalgebra/petsc/VectorBackend.hpp> 11 template <
class DofMap>
25 template <
class BitVector>
28 thread_local std::vector<PetscInt> rows;
30 auto const& dofMap = mat.dofMap_;
31 for (std::size_t i = 0; i < nodes.size(); ++i)
33 rows.push_back(dofMap.global(i));
37 MatIsSymmetricKnown(mat.
matrix(), &
set, &flg);
39 if (
set == PETSC_TRUE && flg == PETSC_TRUE)
40 MatZeroRowsColumns(mat.
matrix(), rows.size(), rows.data(), setDiagonal ? 1.0 : 0.0, x.
vector(), b.
vector());
42 MatZeroRows(mat.
matrix(), rows.size(), rows.data(), setDiagonal ? 1.0 : 0.0, x.
vector(), b.
vector());
50 template <
class BitVector,
class Associations>
51 static void periodicBC(Matrix& mat,
Vector& x,
Vector& b, BitVector
const& left, Associations
const& left2right,
bool setDiagonal =
true)
53 thread_local std::vector<PetscInt> rows, associated_rows;
55 associated_rows.clear();
56 auto const& dofMap = mat.dofMap_;
57 for (PetscInt i = 0; i < PetscInt(left.size()); ++i) {
60 PetscInt row_local[2] = {i, at(left2right,i)};
61 PetscInt row[2] = {dofMap.global(row_local[0]), dofMap.global(row_local[1])};
63 rows.push_back(row[0]);
64 associated_rows.push_back(row[1]);
69 PetscScalar
const* vals;
70 MatGetRow(mat.
matrix(), row[0], &ncols, &cols, &vals);
71 MatSetValues(mat.
matrix(), 1, &row[1], ncols, cols, vals, ADD_VALUES);
72 MatRestoreRow(mat.
matrix(), row[0], &ncols, &cols, &vals);
76 MatAssemblyBegin(mat.
matrix(), MAT_FLUSH_ASSEMBLY);
77 MatAssemblyEnd(mat.
matrix(), MAT_FLUSH_ASSEMBLY);
80 MatSetOption(mat.
matrix(), MAT_NO_OFF_PROC_ZERO_ROWS, PETSC_TRUE);
81 MatZeroRows(mat.
matrix(), rows.size(), rows.data(), setDiagonal ? 1 : 0, PETSC_NULL, PETSC_NULL);
85 for (std::size_t i = 0; i < rows.size(); ++i)
86 MatSetValue(mat.
matrix(), rows[i], associated_rows[i], -1, ADD_VALUES);
87 MatAssemblyBegin(mat.
matrix(), MAT_FINAL_ASSEMBLY);
88 MatAssemblyEnd(mat.
matrix(), MAT_FINAL_ASSEMBLY);
90 std::vector<PetscScalar> data(rows.size());
93 VecGetValues(b.
vector(), rows.size(), rows.data(), data.data());
94 VecSetValues(b.
vector(), rows.size(), associated_rows.data(), data.data(), ADD_VALUES);
95 VecAssemblyBegin(b.
vector());
96 VecAssemblyEnd(b.
vector());
99 data.assign(rows.size(), 0);
100 VecSetValues(b.
vector(), rows.size(), rows.data(), data.data(), INSERT_VALUES);
101 VecAssemblyBegin(b.
vector());
102 VecAssemblyEnd(b.
vector());
105 VecGetValues(x.
vector(), rows.size(), rows.data(), data.data());
106 VecSetValues(x.
vector(), rows.size(), associated_rows.data(), data.data(), INSERT_VALUES);
107 VecAssemblyBegin(x.
vector());
108 VecAssemblyEnd(x.
vector());
112 template <
class Associations>
113 static PetscInt at(Associations
const& m, std::size_t idx)
115 auto it = m.find(idx);
116 assert(it != m.end());
BaseVector const & vector() const
Return the data-vector vector_.
Definition: VectorBackend.hpp:81
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
static void periodicBC(Matrix &mat, Vector &x, Vector &b, BitVector const &left, Associations const &left2right, bool setDiagonal=true)
Periodic boundary conditions.
Definition: Constraints.hpp:51
The basic container that stores a base matrix.
Definition: MatrixBackend.hpp:20
BaseMatrix & matrix()
Return a reference to the data-matrix matrix.
Definition: MatrixBackend.hpp:53
Definition: Constraints.hpp:13
static void dirichletBC(Matrix &mat, Vector &x, Vector &b, BitVector const &nodes, bool setDiagonal=true)
Dirichlet boundary conditions.
Definition: Constraints.hpp:26
The basic container that stores a base vector data.
Definition: VectorBackend.hpp:19