1#ifndef DUNE_PDELAB_BACKEND_ISTL_GENEO_MULTICOMMDATAHANDLE_HH
2#define DUNE_PDELAB_BACKEND_ISTL_GENEO_MULTICOMMDATAHANDLE_HH
13 template<
typename GFS,
typename RankIndex,
typename V>
18 typedef typename V::template LocalView<IndexCache> LocalView;
24 template<
typename MessageBuffer,
typename Entity,
typename LocalView>
25 bool gather(MessageBuffer& buff,
const Entity& e, LocalView& local_view)
const
28 for (std::size_t i = 0; i < local_view.size(); ++i) {
29 buff.write (local_view[i]);
34 template<
typename MessageBuffer,
typename Entity,
typename LocalView>
35 bool scatter(MessageBuffer& buff, std::size_t n,
const Entity& e, LocalView& local_view)
const
37 RankIndex remote_rank = buff.senderRank();
40 int remote_index = std::distance(_neighbor_ranks.begin(), std::find(_neighbor_ranks.begin(), _neighbor_ranks.end(), remote_rank));
41 if (remote_index ==
static_cast<RankIndex
>(_neighbor_ranks.size()))
42 DUNE_THROW(
Exception,
"Received remote rank " << remote_rank <<
", but it's not in the given neighbor set!");
45 auto target_view = LocalView(*_target_vectors[remote_index]);
47 target_view.bind(_index_cache);
50 if (target_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType())) {
51 if (target_view.size() != n)
52 DUNE_THROW(
Exception,
"size mismatch in GridFunctionSpace data handle, have " << target_view.size() <<
"DOFs, but received " << n);
54 for (
size_type i = 0; i < target_view.size(); ++i) {
55 typename LocalView::ElementType x;
64 if (target_view.size() != 0)
65 DUNE_THROW(
Exception,
"expected no DOFs in partition '" << e.partitionType() <<
"', but have " << target_view.size());
67 for (
size_type i = 0; i < target_view.size(); ++i) {
68 typename LocalView::ElementType dummy;
77 template<
typename MessageBuffer,
typename Offsets,
typename Entity,
typename LocalView>
78 bool scatter(MessageBuffer& buff,
const Offsets& remote_offsets,
const Offsets& local_offsets,
const Entity& e, LocalView& local_view)
const
80 if (local_view.cache().gridFunctionSpace().entitySet().partitions().contains(e.partitionType()))
90 RankIndex remote_rank = buff.senderRank();
93 int remote_index = std::distance(_neighbor_ranks.begin(), std::find(_neighbor_ranks.begin(), _neighbor_ranks.end(), remote_rank));
94 if (remote_index ==
static_cast<RankIndex
>(_neighbor_ranks.size()))
95 DUNE_THROW(
Exception,
"Received remote rank " << remote_rank <<
", but it's not in the given neighbor set!");
98 auto target_view = LocalView(*_target_vectors[remote_index]);
100 target_view.bind(_index_cache);
105 bool needs_commit =
false;
106 for (
size_type block = 1; block < local_offsets.size(); ++block)
109 if (remote_offsets[block] == remote_i)
111 local_i = local_offsets[block];
116 if (local_offsets[block] == local_i)
118 for (; remote_i < remote_offsets[block]; ++remote_i)
120 typename LocalView::ElementType dummy;
126 if (remote_offsets[block] - remote_i != local_offsets[block] - local_i)
127 DUNE_THROW(
Exception,
"size mismatch in GridFunctionSpace data handle block " << block <<
", have " << local_offsets[block] - local_i <<
"DOFs, but received " << remote_offsets[block] - remote_i);
129 for (; local_i < local_offsets[block]; ++local_i) {
130 typename LocalView::ElementType x;
132 target_view[local_i] = x;
135 remote_i = remote_offsets[block];
138 target_view.unbind();
143 if (local_view.size() != 0)
144 DUNE_THROW(
Exception,
"expected no DOFs in partition '" << e.partitionType() <<
"', but have " << local_view.size());
146 for (std::size_t i = 0; i < remote_offsets.back(); ++i)
148 typename LocalView::ElementType dummy;
162 std::vector<std::shared_ptr<V> > target_vectors,
163 std::vector<RankIndex> neighbor_ranks)
166 , _target_vectors(target_vectors)
167 , _neighbor_ranks(neighbor_ranks)
172 mutable IndexCache _index_cache;
174 const RankIndex _rank;
175 std::vector<std::shared_ptr<V> > _target_vectors;
176 std::vector<RankIndex> _neighbor_ranks;
184 template<
class GFS,
class V,
typename RankIndex>
187 MultiCommGatherScatter<GFS, RankIndex, V>,
188 Dune::PDELab::DOFDataCommunicationDescriptor<typename V::ElementType,true>>
202 MultiCommDataHandle(
const GFS& gfs_, V& v_, std::vector<std::shared_ptr<V> > target_vectors, std::vector<RankIndex> neighbor_ranks)
For backward compatibility – Do not use this!
Definition: adaptivity.hh:28
Gather/scatter communication that passes a single function from each subdomain to all its neighbors.
Definition: multicommdatahandle.hh:15
MultiCommGatherScatter(const GFS &gfs, RankIndex rank, std::vector< std::shared_ptr< V > > target_vectors, std::vector< RankIndex > neighbor_ranks)
Definition: multicommdatahandle.hh:161
bool scatter(MessageBuffer &buff, const Offsets &remote_offsets, const Offsets &local_offsets, const Entity &e, LocalView &local_view) const
Definition: multicommdatahandle.hh:78
bool scatter(MessageBuffer &buff, std::size_t n, const Entity &e, LocalView &local_view) const
Definition: multicommdatahandle.hh:35
bool gather(MessageBuffer &buff, const Entity &e, LocalView &local_view) const
Definition: multicommdatahandle.hh:25
std::size_t size_type
Definition: multicommdatahandle.hh:22
Gather/scatter communication that passes a single function from each subdomain to all its neighbors....
Definition: multicommdatahandle.hh:189
MultiCommDataHandle(const GFS &gfs_, V &v_, std::vector< std::shared_ptr< V > > target_vectors, std::vector< RankIndex > neighbor_ranks)
Definition: multicommdatahandle.hh:202
Base class for all PDELab exceptions.
Definition: exceptions.hh:19
Definition: entityindexcache.hh:18
void update(const Entity &e)
Definition: entityindexcache.hh:50
Communication descriptor for sending one item of type E per DOF.
Definition: genericdatahandle.hh:25
Implement a data handle with a grid function space.
Definition: genericdatahandle.hh:132