1 // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
2 // vi: set et ts=4 sw=2 sts=2:
3
4 #include "config.h"
5
6 #include <cassert>
7 #include <cstdlib>
8
9 #include <dune/grid/yaspgrid.hh>
10
11 #include <dune/pdelab.hh>
12
13
14 /*
15 * This test makes sure that the GenericDataHandle infrastructure correctly transmits
16 * the sender rank if requested, both for per-DOF and per-entity communications, and
17 * that the reported communication size is still correct.
18 */
19
20 // This custom gather/scatter handler sends the rank and checks its value when receiving
21 struct GatherScatter
22 {
23
24 template<typename Buffer, typename Entity, typename LocalView>
gatherGatherScatter25 bool gather(Buffer& buffer, const Entity& e, const LocalView& local_view) const
26 {
27 std::size_t size = _expected_size >= 0 ? _expected_size : local_view.size();
28 for (std::size_t i = 0 ; i < size ; ++i)
29 buffer.write(_rank);
30 return false;
31 }
32
33 template<typename Buffer, typename Entity, typename LocalView>
scatterGatherScatter34 bool scatter(Buffer& buffer, std::size_t n, const Entity& e, const LocalView& local_view) const
35 {
36 std::size_t size = 0;
37 if (_expected_size >= 0)
38 size = _expected_size;
39 else
40 size = local_view.size();
41
42 assert(size == n);
43
44 for (std::size_t i = 0 ; i < local_view.size() ; ++i)
45 {
46 int rank;
47 buffer.read(rank);
48 assert(rank == buffer.senderRank());
49 }
50 return false;
51 }
52
53 template<typename Buffer, typename O1, typename O2, typename Entity, typename LocalView>
scatterGatherScatter54 bool scatter(Buffer& buffer, const O1& o1, const O2& o2, const Entity& e, const LocalView& local_view) const
55 {
56 assert(false && "Should never get here - this only exists to make the test compile");
57 return false;
58 }
59
60 // If expected_size is < 0, take the size from the LocalView, otherwise use the prescribed value
GatherScatterGatherScatter61 GatherScatter(int rank, int expected_size)
62 : _rank(rank)
63 , _expected_size(expected_size)
64 {}
65
66 private:
67
68 int _rank;
69 int _expected_size;
70
71 };
72
73
74 template<typename GFS, typename V>
75 class DOFDataHandle
76 : public Dune::PDELab::GFSDataHandle<
77 GFS,V,GatherScatter,Dune::PDELab::DOFDataCommunicationDescriptor<typename V::ElementType,true>>
78 {
79
80 using Base = Dune::PDELab::GFSDataHandle<
81 GFS,
82 V,
83 GatherScatter,
84 Dune::PDELab::DOFDataCommunicationDescriptor<typename V::ElementType,true>
85 >;
86
87 public:
88
DOFDataHandle(const GFS & gfs,V & v)89 DOFDataHandle(const GFS& gfs, V& v)
90 : Base(gfs,v,GatherScatter(gfs.gridView().comm().rank(),-1))
91 {}
92
93 };
94
95
96 template<typename GFS, typename V>
97 class EntityDataHandle
98 : public Dune::PDELab::GFSDataHandle<
99 GFS,V,GatherScatter,Dune::PDELab::EntityDataCommunicationDescriptor<typename V::ElementType,true>>
100 {
101
102 using Base = Dune::PDELab::GFSDataHandle<
103 GFS,
104 V,
105 GatherScatter,
106 Dune::PDELab::DOFDataCommunicationDescriptor<typename V::ElementType,true>
107 >;
108
109 public:
110
EntityDataHandle(const GFS & gfs,V & v,std::size_t count)111 EntityDataHandle(const GFS& gfs, V& v, std::size_t count)
112 : Base(gfs,v,GatherScatter(gfs.gridView().comm().rank(),count))
113 {}
114
115 };
116
117
118
main(int argc,char ** argv)119 int main(int argc, char** argv)
120 {
121
122 try {
123 auto& helper = Dune::MPIHelper::instance(argc,argv);
124
125 if (helper.getCollectiveCommunication().size() == 1)
126 {
127 std::cerr << "Communication test requires at least 2 MPI ranks" << std::endl;
128 std::exit(77);
129 }
130
131 Dune::YaspGrid<2> grid({{1.0,1.0}},{{16,16}},0,1);
132 auto gv = grid.leafGridView();
133 using GV = decltype(gv);
134 using F = GV::ctype;
135
136 using FEM = Dune::PDELab::QkLocalFiniteElementMap<GV,F,F,2>;
137 FEM fem(gv);
138
139 using GFS = Dune::PDELab::GridFunctionSpace<GV,FEM>;
140 GFS gfs(gv,fem);
141
142 using V = Dune::PDELab::Backend::Vector<GFS,int>;
143 auto v = V(gfs);
144
145 {
146 // start with a standard DOF communication - no rank transmission
147 auto handle = Dune::PDELab::CopyDataHandle<GFS,V>(gfs,v);
148 gv.communicate(handle,Dune::All_All_Interface,Dune::ForwardCommunication);
149 }
150
151 {
152 // DOF communication with rank transmission
153 auto handle = DOFDataHandle<GFS,V>(gfs,v);
154 gv.communicate(handle,Dune::All_All_Interface,Dune::ForwardCommunication);
155 }
156
157 {
158 // standard per-entity communication - no rank transmission
159 auto handle = Dune::PDELab::CopyDataHandle<GFS,V>(gfs,v);
160 gv.communicate(handle,Dune::All_All_Interface,Dune::ForwardCommunication);
161 }
162
163 {
164 // DOF communication with rank transmission
165 auto handle = DOFDataHandle<GFS,V>(gfs,v);
166 gv.communicate(handle,Dune::All_All_Interface,Dune::ForwardCommunication);
167 }
168
169 return 0;
170 } catch (std::exception& e) {
171 std::cerr << "error: caught exception: " << e.what() << std::endl;
172 std::abort();
173 } catch (...) {
174 std::cerr << "error: caught unknown exception" << std::endl;
175 std::abort();
176 }
177 }
178