Line data Source code
1 : // Include header
2 : #include "IndexingService.h"
3 :
4 : #include <algorithm>
5 : #include <iterator>
6 : #include <numeric>
7 : #include <sstream>
8 :
9 : /*
10 : * Define specialisations of CellIndex.
11 : * Differentiate between dim=2 and dim=3.
12 : *
13 : * Only declarations and default inits should be made here, lower/upper
14 : * boundaries must be determined at runtime using IndexingService
15 : */
16 :
17 : // #DEFINE INDEXING_ENABLE_DIM2
18 : #ifdef INDEXING_ENABLE_DIM2
19 : // We must compile both, dim2 and dim3, so that both can be used (also in the same executable e.g. main_lammps.cpp)
20 : // Dim2 ///////////
21 : namespace coupling {
22 : namespace indexing {
23 :
24 : /*
25 : * Declare specialisations of CellIndex.
26 : * Define their static members.
27 : */
28 :
29 : /*
30 : * NON-MD-TO-MACRO aka MAMICO INDEXING, INCL GHOST LAYER
31 : */
32 :
33 : // scalar, global, !md2macro, !noGL
34 : template <> BaseIndex<2> CellIndex<2>::lowerBoundary{};
35 : template <> BaseIndex<2> CellIndex<2>::upperBoundary{};
36 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2>::numberCellsInDomain{};
37 : template <> unsigned int CellIndex<2>::linearNumberCellsInDomain{};
38 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2>::divisionFactor{};
39 :
40 : // BaseIndex
41 : template <> BaseIndex<2> BaseIndex<2>::lowerBoundary{};
42 : template <> BaseIndex<2> BaseIndex<2>::upperBoundary{};
43 : template <> tarch::la::Vector<2, unsigned int> BaseIndex<2>::numberCellsInDomain{};
44 : template <> unsigned int BaseIndex<2>::linearNumberCellsInDomain{};
45 : template <> tarch::la::Vector<2, unsigned int> BaseIndex<2>::divisionFactor{};
46 :
47 : // scalar, local, !md2macro, !noGL
48 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local>::lowerBoundary{};
49 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local>::upperBoundary{};
50 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local>::numberCellsInDomain{};
51 : template <> unsigned int CellIndex<2, IndexTrait::local>::linearNumberCellsInDomain{};
52 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local>::divisionFactor{};
53 :
54 : // vector, local, !md2macro, !noGL
55 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local>::lowerBoundary{};
56 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local>::upperBoundary{};
57 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local>::numberCellsInDomain{};
58 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::local>::linearNumberCellsInDomain{};
59 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local>::divisionFactor{};
60 :
61 : /*
62 : * MD TO MACRO, INCL GHOST LAYER
63 : */
64 :
65 : // scalar, global, md2macro, !noGL
66 : template <> BaseIndex<2> CellIndex<2, IndexTrait::md2macro>::lowerBoundary{};
67 : template <> BaseIndex<2> CellIndex<2, IndexTrait::md2macro>::upperBoundary{};
68 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::md2macro>::numberCellsInDomain{};
69 : template <> unsigned int CellIndex<2, IndexTrait::md2macro>::linearNumberCellsInDomain{};
70 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::md2macro>::divisionFactor{};
71 :
72 : // vector, global, md2macro, !noGL
73 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro>::lowerBoundary{};
74 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro>::upperBoundary{};
75 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro>::numberCellsInDomain{};
76 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::md2macro>::linearNumberCellsInDomain{};
77 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro>::divisionFactor{};
78 :
79 : // scalar, local, md2macro, !noGL
80 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::md2macro>::lowerBoundary{};
81 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::md2macro>::upperBoundary{};
82 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::md2macro>::numberCellsInDomain{};
83 : template <> unsigned int CellIndex<2, IndexTrait::local, IndexTrait::md2macro>::linearNumberCellsInDomain{};
84 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::md2macro>::divisionFactor{};
85 :
86 : // vector, local, md2macro, !noGL
87 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::lowerBoundary{};
88 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::upperBoundary{};
89 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::numberCellsInDomain{};
90 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::linearNumberCellsInDomain{};
91 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::divisionFactor{};
92 :
93 : /*
94 : * !MD TO MACRO aka MAMICO INDEXING, EXCL GHOST LAYER
95 : */
96 :
97 : // scalar, global, !md2macro, noGL
98 : template <> BaseIndex<2> CellIndex<2, IndexTrait::noGhost>::lowerBoundary{};
99 : template <> BaseIndex<2> CellIndex<2, IndexTrait::noGhost>::upperBoundary{};
100 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::noGhost>::numberCellsInDomain{};
101 : template <> unsigned int CellIndex<2, IndexTrait::noGhost>::linearNumberCellsInDomain{};
102 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::noGhost>::divisionFactor{};
103 :
104 : // vector, global, !md2macro, noGL
105 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::noGhost>::lowerBoundary{};
106 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::noGhost>::upperBoundary{};
107 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::noGhost>::numberCellsInDomain{};
108 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::noGhost>::linearNumberCellsInDomain{};
109 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::noGhost>::divisionFactor{};
110 :
111 : // scalar, local, !md2macro, noGL
112 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::noGhost>::lowerBoundary{};
113 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::noGhost>::upperBoundary{};
114 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::noGhost>::numberCellsInDomain{};
115 : template <> unsigned int CellIndex<2, IndexTrait::local, IndexTrait::noGhost>::linearNumberCellsInDomain{};
116 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::noGhost>::divisionFactor{};
117 :
118 : // vector, local, !md2macro, noGL
119 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::lowerBoundary{};
120 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::upperBoundary{};
121 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::numberCellsInDomain{};
122 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::linearNumberCellsInDomain{};
123 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::divisionFactor{};
124 :
125 : /*
126 : * MD TO MACRO, EXCL GHOST LAYER
127 : */
128 :
129 : // scalar, global, md2macro, noGL
130 : template <> BaseIndex<2> CellIndex<2, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary{};
131 : template <> BaseIndex<2> CellIndex<2, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary{};
132 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::md2macro, IndexTrait::noGhost>::numberCellsInDomain{};
133 : template <> unsigned int CellIndex<2, IndexTrait::md2macro, IndexTrait::noGhost>::linearNumberCellsInDomain{};
134 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::md2macro, IndexTrait::noGhost>::divisionFactor{};
135 :
136 : // vector, global, md2macro, noGL
137 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary{};
138 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary{};
139 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::numberCellsInDomain{};
140 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::linearNumberCellsInDomain{};
141 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::divisionFactor{};
142 :
143 : // scalar, local, md2macro, noGL
144 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary{};
145 : template <> BaseIndex<2> CellIndex<2, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary{};
146 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::numberCellsInDomain{};
147 : template <> unsigned int CellIndex<2, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::linearNumberCellsInDomain{};
148 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::divisionFactor{};
149 :
150 : // vector, local, md2macro, noGL
151 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary{};
152 : template <> BaseIndex<2> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary{};
153 : template <>
154 : tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::numberCellsInDomain{};
155 : template <> unsigned int CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::linearNumberCellsInDomain{};
156 : template <> tarch::la::Vector<2, unsigned int> CellIndex<2, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::divisionFactor{};
157 : } // namespace indexing
158 : } // namespace coupling
159 : #endif
160 :
161 : // Dim3 //////////////////////////
162 : namespace coupling {
163 : namespace indexing {
164 :
165 : /*
166 : * Declare specialisations of CellIndex.
167 : * Define their static members.
168 : */
169 :
170 : /*
171 : * NON-MD-TO-MACRO aka MAMICO INDEXING, INCL GHOST LAYER
172 : */
173 :
174 : // scalar, global, !md2macro, !noGL
175 : template <> BaseIndex<3> I00::lowerBoundary{};
176 : template <> BaseIndex<3> I00::upperBoundary{};
177 : template <> tarch::la::Vector<3, unsigned int> I00::numberCellsInDomain{};
178 : template <> unsigned int I00::linearNumberCellsInDomain{};
179 : template <> tarch::la::Vector<3, unsigned int> I00::divisionFactor{};
180 : template <> const char I00::TNAME[] = "CellIndex<3>";
181 :
182 : // BaseIndex
183 : template <> BaseIndex<3> BaseIndex<3>::lowerBoundary{};
184 : template <> BaseIndex<3> BaseIndex<3>::upperBoundary{};
185 : template <> tarch::la::Vector<3, unsigned int> BaseIndex<3>::numberCellsInDomain{};
186 : template <> unsigned int BaseIndex<3>::linearNumberCellsInDomain{};
187 : template <> tarch::la::Vector<3, unsigned int> BaseIndex<3>::divisionFactor{};
188 : template <> const char CellIndex<3, IndexTrait::vector>::TNAME[] = "CellIndex<3, vector>";
189 :
190 : // scalar, local, !md2macro, !noGL
191 : template <> BaseIndex<3> I02::lowerBoundary{};
192 : template <> BaseIndex<3> I02::upperBoundary{};
193 : template <> tarch::la::Vector<3, unsigned int> I02::numberCellsInDomain{};
194 : template <> unsigned int I02::linearNumberCellsInDomain{};
195 : template <> tarch::la::Vector<3, unsigned int> I02::divisionFactor{};
196 : template <> const char I02::TNAME[] = "CellIndex<3, local>";
197 :
198 : // vector, local, !md2macro, !noGL
199 : template <> BaseIndex<3> I03::lowerBoundary{};
200 : template <> BaseIndex<3> I03::upperBoundary{};
201 : template <> tarch::la::Vector<3, unsigned int> I03::numberCellsInDomain{};
202 : template <> unsigned int I03::linearNumberCellsInDomain{};
203 : template <> tarch::la::Vector<3, unsigned int> I03::divisionFactor{};
204 : template <> const char I03::TNAME[] = "CellIndex<3, vector, local>";
205 :
206 : /*
207 : * MD TO MACRO, INCL GHOST LAYER
208 : */
209 :
210 : // scalar, global, md2macro, !noGL
211 : template <> BaseIndex<3> I04::lowerBoundary{};
212 : template <> BaseIndex<3> I04::upperBoundary{};
213 : template <> tarch::la::Vector<3, unsigned int> I04::numberCellsInDomain{};
214 : template <> unsigned int I04::linearNumberCellsInDomain{};
215 : template <> tarch::la::Vector<3, unsigned int> I04::divisionFactor{};
216 : template <> const char I04::TNAME[] = "CellIndex<3, md2macro>";
217 :
218 : // vector, global, md2macro, !noGL
219 : template <> BaseIndex<3> I05::lowerBoundary{};
220 : template <> BaseIndex<3> I05::upperBoundary{};
221 : template <> tarch::la::Vector<3, unsigned int> I05::numberCellsInDomain{};
222 : template <> unsigned int I05::linearNumberCellsInDomain{};
223 : template <> tarch::la::Vector<3, unsigned int> I05::divisionFactor{};
224 : template <> const char I05::TNAME[] = "CellIndex<3, vector, md2macro>";
225 :
226 : // scalar, local, md2macro, !noGL
227 : template <> BaseIndex<3> I06::lowerBoundary{};
228 : template <> BaseIndex<3> I06::upperBoundary{};
229 : template <> tarch::la::Vector<3, unsigned int> I06::numberCellsInDomain{};
230 : template <> unsigned int I06 ::linearNumberCellsInDomain{};
231 : template <> tarch::la::Vector<3, unsigned int> I06::divisionFactor{};
232 : template <> const char I06::TNAME[] = "CellIndex<3, local, md2macro>";
233 :
234 : // vector, local, md2macro, !noGL
235 : template <> BaseIndex<3> I07::lowerBoundary{};
236 : template <> BaseIndex<3> I07::upperBoundary{};
237 : template <> tarch::la::Vector<3, unsigned int> I07::numberCellsInDomain{};
238 : template <> unsigned int I07::linearNumberCellsInDomain{};
239 : template <> tarch::la::Vector<3, unsigned int> I07::divisionFactor{};
240 : template <> const char I07::TNAME[] = "CellIndex<3, vector, local, md2macro>";
241 :
242 : /*
243 : * !MD TO MACRO aka MAMICO INDEXING, EXCL GHOST LAYER
244 : */
245 :
246 : // scalar, global, !md2macro, noGL
247 : template <> BaseIndex<3> I08::lowerBoundary{};
248 : template <> BaseIndex<3> I08::upperBoundary{};
249 : template <> tarch::la::Vector<3, unsigned int> I08::numberCellsInDomain{};
250 : template <> unsigned int I08::linearNumberCellsInDomain{};
251 : template <> tarch::la::Vector<3, unsigned int> I08::divisionFactor{};
252 : template <> const char I08::TNAME[] = "CellIndex<3, noGhost>";
253 :
254 : // vector, global, !md2macro, noGL
255 : template <> BaseIndex<3> I09::lowerBoundary{};
256 : template <> BaseIndex<3> I09::upperBoundary{};
257 : template <> tarch::la::Vector<3, unsigned int> I09::numberCellsInDomain{};
258 : template <> unsigned int I09::linearNumberCellsInDomain{};
259 : template <> tarch::la::Vector<3, unsigned int> I09::divisionFactor{};
260 : template <> const char I09::TNAME[] = "CellIndex<3, vector, noGhost>";
261 :
262 : // scalar, local, !md2macro, noGL
263 : template <> BaseIndex<3> I10::lowerBoundary{};
264 : template <> BaseIndex<3> I10::upperBoundary{};
265 : template <> tarch::la::Vector<3, unsigned int> I10::numberCellsInDomain{};
266 : template <> unsigned int I10::linearNumberCellsInDomain{};
267 : template <> tarch::la::Vector<3, unsigned int> I10::divisionFactor{};
268 : template <> const char I10::TNAME[] = "CellIndex<3, local, noGhost>";
269 :
270 : // vector, local, !md2macro, noGL
271 : template <> BaseIndex<3> I11::lowerBoundary{};
272 : template <> BaseIndex<3> I11::upperBoundary{};
273 : template <> tarch::la::Vector<3, unsigned int> I11::numberCellsInDomain{};
274 : template <> unsigned int I11::linearNumberCellsInDomain{};
275 : template <> tarch::la::Vector<3, unsigned int> I11::divisionFactor{};
276 : template <> const char I11::TNAME[] = "CellIndex<3, vector, local, noGhost>";
277 :
278 : /*
279 : * MD TO MACRO, EXCL GHOST LAYER
280 : */
281 :
282 : // scalar, global, md2macro, noGL
283 : template <> BaseIndex<3> I12::lowerBoundary{};
284 : template <> BaseIndex<3> I12::upperBoundary{};
285 : template <> tarch::la::Vector<3, unsigned int> I12::numberCellsInDomain{};
286 : template <> unsigned int I12::linearNumberCellsInDomain{};
287 : template <> tarch::la::Vector<3, unsigned int> I12::divisionFactor{};
288 : template <> const char I12::TNAME[] = "CellIndex<3, md2macro, noGhost>";
289 :
290 : // vector, global, md2macro, noGL
291 : template <> BaseIndex<3> I13::lowerBoundary{};
292 : template <> BaseIndex<3> I13::upperBoundary{};
293 : template <> tarch::la::Vector<3, unsigned int> I13::numberCellsInDomain{};
294 : template <> unsigned int I13::linearNumberCellsInDomain{};
295 : template <> tarch::la::Vector<3, unsigned int> I13::divisionFactor{};
296 : template <> const char I13::TNAME[] = "CellIndex<3, vector, md2macro, noGhost>";
297 :
298 : // scalar, local, md2macro, noGL
299 : template <> BaseIndex<3> I14::lowerBoundary{};
300 : template <> BaseIndex<3> I14::upperBoundary{};
301 : template <> tarch::la::Vector<3, unsigned int> I14::numberCellsInDomain{};
302 : template <> unsigned int I14::linearNumberCellsInDomain{};
303 : template <> tarch::la::Vector<3, unsigned int> I14::divisionFactor{};
304 : template <> const char I14::TNAME[] = "CellIndex<3, local, md2macro, noGhost>";
305 :
306 : // vector, local, md2macro, noGL
307 : template <> BaseIndex<3> I15::lowerBoundary{};
308 : template <> BaseIndex<3> I15::upperBoundary{};
309 : template <> tarch::la::Vector<3, unsigned int> I15::numberCellsInDomain{};
310 : template <> unsigned int I15::linearNumberCellsInDomain{};
311 : template <> tarch::la::Vector<3, unsigned int> I15::divisionFactor{};
312 : template <> const char I15::TNAME[] = "CellIndex<3, vector, local, md2macro, noGhost>";
313 : } // namespace indexing
314 : } // namespace coupling
315 :
316 : // delegated init, this does the main work
317 : template <unsigned int dim>
318 780 : void coupling::indexing::IndexingService<dim>::initWithCells(const tarch::la::Vector<dim, std::vector<unsigned int>>& subdomainWeights,
319 : tarch::la::Vector<dim, unsigned int> globalNumberCouplingCells,
320 : tarch::la::Vector<dim, unsigned int> numberProcesses,
321 : coupling::paralleltopology::ParallelTopologyType parallelTopologyType, unsigned int outerRegion,
322 : const unsigned int rank
323 : #if (COUPLING_MD_PARALLEL == COUPLING_MD_YES)
324 : ,
325 : MPI_Comm comm
326 : #endif
327 : ) {
328 :
329 780 : _rank = rank;
330 780 : _parallelTopologyType = parallelTopologyType;
331 : #if (COUPLING_MD_PARALLEL == COUPLING_MD_YES)
332 780 : _comm = comm;
333 : #endif
334 :
335 : #if (COUPLING_MD_ERROR == COUPLING_MD_YES)
336 780 : if (_isInitialized) {
337 48 : std::cout << "IndexingService: WARNING: Initializing twice! " << std::endl;
338 : }
339 : #endif
340 :
341 3080 : for (unsigned int d = 0; d < dim; d++) {
342 2316 : if (globalNumberCouplingCells[d] % numberProcesses[d] != 0) {
343 8 : std::stringstream ss;
344 8 : ss << "IndexingService: initWithCells(): ERROR: Number "
345 : "of coupling cells must be divisible by number of processes! ";
346 8 : ss << "globalNumberCouplingCells = " << globalNumberCouplingCells;
347 8 : ss << ", numberProcesses = " << numberProcesses;
348 16 : throw std::runtime_error(ss.str());
349 8 : }
350 :
351 2308 : const auto totalWeight = std::accumulate(subdomainWeights[d].begin(), subdomainWeights[d].end(), 0u);
352 2308 : if (globalNumberCouplingCells[d] % totalWeight != 0) {
353 0 : std::stringstream ss;
354 0 : ss << "IndexingService: initWithCells(): ERROR: Number "
355 : "of coupling cells must be divisible by total subdomain weights! ";
356 0 : ss << "globalNumberCouplingCells = " << globalNumberCouplingCells;
357 0 : ss << ", total weights for axis " << d << " = " << totalWeight;
358 0 : throw std::runtime_error(ss.str());
359 0 : }
360 :
361 2308 : if (numberProcesses[d] != subdomainWeights[d].size()) {
362 8 : std::stringstream ss;
363 8 : ss << "IndexingService: initWithCells(): ERROR: Number "
364 : "of process not equal to number of subdomain weights given! ";
365 8 : ss << "number of processes for axis " << d << " = " << numberProcesses[d];
366 8 : ss << ", number of subdomain weights for axis " << d << " = " << subdomainWeights[d].size();
367 16 : throw std::runtime_error(ss.str());
368 8 : }
369 : }
370 :
371 : // populate the _subdomainOwnership datastructure (does not include ghost cells)
372 : // after populating, can be used for lookup to find processes for specific cell coords
373 : // ex: if axis 0 has 8 cells, and weights 1,2,1, _subdomainOwnership[0] = {0,0,1,1,1,1,2,2}
374 : // hence cell 3 will be on x axis coord 1
375 3056 : for (unsigned int d = 0; d < dim; d++) {
376 : // with weights 1,2,1, totalWeight = 4
377 2292 : const int totalWeight = std::accumulate(subdomainWeights[d].begin(), subdomainWeights[d].end(), 0u);
378 : // with 8 cells on an axis and totalWeight = 4, multiplier = 2, thus weight 1 will have 2 cells
379 2292 : const int multiplier = globalNumberCouplingCells[d] / totalWeight;
380 2292 : _subdomainOwnership[d].clear();
381 2292 : _subdomainOwnership[d].reserve(globalNumberCouplingCells[d]);
382 6272 : for (unsigned int i = 0; i < numberProcesses[d]; i++) {
383 : // push in the appropriate number of coords into the lookup table
384 30860 : for (unsigned int j = 0; j < multiplier * subdomainWeights[d][i]; j++) {
385 26880 : _subdomainOwnership[d].push_back(i);
386 : }
387 : }
388 : }
389 :
390 : // TODO: make this globalNumberCouplingCells and remove all usages of the
391 : // old meaning (seen above)
392 1528 : const auto globalNumberCouplingCellsInclGL{globalNumberCouplingCells + tarch::la::Vector<dim, unsigned int>{2}};
393 :
394 : // init boundaries of all global, non-m2m, GL including indexing types
395 4584 : CellIndex<dim>::lowerBoundary = I01{0};
396 5348 : CellIndex<dim>::upperBoundary = I01{tarch::la::Vector<dim, int>{globalNumberCouplingCellsInclGL - tarch::la::Vector<dim, unsigned int>{1}}};
397 764 : CellIndex<dim>::setDomainParameters();
398 :
399 3056 : CellIndex<dim, IndexTrait::vector>::lowerBoundary = CellIndex<dim>::lowerBoundary;
400 764 : CellIndex<dim, IndexTrait::vector>::upperBoundary = CellIndex<dim>::upperBoundary;
401 764 : CellIndex<dim, IndexTrait::vector>::setDomainParameters();
402 :
403 : // init boundaries of all global, non-m2m, GL excluding indexing types
404 4584 : CellIndex<dim, IndexTrait::noGhost>::lowerBoundary = I01{1};
405 3056 : CellIndex<dim, IndexTrait::noGhost>::upperBoundary =
406 2292 : I01{tarch::la::Vector<dim, int>{globalNumberCouplingCellsInclGL - tarch::la::Vector<dim, unsigned int>{2}}};
407 764 : CellIndex<dim, IndexTrait::noGhost>::setDomainParameters();
408 :
409 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::noGhost>::lowerBoundary = CellIndex<dim, IndexTrait::noGhost>::lowerBoundary;
410 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::noGhost>::upperBoundary = CellIndex<dim, IndexTrait::noGhost>::upperBoundary;
411 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::noGhost>::setDomainParameters();
412 :
413 : // CellIndex<dim> and CellIndex<dim, IndexTrait::vector> have been set up, so from here on it is ok to do some basic conversions
414 : #if (COUPLING_MD_ERROR == COUPLING_MD_YES)
415 764 : _isInitialized = true;
416 : #endif
417 :
418 : // init boundaries of all global, m2m, GL excluding indexing types
419 :
420 4584 : CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary = BaseIndex<dim>{tarch::la::Vector<dim, int>{(int)(outerRegion + 1)}};
421 764 : CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary =
422 2292 : BaseIndex<dim>::upperBoundary - BaseIndex<dim>{tarch::la::Vector<dim, int>{(int)(outerRegion + 1)}};
423 764 : CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::setDomainParameters();
424 :
425 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary =
426 : CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary;
427 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary =
428 : CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary;
429 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro, IndexTrait::noGhost>::setDomainParameters();
430 :
431 : // init boundaries of all global, m2m, GL including indexing types
432 3056 : CellIndex<dim, IndexTrait::md2macro>::lowerBoundary =
433 3820 : I01{CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary.get() - tarch::la::Vector<dim, int>{1}};
434 3056 : CellIndex<dim, IndexTrait::md2macro>::upperBoundary =
435 3820 : I01{CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary.get() + tarch::la::Vector<dim, int>{1}};
436 764 : CellIndex<dim, IndexTrait::md2macro>::setDomainParameters();
437 :
438 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro>::lowerBoundary = CellIndex<dim, IndexTrait::md2macro>::lowerBoundary;
439 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro>::upperBoundary = CellIndex<dim, IndexTrait::md2macro>::upperBoundary;
440 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::md2macro>::setDomainParameters();
441 :
442 : // handle all local indexing types
443 :
444 764 : _numberProcesses = numberProcesses;
445 :
446 764 : _scalarNumberProcesses = _numberProcesses[0];
447 2292 : for (unsigned int d = 1; d < dim; d++)
448 1528 : _scalarNumberProcesses *= _numberProcesses[d];
449 :
450 1528 : _parallelTopology = coupling::paralleltopology::ParallelTopologyFactory::getParallelTopology<dim>(parallelTopologyType, _numberProcesses);
451 :
452 764 : const unsigned int topologyOffset = (_rank / _scalarNumberProcesses) * _scalarNumberProcesses;
453 764 : auto coords = _parallelTopology->getProcessCoordinates(_rank, topologyOffset);
454 : // init boundaries of all local, non-m2m, GL including indexing types
455 : {
456 764 : tarch::la::Vector<3, int> boxMin, boxMax;
457 3056 : for (unsigned int i = 0; i < dim; i++) {
458 : // calculate box bounds from cell ownership
459 : // find the first occurence of owned rank
460 2292 : boxMin[i] = std::distance(_subdomainOwnership[i].begin(), std::find(_subdomainOwnership[i].begin(), _subdomainOwnership[i].end(), coords[i]));
461 : // find the last occurence, add one to convert reverse iterator to forward
462 2292 : boxMax[i] =
463 2292 : std::distance(_subdomainOwnership[i].begin(), (std::find(_subdomainOwnership[i].rbegin(), _subdomainOwnership[i].rend(), coords[i]) + 1).base());
464 : }
465 : // _subdomainOwnership does not include ghost, so when we take the first occurence value, it is noGhost
466 : // however directly casting it into baseIndex shifts everything left by 1, since baseIndex expects ghost
467 : // thus we directly use the BaseIndex cast, and consequently boxMax requires +2 for two ghost layers
468 : // the "safer" way to write this would be I08{boxMin} -1 and I08{boxMax} + 1, but the result is the same
469 3056 : CellIndex<dim, IndexTrait::local>::lowerBoundary = BaseIndex<dim>{boxMin};
470 4584 : CellIndex<dim, IndexTrait::local>::upperBoundary = BaseIndex<dim>{boxMax + tarch::la::Vector<dim, int>{2}};
471 764 : CellIndex<dim, IndexTrait::local>::setDomainParameters();
472 : }
473 :
474 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::local>::lowerBoundary = CellIndex<dim, IndexTrait::local>::lowerBoundary;
475 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local>::upperBoundary = CellIndex<dim, IndexTrait::local>::upperBoundary;
476 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local>::setDomainParameters();
477 :
478 : // init boundaries of all local, non-m2m, GL excluding indexing types
479 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>::lowerBoundary =
480 3820 : I01{CellIndex<dim, IndexTrait::local>::lowerBoundary.get() + tarch::la::Vector<dim, int>{1}};
481 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>::upperBoundary =
482 3820 : I01{CellIndex<dim, IndexTrait::local>::upperBoundary.get() - tarch::la::Vector<dim, int>{1}};
483 764 : CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>::setDomainParameters();
484 :
485 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::lowerBoundary =
486 : CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>::lowerBoundary;
487 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::upperBoundary =
488 : CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>::upperBoundary;
489 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::noGhost>::setDomainParameters();
490 :
491 : // init boundaries of all local, m2m, GL excluding indexing types
492 : {
493 : using LocalNoGlIndex = CellIndex<dim, IndexTrait::local, IndexTrait::noGhost>;
494 : using GlobalMD2MacroNoGlIndex = CellIndex<dim, IndexTrait::md2macro, IndexTrait::noGhost>;
495 :
496 764 : tarch::la::Vector<dim, int> lowerBoundary;
497 764 : tarch::la::Vector<dim, int> upperBoundary;
498 3056 : for (unsigned int d = 0; d < dim; d++) {
499 5215 : lowerBoundary[d] = std::max(LocalNoGlIndex::lowerBoundary.get()[d], GlobalMD2MacroNoGlIndex::lowerBoundary.get()[d]);
500 9811 : upperBoundary[d] = std::min(LocalNoGlIndex::upperBoundary.get()[d], GlobalMD2MacroNoGlIndex::upperBoundary.get()[d]);
501 : }
502 :
503 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary = I01{lowerBoundary};
504 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary = I01{upperBoundary};
505 764 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::setDomainParameters();
506 : }
507 :
508 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary =
509 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary;
510 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary =
511 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary;
512 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::setDomainParameters();
513 :
514 : // init boundaries of all local, m2m, GL including indexing types
515 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro>::lowerBoundary =
516 3820 : I01{CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::lowerBoundary.get() - tarch::la::Vector<dim, int>{1}};
517 3056 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro>::upperBoundary =
518 3820 : I01{CellIndex<dim, IndexTrait::local, IndexTrait::md2macro, IndexTrait::noGhost>::upperBoundary.get() + tarch::la::Vector<dim, int>{1}};
519 764 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro>::setDomainParameters();
520 :
521 3056 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::lowerBoundary =
522 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro>::lowerBoundary;
523 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::upperBoundary =
524 : CellIndex<dim, IndexTrait::local, IndexTrait::md2macro>::upperBoundary;
525 764 : CellIndex<dim, IndexTrait::vector, IndexTrait::local, IndexTrait::md2macro>::setDomainParameters();
526 764 : }
527 :
528 : /*
529 : * This was in large parts stolen from IndexConversion.
530 : */
531 : template <unsigned int dim>
532 8 : std::vector<unsigned int> coupling::indexing::IndexingService<dim>::getRanksForGlobalIndex(const BaseIndex<dim>& globalCellIndex,
533 : unsigned int topologyOffset) const {
534 :
535 : #if (COUPLING_MD_ERROR == COUPLING_MD_YES)
536 8 : if (!_isInitialized) {
537 8 : throw std::runtime_error(std::string("coupling::indexing::getRanksForGlobalIndex: IndexingService not initialized! "));
538 : }
539 : #endif
540 :
541 4 : std::vector<unsigned int> ranks;
542 : // using the old meaning of 'globalNumberCouplingCells' from
543 : // IndexConversion
544 4 : const auto globalNumberCouplingCells = I08::numberCellsInDomain;
545 :
546 : // start and end coordinates of neighbouring cells.
547 4 : tarch::la::Vector<3, unsigned int> start(0);
548 4 : tarch::la::Vector<3, unsigned int> end(0);
549 4 : tarch::la::Vector<3, unsigned int> loopIndex(0);
550 :
551 : // determine up to 3^dim neighbouring cells in the surrounding of
552 : // globalCellIndex; reduce this number if globalCellIndex lies on the global
553 : // boundary
554 16 : for (unsigned int d = 0; d < dim; d++) {
555 12 : if ((unsigned int)globalCellIndex.get()[d] > 0) {
556 12 : start[d] = (unsigned int)globalCellIndex.get()[d] - 1;
557 : }
558 12 : end[d] = globalNumberCouplingCells[d] + 1;
559 12 : if ((unsigned int)globalCellIndex.get()[d] < end[d]) {
560 12 : end[d] = (unsigned int)globalCellIndex.get()[d] + 1;
561 : }
562 : }
563 :
564 : // loop over neighbouring regions
565 16 : for (loopIndex[2] = start[2]; loopIndex[2] <= end[2]; loopIndex[2]++) {
566 48 : for (loopIndex[1] = start[1]; loopIndex[1] <= end[1]; loopIndex[1]++) {
567 144 : for (loopIndex[0] = start[0]; loopIndex[0] <= end[0]; loopIndex[0]++) {
568 :
569 : // determine the global cell index of this particular grid cell
570 108 : tarch::la::Vector<dim, unsigned int> thisGlobalCellIndex(0);
571 432 : for (unsigned int d = 0; d < dim; d++) {
572 324 : thisGlobalCellIndex[d] = loopIndex[d];
573 : }
574 :
575 : // determine the unique rank for this cell
576 108 : const unsigned int rank = getUniqueRankForCouplingCell(thisGlobalCellIndex, globalNumberCouplingCells, topologyOffset);
577 :
578 : // add this rank to the vector with all ranks if we did not add this one
579 : // before
580 108 : bool isContained = false;
581 108 : const unsigned int thisSize = (unsigned int)ranks.size();
582 108 : for (unsigned int i = 0; i < thisSize; i++) {
583 104 : if (ranks[i] == rank) {
584 : isContained = true;
585 : break;
586 : }
587 : }
588 108 : if (!isContained) {
589 4 : ranks.push_back(rank);
590 : }
591 : }
592 : }
593 : }
594 :
595 4 : return ranks;
596 0 : }
597 :
598 : /*
599 : * This was in large parts stolen from IndexConversion.
600 : * Note that this uses the globalNumberCouplingCells definition excl. the
601 : * ghost layer.
602 : */
603 : template <unsigned int dim>
604 120 : unsigned int coupling::indexing::IndexingService<dim>::getUniqueRankForCouplingCell(tarch::la::Vector<dim, unsigned int> globalCellIndex,
605 : const tarch::la::Vector<dim, unsigned int>& globalNumberCouplingCells,
606 : unsigned int topologyOffset) const {
607 : // vector containing avg number of macro cells, not counting global GL.
608 472 : for (unsigned int d = 0; d < dim; d++) {
609 356 : if (globalCellIndex[d] >= globalNumberCouplingCells[d] + 2) { // greater or equal to the total global number incl GL (+2)
610 : using namespace std::string_literals;
611 8 : throw std::runtime_error("IndexingService: getUniqueRankForCouplingCell(): Global cell index greater than global size in dim "s + std::to_string(d));
612 : }
613 352 : if (globalNumberCouplingCells[d] % _numberProcesses[d] != 0) {
614 0 : std::stringstream ss;
615 0 : ss << "IndexingService: getUniqueRankForCouplingCell(): ERROR: Number "
616 : "of coupling cells must be divisible by number of processes! ";
617 0 : ss << "globalNumberCouplingCells = " << globalNumberCouplingCells;
618 0 : ss << ", numberProcesses = " << _numberProcesses;
619 0 : throw std::runtime_error(ss.str());
620 0 : }
621 : }
622 :
623 116 : tarch::la::Vector<dim, unsigned int> processCoords(0);
624 464 : for (unsigned int d = 0; d < dim; d++) {
625 : // the passed index includes ghost cells, while _subdomainOwnership is noGhost
626 : // thus, we have two edge cases (the two ghost cells, number 0 and size-1)
627 : // we handle these two cases first, and then use the lookup table for everything else
628 348 : if (globalCellIndex[d] == 0) // first ghost, return first coord in this axis (usually 0)
629 108 : processCoords[d] = _subdomainOwnership[d][0];
630 240 : else if (globalCellIndex[d] == globalNumberCouplingCells[d] + 1) // last ghost, return last coord
631 0 : processCoords[d] = _subdomainOwnership[d][_subdomainOwnership[d].size() - 1];
632 : else // all internal cells
633 240 : processCoords[d] = _subdomainOwnership[d][globalCellIndex[d] - 1];
634 : }
635 232 : return _parallelTopology->getRank(processCoords, topologyOffset);
636 : }
637 :
638 : template <unsigned int dim>
639 0 : unsigned int coupling::indexing::IndexingService<dim>::getUniqueRankForCouplingCell(const BaseIndex<dim>& globalCellIndex, unsigned int topologyOffset) const {
640 :
641 : #if (COUPLING_MD_ERROR == COUPLING_MD_YES)
642 0 : if (!_isInitialized) {
643 0 : throw std::runtime_error(std::string("coupling::indexing::getUniqueRankForCouplingCell: IndexingService not initialized! "));
644 : }
645 : #endif
646 :
647 0 : return getUniqueRankForCouplingCell((tarch::la::Vector<dim, unsigned int>)(globalCellIndex.get()), I09::numberCellsInDomain, topologyOffset);
648 : }
649 :
650 : // declare specialisation of IndexingService
651 : #ifdef INDEXING_ENABLE_DIM2
652 : template class coupling::indexing::IndexingService<2>;
653 : #endif
654 : template class coupling::indexing::IndexingService<3>;
|
