5 #include <dune/common/ftraits.hh> 6 #include <dune/common/rangeutilities.hh> 7 #include <dune/common/shared_ptr.hh> 8 #include <dune/functions/functionspacebases/flatmultiindex.hh> 10 #include <amdis/Output.hpp> 11 #include <amdis/common/FakeContainer.hpp> 12 #include <amdis/operations/Assigner.hpp> 13 #include <amdis/typetree/MultiIndex.hpp> 18 template <
class DofMap>
22 using AssignMode = std::conditional_t<std::is_same_v<A,Assigner::plus_assign>,
23 std::integral_constant<::InsertMode, ADD_VALUES>,
24 std::integral_constant<::InsertMode, INSERT_VALUES>>;
38 template <
class Basis>
40 : dofMap_(basis.communicator().dofMap())
45 : dofMap_(other.dofMap_)
52 : dofMap_(other.dofMap_)
53 , initialized_(other.initialized_)
55 if (other.initialized_) {
56 VecDuplicate(other.vector_, &vector_);
57 VecCopy(other.vector_, vector_);
95 return dofMap_.localSize();
101 return dofMap_.globalSize();
106 template <
class SizeInfo>
107 void init(SizeInfo
const&,
bool clear)
113 info(3,
" create new vector needed {} seconds", t.elapsed());
119 VecGetLocalSize(vector_, &localSize);
122 if (std::size_t(localSize) != dofMap_.localSize() ||
123 std::size_t(
globalSize) != dofMap_.globalSize()) {
126 info(3,
" destroy old and create new vector needed {} seconds", t.elapsed());
142 VecAssemblyBegin(vector_);
143 VecAssemblyEnd(vector_);
144 info(3,
" finish vector assembly needed {} seconds", t.elapsed());
152 VecGhostUpdateBegin(vector_, INSERT_VALUES, SCATTER_FORWARD);
153 VecGhostUpdateEnd(vector_, INSERT_VALUES, SCATTER_FORWARD);
154 info(3,
" synchronize vector needed {} seconds", t.elapsed());
159 template <
class MultiIndex>
163 gather(std::array<MultiIndex,1>{idx}, &y);
169 template <
class MultiIndex,
class Assign>
172 assert(initialized_);
173 VecSetValue(vector_, dofMap_.global(flatMultiIndex(dof)), value, AssignMode<Assign>::value);
178 template <
class IndexRange,
class OutputIterator>
179 void gather(IndexRange
const& localInd, OutputIterator buffer)
const 181 forEach(localInd, [&buffer](
auto&&,
auto const& value) { *buffer = value; ++buffer; });
190 template <
class IndexRange,
class LocalValues,
class MaskRange,
class Assign>
191 void scatter(IndexRange
const& localInd, LocalValues
const& localVal, MaskRange
const& mask, Assign)
193 static_assert(std::is_same_v<typename LocalValues::value_type, PetscScalar>,
"");
194 assert(initialized_);
195 assert(localInd.size() == localVal.size());
198 std::vector<PetscInt> globalInd;
199 globalInd.reserve(localInd.size());
201 auto ind_it = std::begin(localInd);
202 auto mask_it = std::begin(mask);
203 for (; ind_it != std::end(localInd); ++mask_it, ++ind_it)
204 globalInd.push_back(*mask_it ? dofMap_.global(flatMultiIndex(*ind_it)) : -1);
207 VecSetOption(vector_, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
209 VecSetValues(vector_, globalInd.size(), globalInd.data(), localVal.data(), AssignMode<Assign>::value);
210 VecSetOption(vector_, VEC_IGNORE_NEGATIVE_INDICES, PETSC_FALSE);
215 template <
class IndexRange,
class Func>
216 void forEach(IndexRange
const& localInd, Func&& f)
const 219 Vec localGhostedVec =
nullptr;
220 VecGhostGetLocalForm(vector_, &localGhostedVec);
221 assert(localGhostedVec !=
nullptr);
225 VecGetArray(localGhostedVec, &ptr);
227 for (
auto const& dof : localInd)
230 if (dofMap_.owner(i)) {
232 f(dof, ptr[dofMap_.dofToLocal(i)]);
236 f(dof, ptr[dofMap_.localSize() + dofMap_.dofToGhost(i)]);
240 VecRestoreArray(localGhostedVec, &ptr);
241 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
244 template <
class Func>
248 [](
auto i) {
return Dune::Functions::FlatMultiIndex<decltype(i)>{i}; }), FWD(f));
251 template <
class Func>
255 [](
auto i) {
return Dune::Functions::FlatMultiIndex<decltype(i)>{i}; }), FWD(f));
261 void set(PetscScalar value)
263 Vec localGhostedVec =
nullptr;
264 VecGhostGetLocalForm(vector_, &localGhostedVec);
265 assert(localGhostedVec !=
nullptr);
267 VecSet(localGhostedVec, value);
268 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
275 Vec localGhostedVec =
nullptr;
276 VecGhostGetLocalForm(vector_, &localGhostedVec);
277 assert(localGhostedVec !=
nullptr);
279 VecZeroEntries(localGhostedVec);
280 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
285 assert(&lhs.dofMap_ == &rhs.dofMap_);
286 std::swap(lhs.vector_, rhs.vector_);
287 std::swap(lhs.initialized_, rhs.initialized_);
292 MPI_Comm comm()
const 297 return PETSC_COMM_SELF;
304 VecCreateGhost(comm(),
305 dofMap_.localSize(), dofMap_.globalSize(), dofMap_.ghostSize(), dofMap_.ghostIndices().data(),
313 VecDestroy(&vector_);
315 initialized_ =
false;
320 DofMap
const& dofMap_;
323 mutable Vec vector_ =
nullptr;
324 bool initialized_ =
false;
constexpr bool MultiIndex
A multi-index type.
Definition: Concepts.hpp:150
BaseVector const & vector() const
Return the data-vector vector_.
Definition: VectorBackend.hpp:81
PetscVector(PetscVector &&other)
Move constructor.
Definition: VectorBackend.hpp:44
BaseVector & vector()
Return the data-vector vector_.
Definition: VectorBackend.hpp:87
PetscVector(PetscVector const &other)
Copy constructor.
Definition: VectorBackend.hpp:51
~PetscVector()
Destructor.
Definition: VectorBackend.hpp:62
std::size_t globalSize() const
Return the number of entries in the global vector.
Definition: VectorBackend.hpp:99
void forEach(IndexRange const &localInd, Func &&f) const
Apply the functor f to each vector entry to local index in localInd range.
Definition: VectorBackend.hpp:216
PetscScalar value_type
The type of the elements of the DOFVector.
Definition: VectorBackend.hpp:31
void insert(MultiIndex const &dof, PetscScalar value, Assign)
Access the entry i of the vector with write-access.
Definition: VectorBackend.hpp:170
Contains all classes needed for solving linear and non linear equation systems.
Definition: AdaptBase.hpp:6
PetscVector & operator=(PetscVector &&other)
Move assignment operator.
Definition: VectorBackend.hpp:68
auto mappedRangeView(R &&range, F const &f)
Create a MappedRangeView.
Definition: MappedRangeView.hpp:445
PetscVector & operator=(PetscVector const &other)
Copy assignment operator.
Definition: VectorBackend.hpp:75
void init(SizeInfo const &, bool clear)
Resize the vector_ to the size s.
Definition: VectorBackend.hpp:107
void setZero()
Zero all entries, including the ghost entries.
Definition: VectorBackend.hpp:273
void scatter(IndexRange const &localInd, LocalValues const &localVal, MaskRange const &mask, Assign)
Store all values given by localVal with indices localInd in the vector.
Definition: VectorBackend.hpp:191
std::size_t localSize() const
Return the number of entries in the local part of the vector.
Definition: VectorBackend.hpp:93
void gather(IndexRange const &localInd, OutputIterator buffer) const
Collect all values from this vector to indices given by localInd and store it into the output buffer...
Definition: VectorBackend.hpp:179
A container-like data-structure not storing anything and with empty implementations in many container...
Definition: FakeContainer.hpp:34
void info(int noInfoLevel, std::string const &str, Args &&... args)
prints a message, if Environment::infoLevel() >= noInfoLevel
Definition: Output.hpp:105
PetscVector(Basis const &basis)
Constructor. Constructs new BaseVector.
Definition: VectorBackend.hpp:39
PetscInt size_type
The index/size - type.
Definition: VectorBackend.hpp:34
void synchronize()
Update the ghost regions with values from the owning process.
Definition: VectorBackend.hpp:149
::Vec BaseVector
The type of the base vector.
Definition: VectorBackend.hpp:28
void finish()
Finish assembly. Must be called before vector can be used in a KSP.
Definition: VectorBackend.hpp:139
PetscScalar at(MultiIndex const &idx) const
Access the entry i of the vector with read-only-access.
Definition: VectorBackend.hpp:160
The basic container that stores a base vector data.
Definition: VectorBackend.hpp:19
static Dune::MPIHelper::MPICommunicator comm()
Return the MPI_Comm object (or a fake communicator)
Definition: Environment.hpp:86