9 #include <dune/common/timer.hh> 11 #include <amdis/Output.hpp> 12 #include <amdis/linearalgebra/SymmetryStructure.hpp> 16 template <
class>
class Constraints;
19 template <
class DofMap>
36 template <
class Basis>
38 : dofMap_(rowBasis.communicator().dofMap())
65 template <
class RowIndex,
class ColIndex>
66 void insert(RowIndex
const& r, ColIndex
const& c, PetscScalar value)
68 PetscInt row = dofMap_.global(r);
69 PetscInt col = dofMap_.global(c);
70 MatSetValue(matrix_, row, col, value, ADD_VALUES);
74 template <
class LocalInd,
class LocalMatrix>
75 void scatter(LocalInd
const& localInd, LocalMatrix
const& localMat)
77 thread_local std::vector<PetscInt> idx;
78 idx.resize(localInd.size());
81 std::transform(localInd.begin(), localInd.end(), idx.begin(),
82 [
this](
auto const& mi) {
return dofMap_.global(mi); });
84 MatSetValues(matrix_, idx.size(), idx.data(), idx.size(), idx.data(), localMat.data(), ADD_VALUES);
88 template <
class RowLocalInd,
class ColLocalInd,
class LocalMatrix>
89 void scatter(RowLocalInd
const& rowLocalInd, ColLocalInd
const& colLocalInd, LocalMatrix
const& localMat)
91 thread_local std::vector<PetscInt> ri;
92 thread_local std::vector<PetscInt> ci;
93 ri.resize(rowLocalInd.size());
94 ci.resize(colLocalInd.size());
97 std::transform(rowLocalInd.begin(), rowLocalInd.end(), ri.begin(),
98 [
this](
auto const& mi) {
return dofMap_.global(mi); });
99 std::transform(colLocalInd.begin(), colLocalInd.end(), ci.begin(),
100 [
this](
auto const& mi) {
return dofMap_.global(mi); });
102 MatSetValues(matrix_, ri.size(), ri.data(), ci.size(), ci.data(), localMat.data(), ADD_VALUES);
106 template <
class Pattern>
107 void init(Pattern
const& pattern)
113 info(3,
" destroy old matrix needed {} seconds", t.elapsed());
118 dofMap_.localSize(), dofMap_.localSize(),
119 dofMap_.globalSize(), dofMap_.globalSize(),
120 0, pattern.d_nnz().data(), 0, pattern.o_nnz().data(), &matrix_);
123 MatSetOption(matrix_, MAT_KEEP_NONZERO_PATTERN, PETSC_TRUE);
126 switch (pattern.symmetry()) {
127 case SymmetryStructure::spd:
128 MatSetOption(matrix_, MAT_SPD, PETSC_TRUE);
130 case SymmetryStructure::symmetric:
131 MatSetOption(matrix_, MAT_SYMMETRIC, PETSC_TRUE);
133 case SymmetryStructure::hermitian:
134 MatSetOption(matrix_, MAT_HERMITIAN, PETSC_TRUE);
136 case SymmetryStructure::structurally_symmetric:
137 MatSetOption(matrix_, MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE);
144 info(3,
" create new matrix needed {} seconds", t.elapsed());
153 MatZeroEntries(matrix_);
161 MatAssemblyBegin(matrix_, MAT_FINAL_ASSEMBLY);
162 MatAssemblyEnd(matrix_, MAT_FINAL_ASSEMBLY);
163 info(3,
" finish matrix assembling needed {} seconds", t.elapsed());
170 MatGetInfo(matrix_, MAT_LOCAL, &info);
171 return std::size_t(info.nz_used);
176 MPI_Comm comm()
const 181 return PETSC_COMM_SELF;
189 MatDestroy(&matrix_);
194 DofMap
const& dofMap_;
199 bool initialized_ =
false;
void finish()
Finish assembly. Must be called before matrix can be used in a KSP.
Definition: MatrixBackend.hpp:158
BaseMatrix const & matrix() const
Return a reference to the data-matrix matrix.
Definition: MatrixBackend.hpp:59
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
void scatter(LocalInd const &localInd, LocalMatrix const &localMat)
Insert an element-matrix with row-indices == col-indices.
Definition: MatrixBackend.hpp:75
void init(Pattern const &pattern)
Create and initialize the matrix.
Definition: MatrixBackend.hpp:107
::Mat BaseMatrix
The matrix type of the underlying base matrix.
Definition: MatrixBackend.hpp:26
void scatter(RowLocalInd const &rowLocalInd, ColLocalInd const &colLocalInd, LocalMatrix const &localMat)
Insert an element-matrix.
Definition: MatrixBackend.hpp:89
The basic container that stores a base matrix.
Definition: MatrixBackend.hpp:20
PetscInt size_type
The index/size - type.
Definition: MatrixBackend.hpp:32
PetscScalar value_type
The type of the elements of the DOFMatrix.
Definition: MatrixBackend.hpp:29
void info(int noInfoLevel, std::string const &str, Args &&... args)
prints a message, if Environment::infoLevel() >= noInfoLevel
Definition: Output.hpp:105
PetscMatrix(Basis const &rowBasis, Basis const &)
Constructor. Constructs new BaseMatrix.
Definition: MatrixBackend.hpp:37
void init()
Reuse the matrix pattern and set all entries to zero.
Definition: MatrixBackend.hpp:151
BaseMatrix & matrix()
Return a reference to the data-matrix matrix.
Definition: MatrixBackend.hpp:53
std::size_t nnz() const
Return the local number of nonzeros in the matrix.
Definition: MatrixBackend.hpp:167
Definition: Constraints.hpp:13
void insert(RowIndex const &r, ColIndex const &c, PetscScalar value)
Insert a single value into the matrix.
Definition: MatrixBackend.hpp:66
static Dune::MPIHelper::MPICommunicator comm()
Return the MPI_Comm object (or a fake communicator)
Definition: Environment.hpp:86