5 #include <dune/common/ftraits.hh> 6 #include <dune/common/rangeutilities.hh> 7 #include <dune/common/shared_ptr.hh> 8 #include <dune/functions/common/multiindex.hh> 10 #include <amdis/Output.hpp> 11 #include <amdis/algorithm/ForEach.hpp> 12 #include <amdis/algorithm/InnerProduct.hpp> 13 #include <amdis/algorithm/Transform.hpp> 14 #include <amdis/common/FakeContainer.hpp> 15 #include <amdis/operations/Assigner.hpp> 16 #include <amdis/typetree/MultiIndex.hpp> 21 template <
class DofMap>
25 using AssignMode = std::conditional_t<std::is_same_v<A,Assigner::plus_assign>,
26 std::integral_constant<::InsertMode, ADD_VALUES>,
27 std::integral_constant<::InsertMode, INSERT_VALUES>>;
47 : dofMap_(other.dofMap_)
54 : dofMap_(other.dofMap_)
55 , initialized_(other.initialized_)
57 if (other.initialized_) {
58 VecDuplicate(other.vector_, &vector_);
59 VecCopy(other.vector_, vector_);
97 return dofMap_.localSize();
103 return dofMap_.globalSize();
108 template <
class Basis>
109 void init(Basis
const&,
bool clear)
117 info(3,
" create new vector needed {} seconds", t.elapsed());
131 VecAssemblyBegin(vector_);
132 VecAssemblyEnd(vector_);
133 info(3,
" finish vector assembly needed {} seconds", t.elapsed());
141 VecGhostUpdateBegin(vector_, INSERT_VALUES, SCATTER_FORWARD);
142 VecGhostUpdateEnd(vector_, INSERT_VALUES, SCATTER_FORWARD);
143 info(3,
" synchronize vector needed {} seconds", t.elapsed());
148 template <
class MultiIndex>
152 gather(std::array<MultiIndex,1>{idx}, &y);
158 template <
class MultiIndex,
class Assign>
161 assert(initialized_);
162 VecSetValue(vector_, dofMap_.global(flatMultiIndex(dof)), value, AssignMode<Assign>::value);
167 template <
class IndexRange,
class OutputIterator>
168 void gather(IndexRange
const& localInd, OutputIterator buffer)
const 170 forEach(localInd, [&buffer](
auto&&,
auto const& value) { *buffer = value; ++buffer; });
179 template <
class IndexRange,
class LocalValues,
class MaskRange,
class Assign>
180 void scatter(IndexRange
const& localInd, LocalValues
const& localVal, MaskRange
const& mask, Assign)
182 if constexpr(not std::is_same_v<typename LocalValues::value_type, PetscScalar>)
184 std::vector<PetscScalar> newLocalVal(localVal.begin(), localVal.end());
185 scatter(localInd, newLocalVal, mask, Assign{});
189 assert(initialized_);
190 assert(localInd.size() == localVal.size());
193 std::vector<PetscInt> globalInd;
194 globalInd.reserve(localInd.size());
196 auto ind_it = std::begin(localInd);
197 auto mask_it = std::begin(mask);
198 for (; ind_it != std::end(localInd); ++mask_it, ++ind_it)
199 globalInd.push_back(*mask_it ? dofMap_.global(flatMultiIndex(*ind_it)) : -1);
202 VecSetOption(vector_, VEC_IGNORE_NEGATIVE_INDICES, PETSC_TRUE);
204 VecSetValues(vector_, globalInd.size(), globalInd.data(), localVal.data(), AssignMode<Assign>::value);
205 VecSetOption(vector_, VEC_IGNORE_NEGATIVE_INDICES, PETSC_FALSE);
210 template <
class IndexRange,
class Func>
211 void forEach(IndexRange
const& localInd, Func&& f)
const 214 Vec localGhostedVec =
nullptr;
215 VecGhostGetLocalForm(vector_, &localGhostedVec);
216 assert(localGhostedVec !=
nullptr);
220 VecGetArray(localGhostedVec, &ptr);
222 for (
auto const& dof : localInd)
225 if (dofMap_.owner(i)) {
227 f(dof, ptr[dofMap_.dofToLocal(i)]);
231 f(dof, ptr[dofMap_.localSize() + dofMap_.dofToGhost(i)]);
235 VecRestoreArray(localGhostedVec, &ptr);
236 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
241 template <
class Func>
245 PetscScalar
const* ptr;
246 VecGetArrayRead(vector_, &ptr);
249 f(ptr[dofMap_.dofToLocal(i)]);
251 VecRestoreArrayRead(vector_, &ptr);
257 template <
class Func>
262 VecGetArray(vector_, &ptr);
265 f(ptr[dofMap_.dofToLocal(i)]);
267 VecRestoreArray(vector_, &ptr);
272 template <
class Operation,
class... VecIns>
276 std::array<PetscScalar
const*,
sizeof...(VecIns)> ptrsIn{};
279 VecGetArray(vecOut.vector_, &ptrOut);
280 Ranges::forIndices<
sizeof...(VecIns)>([&](
auto ii) {
281 VecGetArrayRead(std::get<ii>(std::tie(vecIn...)).vector_, &ptrsIn[ii]);
285 auto idx = vecOut.dofMap_.dofToLocal(i);
286 Ranges::applyIndices<
sizeof...(VecIns)>([&](
auto... ii) {
287 ptrOut[idx] = op(ptrsIn[ii][idx]...);
292 VecRestoreArray(vecOut.vector_, &ptrOut);
293 Ranges::forIndices<
sizeof...(VecIns)>([&](
auto ii) {
294 VecRestoreArrayRead(std::get<ii>(std::tie(vecIn...)).vector_, &ptrsIn[ii]);
300 template <
class T,
class BinOp1,
class BinOp2>
303 PetscScalar
const* ptr1;
304 PetscScalar
const* ptr2;
307 VecGetArrayRead(in1.vector_, &ptr1);
308 VecGetArrayRead(in2.vector_, &ptr2);
311 auto idx = in1.dofMap_.dofToLocal(i);
312 init = op1(std::move(init), op2(ptr1[idx], ptr2[idx]));
316 VecRestoreArrayRead(in1.vector_, &ptr1);
317 VecRestoreArrayRead(in2.vector_, &ptr2);
324 void set(PetscScalar value)
326 Vec localGhostedVec =
nullptr;
327 VecGhostGetLocalForm(vector_, &localGhostedVec);
328 assert(localGhostedVec !=
nullptr);
330 VecSet(localGhostedVec, value);
331 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
338 Vec localGhostedVec =
nullptr;
339 VecGhostGetLocalForm(vector_, &localGhostedVec);
340 assert(localGhostedVec !=
nullptr);
342 VecZeroEntries(localGhostedVec);
343 VecGhostRestoreLocalForm(vector_, &localGhostedVec);
348 assert(&lhs.dofMap_ == &rhs.dofMap_);
349 std::swap(lhs.vector_, rhs.vector_);
350 std::swap(lhs.initialized_, rhs.initialized_);
354 MPI_Comm comm()
const 356 return dofMap_.comm();
363 VecCreateGhost(comm(),
364 dofMap_.localSize(), dofMap_.globalSize(), dofMap_.ghostSize(), dofMap_.ghostIndices().data(),
372 VecDestroy(&vector_);
374 initialized_ =
false;
379 DofMap
const& dofMap_;
382 mutable Vec vector_ =
nullptr;
383 bool initialized_ =
false;
389 template <
class DofMap>
392 template <
class Vec,
class UnaryFunction>
393 static void impl(Vec&& vec, UnaryFunction f)
399 template <
class DofMap>
402 template <
class Operation,
class... VecIns>
409 template <
class DofMap>
412 template <
class T,
class BinOp1,
class BinOp2>
void forEachLocal(Func &&f) const
Perform the f(value) on the local elements of this vector.
Definition: VectorBackend.hpp:242
constexpr bool MultiIndex
A multi-index type.
Definition: Concepts.hpp:149
BaseVector const & vector() const
Return the data-vector vector_.
Definition: VectorBackend.hpp:83
PetscVector(PetscVector &&other)
Move constructor.
Definition: VectorBackend.hpp:46
BaseVector & vector()
Return the data-vector vector_.
Definition: VectorBackend.hpp:89
PetscVector(PetscVector const &other)
Copy constructor.
Definition: VectorBackend.hpp:53
~PetscVector()
Destructor.
Definition: VectorBackend.hpp:64
void forEachLocal(Func &&f)
Perform the f(value) on the local elements of this vector.
Definition: VectorBackend.hpp:258
std::size_t globalSize() const
Return the number of entries in the global vector.
Definition: VectorBackend.hpp:101
Definition: ForEach.hpp:17
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:211
PetscScalar value_type
The type of the elements of the DOFVector.
Definition: VectorBackend.hpp:34
void insert(MultiIndex const &dof, PetscScalar value, Assign)
Access the entry i of the vector with write-access.
Definition: VectorBackend.hpp:159
Definition: AdaptBase.hpp:6
PetscVector & operator=(PetscVector &&other)
Move assignment operator.
Definition: VectorBackend.hpp:70
static void transformLocal(PetscVector &vecOut, Operation op, VecIns const &... vecIn)
Perform the valueOut = op(valueIn...) on the local elements of [vecOut] and [vecIn...].
Definition: VectorBackend.hpp:273
PetscVector & operator=(PetscVector const &other)
Copy assignment operator.
Definition: VectorBackend.hpp:77
General implementation of recursive inner-product.
Definition: InnerProduct.hpp:16
void setZero()
Zero all entries, including the ghost entries.
Definition: VectorBackend.hpp:336
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:180
std::size_t localSize() const
Return the number of entries in the local part of the vector.
Definition: VectorBackend.hpp:95
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:168
A container-like data-structure not storing anything and with empty implementations in many container...
Definition: FakeContainer.hpp:34
PetscInt size_type
The index/size - type.
Definition: VectorBackend.hpp:37
static T innerProductLocal(PetscVector const &in1, PetscVector const &in2, T init, BinOp1 op1, BinOp2 op2)
Perform the op1(init, op2(value1, value2)) on the local elements of [in1] and [in2].
Definition: VectorBackend.hpp:301
void synchronize()
Update the ghost regions with values from the owning process.
Definition: VectorBackend.hpp:138
::Vec BaseVector
The type of the base vector.
Definition: VectorBackend.hpp:31
void finish()
Finish assembly. Must be called before vector can be used in a KSP.
Definition: VectorBackend.hpp:128
void init(Basis const &, bool clear)
Resize the vector_ to the size s.
Definition: VectorBackend.hpp:109
PetscScalar at(MultiIndex const &idx) const
Access the entry i of the vector with read-only-access.
Definition: VectorBackend.hpp:149
PetscVector(DofMap const &dofMap)
Constructor. Constructs new BaseVector.
Definition: VectorBackend.hpp:41
The basic container that stores a base vector data.
Definition: VectorBackend.hpp:22