preCICE v3.1.2
Loading...
Searching...
No Matches
PartitionOfUnityMappingTest.cpp
Go to the documentation of this file.
1#include <Eigen/Core>
2#include <algorithm>
3#include <memory>
4#include <ostream>
5#include <string>
6#include <utility>
7#include <vector>
8#include "logging/Logger.hpp"
9#include "mapping/Mapping.hpp"
10#include "mapping/PartitionOfUnityMapping.hpp"
11#include "mapping/impl/BasisFunctions.hpp"
12#include "mesh/Data.hpp"
13#include "mesh/Mesh.hpp"
14#include "mesh/SharedPointer.hpp"
15#include "mesh/Utils.hpp"
16#include "mesh/Vertex.hpp"
17#include "testing/TestContext.hpp"
18#include "testing/Testing.hpp"
19
20using namespace precice;
21using namespace precice::mesh;
22using namespace precice::mapping;
23using namespace precice::testing;
24using precice::testing::TestContext;
25
26BOOST_AUTO_TEST_SUITE(MappingTests)
27BOOST_AUTO_TEST_SUITE(PartitionOfUnityMapping)
28
29void addGlobalIndex(mesh::PtrMesh &mesh, int offset = 0)
30{
31 for (mesh::Vertex &v : mesh->vertices()) {
32 v.setGlobalIndex(v.getID() + offset);
33 }
34}
35
36double sumComponentWise(const Eigen::VectorXd &vec, int component, int dataDimension)
37{
38 PRECICE_ASSERT(component < dataDimension);
39 double sum = 0;
40 for (unsigned int i = 0; i < vec.size() / dataDimension; ++i) {
41 sum += vec[(i * dataDimension) + component];
42 }
43 return sum;
44}
45
46BOOST_AUTO_TEST_SUITE(Serial)
47
48void perform2DTestConsistentMapping(Mapping &mapping)
49{
50 int dimensions = 2;
51 using Eigen::Vector2d;
52
53 // Create mesh to map from
54 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
55 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
56 int inDataID = inData->getID();
57 inMesh->createVertex(Vector2d(0.0, 0.0));
58 inMesh->createVertex(Vector2d(1.0, 0.0));
59 inMesh->createVertex(Vector2d(1.0, 1.0));
60 inMesh->createVertex(Vector2d(0.0, 1.0));
61
62 inMesh->createVertex(Vector2d(2.0, 0.0));
63 inMesh->createVertex(Vector2d(3.0, 0.0));
64 inMesh->createVertex(Vector2d(3.0, 1.0));
65 inMesh->createVertex(Vector2d(2.0, 1.0));
66
67 inMesh->createVertex(Vector2d(4.0, 0.0));
68 inMesh->createVertex(Vector2d(5.0, 0.0));
69 inMesh->createVertex(Vector2d(5.0, 1.0));
70 inMesh->createVertex(Vector2d(4.0, 1.0));
71
72 inMesh->allocateDataValues();
73 addGlobalIndex(inMesh);
74
75 auto &values = inData->values();
76 values << 1.0, 2.0, 2.0, 1.0, 3.0, 4.0, 4.0, 3.0, 5.0, 6.0, 6.0, 5.0;
77
78 // Create mesh to map to
79 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
80 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
81 int outDataID = outData->getID();
82 mesh::Vertex &vertex = outMesh->createVertex(Vector2d(0, 0));
83 mesh::Vertex &vertex1 = outMesh->createVertex(Vector2d(3.5, 0.5));
84 mesh::Vertex &vertex2 = outMesh->createVertex(Vector2d(2.5, 0.5));
85 // vertex will be changed
86 outMesh->createVertex(Vector2d(0, 0));
87 outMesh->createVertex(Vector2d(5, 1));
88 outMesh->allocateDataValues();
89 addGlobalIndex(outMesh);
90
91 // Setup mapping with mapping coordinates and geometry used
92 mapping.setMeshes(inMesh, outMesh);
93 BOOST_TEST(mapping.hasComputedMapping() == false);
94
95 vertex.setCoords(Vector2d(0.0, 0.0));
96 mapping.computeMapping();
97 mapping.map(inDataID, outDataID);
98 double value = outData->values()(0);
99 double value1 = outData->values()(1);
100 double value2 = outData->values()(2);
101 BOOST_TEST(mapping.hasComputedMapping() == true);
102 BOOST_TEST(value == 1.0);
103 BOOST_TEST(value1 == 4.5);
104 BOOST_TEST(value2 == 3.5);
105
106 vertex.setCoords(Vector2d(0.0, 0.5));
107 mapping.clear();
108 outData->values().setZero();
109 BOOST_TEST(mapping.hasComputedMapping() == false);
110 mapping.computeMapping();
111 mapping.map(inDataID, outDataID);
112 value = outData->values()(0);
113 BOOST_TEST(mapping.hasComputedMapping() == true);
114 BOOST_TEST(value == 1.0);
115
116 vertex.setCoords(Vector2d(0.0, 1.0));
117 mapping.clear();
118 outData->values().setZero();
119 mapping.computeMapping();
120 mapping.map(inDataID, outDataID);
121 value = outData->values()(0);
122 BOOST_TEST(mapping.hasComputedMapping() == true);
123 BOOST_TEST(value == 1.0);
124
125 vertex.setCoords(Vector2d(1.0, 0.0));
126 mapping.clear();
127 outData->values().setZero();
128 mapping.computeMapping();
129 mapping.map(inDataID, outDataID);
130 value = outData->values()(0);
131 BOOST_TEST(mapping.hasComputedMapping() == true);
132 BOOST_TEST(value == 2.0);
133
134 vertex.setCoords(Vector2d(1.0, 0.5));
135 mapping.clear();
136 outData->values().setZero();
137 mapping.computeMapping();
138 mapping.map(inDataID, outDataID);
139 value = outData->values()(0);
140 BOOST_TEST(mapping.hasComputedMapping() == true);
141 BOOST_TEST(value == 2.0);
142
143 vertex.setCoords(Vector2d(1.0, 1.0));
144 mapping.clear();
145 outData->values().setZero();
146 mapping.computeMapping();
147 mapping.map(inDataID, outDataID);
148 value = outData->values()(0);
149 BOOST_TEST(mapping.hasComputedMapping() == true);
150 BOOST_TEST(value == 2.0);
151
152 vertex.setCoords(Vector2d(0.5, 0.0));
153 mapping.clear();
154 outData->values().setZero();
155 mapping.computeMapping();
156 mapping.map(inDataID, outDataID);
157 value = outData->values()(0);
158 BOOST_TEST(mapping.hasComputedMapping() == true);
159 BOOST_TEST(value == 1.5);
160
161 vertex.setCoords(Vector2d(0.5, 0.5));
162 mapping.clear();
163 outData->values().setZero();
164 mapping.computeMapping();
165 mapping.map(inDataID, outDataID);
166 value = outData->values()(0);
167 BOOST_TEST(mapping.hasComputedMapping() == true);
168 BOOST_TEST(value == 1.5);
169
170 vertex.setCoords(Vector2d(0.5, 1.0));
171 mapping.clear();
172 outData->values().setZero();
173 mapping.computeMapping();
174 mapping.map(inDataID, outDataID);
175 value = outData->values()(0);
176 BOOST_TEST(mapping.hasComputedMapping() == true);
177 BOOST_TEST(value == 1.5);
178}
179
180void perform2DTestConsistentMappingVector(Mapping &mapping)
181{
182 int dimensions = 2;
183 int dataDimensions = dimensions;
184 using Eigen::Vector2d;
185
186 // Create mesh to map from
187 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
188 mesh::PtrData inData = inMesh->createData("InData", dataDimensions, 0_dataID);
189 int inDataID = inData->getID();
190 inMesh->createVertex(Vector2d(0.0, 0.0));
191 inMesh->createVertex(Vector2d(1.0, 0.0));
192 inMesh->createVertex(Vector2d(1.0, 1.0));
193 inMesh->createVertex(Vector2d(0.0, 1.0));
194
195 inMesh->createVertex(Vector2d(2.0, 0.0));
196 inMesh->createVertex(Vector2d(3.0, 0.0));
197 inMesh->createVertex(Vector2d(3.0, 1.0));
198 inMesh->createVertex(Vector2d(2.0, 1.0));
199
200 inMesh->createVertex(Vector2d(4.0, 0.0));
201 inMesh->createVertex(Vector2d(5.0, 0.0));
202 inMesh->createVertex(Vector2d(5.0, 1.0));
203 inMesh->createVertex(Vector2d(4.0, 1.0));
204
205 inMesh->allocateDataValues();
206 addGlobalIndex(inMesh);
207
208 auto &values = inData->values();
209 // We select the second component = 2 * first component
210 values << 1.0, 2.0, 2.0, 4.0, 2.0, 4.0, 1.0, 2.0, 3.0, 6.0, 4.0, 8.0, 4.0, 8.0, 3.0, 6.0, 5.0, 10.0, 6.0, 12.0, 6.0, 12.0, 5.0, 10.0;
211
212 // Create mesh to map to
213 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
214 mesh::PtrData outData = outMesh->createData("OutData", dataDimensions, 1_dataID);
215 int outDataID = outData->getID();
216 mesh::Vertex &vertex = outMesh->createVertex(Vector2d(0, 0));
217 mesh::Vertex &vertex1 = outMesh->createVertex(Vector2d(3.5, 0.5));
218 mesh::Vertex &vertex2 = outMesh->createVertex(Vector2d(2.5, 0.5));
219 // vertex will be changed
220 outMesh->createVertex(Vector2d(0, 0));
221 outMesh->createVertex(Vector2d(5, 1));
222 outMesh->allocateDataValues();
223 addGlobalIndex(outMesh);
224
225 // Setup mapping with mapping coordinates and geometry used
226 mapping.setMeshes(inMesh, outMesh);
227 BOOST_TEST(mapping.hasComputedMapping() == false);
228
229 vertex.setCoords(Vector2d(0.0, 0.0));
230 mapping.computeMapping();
231 mapping.map(inDataID, outDataID);
232 DataID index = 0;
233 DataID index1 = 1 * dataDimensions;
234 DataID index2 = 2 * dataDimensions;
235 BOOST_TEST(mapping.hasComputedMapping() == true);
236 BOOST_TEST(outData->values()(index) == 1.0);
237 BOOST_TEST(outData->values()(index + 1) == 2.0);
238 BOOST_TEST(outData->values()(index1) == 4.5);
239 BOOST_TEST(outData->values()(index1 + 1) == 9.0);
240 BOOST_TEST(outData->values()(index2) == 3.5);
241 BOOST_TEST(outData->values()(index2 + 1) == 7.0);
242
243 vertex.setCoords(Vector2d(0.0, 0.5));
244 mapping.clear();
245 outData->values().setZero();
246 mapping.computeMapping();
247 mapping.map(inDataID, outDataID);
248 BOOST_TEST(mapping.hasComputedMapping() == true);
249 BOOST_TEST(outData->values()(0) == 1.0);
250 BOOST_TEST(outData->values()(1) == 2.0);
251
252 vertex.setCoords(Vector2d(0.0, 1.0));
253 mapping.clear();
254 outData->values().setZero();
255 mapping.computeMapping();
256 mapping.map(inDataID, outDataID);
257 BOOST_TEST(mapping.hasComputedMapping() == true);
258 BOOST_TEST(outData->values()(0) == 1.0);
259 BOOST_TEST(outData->values()(1) == 2.0);
260
261 vertex.setCoords(Vector2d(1.0, 0.0));
262 mapping.clear();
263 outData->values().setZero();
264 mapping.computeMapping();
265 mapping.map(inDataID, outDataID);
266 BOOST_TEST(mapping.hasComputedMapping() == true);
267 BOOST_TEST(outData->values()(0) == 2.0);
268 BOOST_TEST(outData->values()(1) == 4.0);
269
270 vertex.setCoords(Vector2d(1.0, 0.5));
271 mapping.clear();
272 outData->values().setZero();
273 mapping.computeMapping();
274 mapping.map(inDataID, outDataID);
275 BOOST_TEST(mapping.hasComputedMapping() == true);
276 BOOST_TEST(outData->values()(0) == 2.0);
277 BOOST_TEST(outData->values()(1) == 4.0);
278
279 vertex.setCoords(Vector2d(1.0, 1.0));
280 mapping.clear();
281 outData->values().setZero();
282 mapping.computeMapping();
283 mapping.map(inDataID, outDataID);
284 BOOST_TEST(mapping.hasComputedMapping() == true);
285 BOOST_TEST(outData->values()(0) == 2.0);
286 BOOST_TEST(outData->values()(1) == 4.0);
287
288 vertex.setCoords(Vector2d(0.5, 0.0));
289 mapping.clear();
290 outData->values().setZero();
291 mapping.computeMapping();
292 mapping.map(inDataID, outDataID);
293 BOOST_TEST(mapping.hasComputedMapping() == true);
294 BOOST_TEST(outData->values()(0) == 1.5);
295 BOOST_TEST(outData->values()(1) == 3.0);
296
297 vertex.setCoords(Vector2d(0.5, 0.5));
298 mapping.clear();
299 outData->values().setZero();
300 mapping.computeMapping();
301 mapping.map(inDataID, outDataID);
302 BOOST_TEST(mapping.hasComputedMapping() == true);
303 BOOST_TEST(outData->values()(0) == 1.5);
304 BOOST_TEST(outData->values()(1) == 3.0);
305
306 vertex.setCoords(Vector2d(0.5, 1.0));
307 mapping.clear();
308 outData->values().setZero();
309 mapping.computeMapping();
310 mapping.map(inDataID, outDataID);
311 BOOST_TEST(mapping.hasComputedMapping() == true);
312 BOOST_TEST(outData->values()(0) == 1.5);
313 BOOST_TEST(outData->values()(1) == 3.0);
314}
315
316// Tests the automatic reduction of dead axis
317void performTestConsistentMapDeadAxis(Mapping &mapping, int dim)
318{
319 // Create mesh to map from
320 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dim, testing::nextMeshID()));
321 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
322 int inDataID = inData->getID();
323
324 // create 20 vertices with two very close lines
325 for (unsigned int i = 0; i < 20; ++i)
326 if (dim == 2) {
327 inMesh->createVertex(Eigen::Vector2d(0.0, 0.0 + i * dim));
328 inMesh->createVertex(Eigen::Vector2d(1e-5, 0.0 + i * dim));
329 } else {
330 inMesh->createVertex(Eigen::Vector3d(7.0, 7.0, 40 + i * dim));
331 inMesh->createVertex(Eigen::Vector3d(7.001, 7.001, 40 + i * dim));
332 }
333
334 inMesh->allocateDataValues();
335 addGlobalIndex(inMesh);
336
337 auto &values = inData->values();
338 for (std::size_t v = 0; v < 20; ++v) {
339 values[v * 2] = v * dim * 1e-5;
340 values[v * 2 + 1] = v * dim * 1e-5 + 2;
341 }
342
343 // Create mesh to map to
344 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dim, testing::nextMeshID()));
345 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
346 int outDataID = outData->getID();
347
348 for (unsigned int i = 0; i < 15; ++i) {
349 if (dim == 2)
350 outMesh->createVertex(Eigen::Vector2d(.1, 1 + 1 + i * dim));
351 else {
352 outMesh->createVertex(Eigen::Vector3d(7.4, 7.4, 41 + i * dim));
353 }
354 }
355 outMesh->allocateDataValues();
356 addGlobalIndex(outMesh);
357
358 // Setup mapping with mapping coordinates and geometry used
359 mapping.setMeshes(inMesh, outMesh);
360 BOOST_TEST(mapping.hasComputedMapping() == false);
361
362 mapping.computeMapping();
363 mapping.map(inDataID, outDataID);
364 double value = outData->values()(0);
365 double value1 = outData->values()(1);
366 double value2 = outData->values()(2);
367 BOOST_TEST(mapping.hasComputedMapping() == true);
368 if (dim == 2) {
369 BOOST_TEST(math::equals(value, 19935.374757794027, 0.2));
370 BOOST_TEST(math::equals(value1, 19935.37795225389, 0.2));
371 BOOST_TEST(math::equals(value2, 19935.35164938603, 0.2));
372 } else {
373 BOOST_TEST(value == 829.28063055069435);
374 BOOST_TEST(value1 == 819.35577218983303);
375 BOOST_TEST(value2 == 829.38388713302811);
376 }
377}
378
379void perform2DTestConservativeMapping(Mapping &mapping)
380{
381 int dimensions = 2;
382 using Eigen::Vector2d;
383
384 // Create mesh to map from
385 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
386 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
387 int inDataID = inData->getID();
388
389 // vertices will be changed
390 mesh::Vertex &vertex = inMesh->createVertex(Vector2d(1.0, 0.0));
391 mesh::Vertex &vertex1 = inMesh->createVertex(Vector2d(3.5, 0.0));
392 mesh::Vertex &vertex2 = inMesh->createVertex(Vector2d(2.5, 0.5));
393
394 inMesh->createVertex(Vector2d(0, 0));
395 inMesh->createVertex(Vector2d(5, 1));
396 inMesh->allocateDataValues();
397 addGlobalIndex(inMesh);
398 inMesh->setGlobalNumberOfVertices(inMesh->nVertices());
399
400 auto &values = inData->values();
401 values << 10.0, 0.0, 0.0, 10.0, 0.0;
402 // Create mesh to map to
403 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
404 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
405 int outDataID = outData->getID();
406
407 outMesh->createVertex(Vector2d(0.0, 0.0));
408 outMesh->createVertex(Vector2d(1.0, 0.0));
409 outMesh->createVertex(Vector2d(1.0, 1.0));
410 outMesh->createVertex(Vector2d(0.0, 1.0));
411
412 outMesh->createVertex(Vector2d(2.0, 0.0));
413 outMesh->createVertex(Vector2d(3.0, 0.0));
414 outMesh->createVertex(Vector2d(3.0, 1.0));
415 outMesh->createVertex(Vector2d(2.0, 1.0));
416
417 outMesh->createVertex(Vector2d(4.0, 0.0));
418 outMesh->createVertex(Vector2d(5.0, 0.0));
419 outMesh->createVertex(Vector2d(5.0, 1.0));
420 outMesh->createVertex(Vector2d(4.0, 1.0));
421
422 outMesh->allocateDataValues();
423 addGlobalIndex(outMesh);
424 outMesh->setGlobalNumberOfVertices(outMesh->nVertices());
425
426 // Setup mapping with mapping coordinates and geometry used
427 mapping.setMeshes(inMesh, outMesh);
428 BOOST_TEST(mapping.hasComputedMapping() == false);
429
430 // Test that we get point-wise exact data values on output vertices 1 and 2
431 // remaining vertices are zero
432 mapping.computeMapping();
433 mapping.map(inDataID, outDataID);
434
435 BOOST_TEST(outData->values().sum() == inData->values().sum());
436 BOOST_TEST(outData->values()(0) == 10);
437 BOOST_TEST(outData->values()(1) == 10);
438 BOOST_TEST(outData->values()(2) == 0);
439 BOOST_TEST(outData->values()(3) == 0);
440
441 // Test that we have a symmetric solution if we have a force acting
442 // in between two output vertices (here 5 and 8), there is actually
443 // no guarantee, but it holds for this setup
444 values << 0.0, 10.0, 0.0, 0.0, 0.0;
445
446 mapping.clear();
447 outData->values().setZero();
448 mapping.computeMapping();
449 mapping.map(inDataID, outDataID);
450 BOOST_TEST(outData->values().sum() == inData->values().sum());
451 BOOST_TEST(outData->values()(5) == outData->values()(8));
452
453 // Test the conservation property if we have everywhere non-zero input data
454 values << 5.0, 10.0, 7.0, 3.0, 4.0;
455
456 mapping.clear();
457 outData->values().setZero();
458 mapping.computeMapping();
459 mapping.map(inDataID, outDataID);
460
461 BOOST_TEST(outData->values().sum() == inData->values().sum());
462 BOOST_TEST(mapping.hasComputedMapping() == true);
463
464 // Test the conservation property if we have everywhere non-zero input data
465 values << 3.0, 4.0, 5.0, 7.0, 9.0;
466
467 mapping.clear();
468 outData->values().setZero();
469 mapping.computeMapping();
470 mapping.map(inDataID, outDataID);
471 BOOST_TEST(outData->values().sum() == inData->values().sum());
472 BOOST_TEST(mapping.hasComputedMapping() == true);
473 BOOST_TEST(outData->values()(5) == 3.2931869752057001);
474}
475
476void perform2DTestConservativeMappingVector(Mapping &mapping)
477{
478 int dimensions = 2;
479 int dataDimensions = dimensions;
480 using Eigen::Vector2d;
481
482 // Create mesh to map to
483 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
484 mesh::PtrData inData = inMesh->createData("InData", dataDimensions, 0_dataID);
485 DataID inDataID = inData->getID();
486 mesh::Vertex &vertex = inMesh->createVertex(Vector2d(0, 0));
487 mesh::Vertex &vertex1 = inMesh->createVertex(Vector2d(3.5, 0.5));
488 mesh::Vertex &vertex2 = inMesh->createVertex(Vector2d(2.5, 0.5));
489 // vertex will be changed
490 inMesh->createVertex(Vector2d(0, 0));
491 inMesh->createVertex(Vector2d(5, 1));
492
493 inMesh->allocateDataValues();
494 addGlobalIndex(inMesh);
495
496 auto &values = inData->values();
497 // We select the second component = 2 * first component
498 values << 1.0, 2.0, 2.0, 4.0, 2.0, 4.0, 1.0, 2.0, 3.0, 6.0;
499
500 // Create mesh to map from
501 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
502 mesh::PtrData outData = outMesh->createData("OutData", dataDimensions, 1_dataID);
503 DataID outDataID = outData->getID();
504
505 outMesh->createVertex(Vector2d(0.0, 0.0));
506 outMesh->createVertex(Vector2d(1.0, 0.0));
507 outMesh->createVertex(Vector2d(1.0, 1.0));
508 outMesh->createVertex(Vector2d(0.0, 1.0));
509
510 outMesh->createVertex(Vector2d(2.0, 0.0));
511 outMesh->createVertex(Vector2d(3.0, 0.0));
512 outMesh->createVertex(Vector2d(3.0, 1.0));
513 outMesh->createVertex(Vector2d(2.0, 1.0));
514
515 outMesh->createVertex(Vector2d(4.0, 0.0));
516 outMesh->createVertex(Vector2d(5.0, 0.0));
517 outMesh->createVertex(Vector2d(5.0, 1.0));
518 outMesh->createVertex(Vector2d(4.0, 1.0));
519
520 outMesh->allocateDataValues();
521 addGlobalIndex(outMesh);
522
523 // Setup mapping with mapping coordinates and geometry used
524 mapping.setMeshes(inMesh, outMesh);
525 BOOST_TEST(mapping.hasComputedMapping() == false);
526
527 // Check that we preserve data component-wise
528 vertex.setCoords(Vector2d(1.0, 0.0));
529 mapping.computeMapping();
530 mapping.map(inDataID, outDataID);
531
532 // component 0
533 double insumc0 = sumComponentWise(inData->values(), 0, dataDimensions);
534 double outsumc0 = sumComponentWise(outData->values(), 0, dataDimensions);
535 BOOST_TEST(insumc0 == outsumc0);
536
537 // component 1
538 double insumc1 = sumComponentWise(inData->values(), 1, dataDimensions);
539 double outsumc1 = sumComponentWise(outData->values(), 1, dataDimensions);
540
541 BOOST_TEST(insumc0 * 2 == insumc1);
542 BOOST_TEST(outsumc0 * 2 == outsumc1);
543 BOOST_TEST(insumc1 == outsumc1);
544
545 BOOST_TEST(mapping.hasComputedMapping() == true);
546
547 // Check for specific sum values for each component
548 values << 1.0, 0.0, 12.0, 0.0, 2.0, 0.0, 1.0, 0.0, 3.0, 0.0;
549
550 mapping.clear();
551 outData->values().setZero();
552 mapping.computeMapping();
553 mapping.map(inDataID, outDataID);
554
555 // component 0
556 outsumc0 = sumComponentWise(outData->values(), 0, dataDimensions);
557 BOOST_TEST(19.0 == outsumc0);
558
559 // we expect no contribution for the zeroth component
560 // component 1
561 outsumc1 = sumComponentWise(outData->values(), 1, dataDimensions);
562 BOOST_TEST(0.0 == outsumc1);
563
564 BOOST_TEST(mapping.hasComputedMapping() == true);
565
566 // Check for specific sum values for each component
567 values << 0.0, 2.0, 0.0, 4.0, 0.0, 8.0, 0.0, 7.0, 0.0, 0.0;
568
569 mapping.clear();
570 outData->values().setZero();
571 mapping.computeMapping();
572 mapping.map(inDataID, outDataID);
573
574 // component 0
575 outsumc0 = sumComponentWise(outData->values(), 0, dataDimensions);
576 BOOST_TEST(0.0 == outsumc0);
577
578 // we expect no contribution for the zeroth component
579 // component 1
580 outsumc1 = sumComponentWise(outData->values(), 1, dataDimensions);
581 BOOST_TEST(21.0 == outsumc1);
582
583 BOOST_TEST(mapping.hasComputedMapping() == true);
584
585 // Check for the exact reproduction of matching vertices
586 values << 10.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 27.0;
587
588 mapping.clear();
589 outData->values().setZero();
590 mapping.computeMapping();
591 mapping.map(inDataID, outDataID);
592
593 // component 0
594 outsumc0 = sumComponentWise(outData->values(), 0, dataDimensions);
595 BOOST_TEST(10.0 == outsumc0);
596 BOOST_TEST(outData->values()(2) == 10);
597
598 // component 1
599 outsumc1 = sumComponentWise(outData->values(), 1, dataDimensions);
600 BOOST_TEST(27.0 == outsumc1);
601 BOOST_TEST(outData->values()(21) == 27);
602
603 BOOST_TEST(mapping.hasComputedMapping() == true);
604
605 // Check for the interpolation at some vertex
606 values << 3.0, 6.0, 2.0, 4.0, 12.0, 24.0, 1.0, 2.0, 3.0, 6.0;
607
608 mapping.clear();
609 outData->values().setZero();
610 mapping.computeMapping();
611 mapping.map(inDataID, outDataID);
612
613 BOOST_TEST(outData->values().sum() == inData->values().sum());
614
615 double expectedValue = 3.5933508322619825;
616 // component 0
617 BOOST_TEST(outData->values()(12) == expectedValue);
618
619 // component 1
620 BOOST_TEST(outData->values()(13) == 2 * expectedValue);
621}
622
623void perform3DTestConsistentMapping(Mapping &mapping)
624{
625 int dimensions = 3;
626 using Eigen::Vector3d;
627
628 // Create mesh to map from
629 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
630 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
631 int inDataID = inData->getID();
632 inMesh->createVertex(Vector3d(0.0, 0.0, 0.0));
633 inMesh->createVertex(Vector3d(1.0, 0.0, 0.0));
634 inMesh->createVertex(Vector3d(1.0, 1.0, 0.0));
635 inMesh->createVertex(Vector3d(0.0, 1.0, 0.0));
636
637 inMesh->createVertex(Vector3d(0.0, 0.0, 1.0));
638 inMesh->createVertex(Vector3d(1.0, 0.0, 1.0));
639 inMesh->createVertex(Vector3d(1.0, 1.0, 1.0));
640 inMesh->createVertex(Vector3d(0.0, 1.0, 1.0));
641
642 inMesh->createVertex(Vector3d(2.0, 0.0, 0.0));
643 inMesh->createVertex(Vector3d(3.0, 0.0, 0.0));
644 inMesh->createVertex(Vector3d(3.0, 1.0, 0.0));
645 inMesh->createVertex(Vector3d(2.0, 1.0, 0.0));
646
647 inMesh->createVertex(Vector3d(2.0, 0.0, 1.0));
648 inMesh->createVertex(Vector3d(3.0, 0.0, 1.0));
649 inMesh->createVertex(Vector3d(3.0, 1.0, 1.0));
650 inMesh->createVertex(Vector3d(2.0, 1.0, 1.0));
651
652 inMesh->createVertex(Vector3d(4.0, 0.0, 0.0));
653 inMesh->createVertex(Vector3d(5.0, 0.0, 0.0));
654 inMesh->createVertex(Vector3d(5.0, 1.0, 0.0));
655 inMesh->createVertex(Vector3d(4.0, 1.0, 0.0));
656
657 inMesh->createVertex(Vector3d(4.0, 0.0, 1.0));
658 inMesh->createVertex(Vector3d(5.0, 0.0, 1.0));
659 inMesh->createVertex(Vector3d(5.0, 1.0, 1.0));
660 inMesh->createVertex(Vector3d(4.0, 1.0, 1.0));
661
662 inMesh->allocateDataValues();
663 addGlobalIndex(inMesh);
664
665 auto &values = inData->values();
666 // Set the values in the parallel (z) plane 3*values
667 values << 1.0, 2.0, 2.0, 1.0, 3.0, 6.0, 6.0, 3.0, 3.0, 4.0, 4.0, 3.0, 9.0, 12.0, 12.0, 9.0, 5.0, 6.0, 6.0, 5.0, 15.0, 18.0, 18.0, 15.0;
668
669 // Create mesh to map to
670 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
671 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
672 int outDataID = outData->getID();
673 mesh::Vertex &vertex = outMesh->createVertex(Vector3d(0, 0, 0));
674 mesh::Vertex &vertex1 = outMesh->createVertex(Vector3d(3.5, 0.5, 0.5));
675 mesh::Vertex &vertex2 = outMesh->createVertex(Vector3d(2.5, 0.5, 1.0));
676 outMesh->createVertex(Vector3d(0, 0, 0.5));
677 outMesh->createVertex(Vector3d(5, 1, 0.5));
678 outMesh->allocateDataValues();
679 addGlobalIndex(outMesh);
680
681 // Setup mapping with mapping coordinates and geometry used
682 mapping.setMeshes(inMesh, outMesh);
683 BOOST_TEST(mapping.hasComputedMapping() == false);
684
685 vertex.setCoords(Vector3d(0.0, 0.0, 0.0));
686 mapping.computeMapping();
687 mapping.map(inDataID, outDataID);
688 double value = outData->values()(0);
689 double value1 = outData->values()(1);
690 double value2 = outData->values()(2);
691 BOOST_TEST(mapping.hasComputedMapping() == true);
692 BOOST_TEST(value == 1.0);
693 BOOST_TEST(value1 == 9.0);
694 BOOST_TEST(value2 == 10.5);
695
696 vertex.setCoords(Vector3d(0.0, 0.5, 0.5));
697 mapping.clear();
698 outData->values().setZero();
699 mapping.computeMapping();
700 mapping.map(inDataID, outDataID);
701 value = outData->values()(0);
702 BOOST_TEST(mapping.hasComputedMapping() == true);
703 BOOST_TEST(value == 2.0);
704
705 vertex.setCoords(Vector3d(0.0, 1.0, 1.0));
706 mapping.clear();
707 outData->values().setZero();
708 mapping.computeMapping();
709 mapping.map(inDataID, outDataID);
710 value = outData->values()(0);
711 BOOST_TEST(mapping.hasComputedMapping() == true);
712 BOOST_TEST(value == 3.0);
713
714 vertex.setCoords(Vector3d(1.0, 0.0, 0.0));
715 mapping.clear();
716 outData->values().setZero();
717 mapping.computeMapping();
718 mapping.map(inDataID, outDataID);
719 value = outData->values()(0);
720 BOOST_TEST(mapping.hasComputedMapping() == true);
721 BOOST_TEST(value == 2.0);
722
723 vertex.setCoords(Vector3d(1.0, 0.5, 0.5));
724 mapping.clear();
725 outData->values().setZero();
726 mapping.computeMapping();
727 mapping.map(inDataID, outDataID);
728 value = outData->values()(0);
729 BOOST_TEST(mapping.hasComputedMapping() == true);
730 BOOST_TEST(value == 4.0);
731
732 vertex.setCoords(Vector3d(1.0, 1.0, 1.0));
733 mapping.clear();
734 outData->values().setZero();
735 mapping.computeMapping();
736 mapping.map(inDataID, outDataID);
737 value = outData->values()(0);
738 BOOST_TEST(mapping.hasComputedMapping() == true);
739 BOOST_TEST(value == 6.0);
740
741 vertex.setCoords(Vector3d(0.5, 0.0, 0.5));
742 mapping.clear();
743 outData->values().setZero();
744 mapping.computeMapping();
745 mapping.map(inDataID, outDataID);
746 value = outData->values()(0);
747 BOOST_TEST(mapping.hasComputedMapping() == true);
748 BOOST_TEST(value == 3.0);
749
750 vertex.setCoords(Vector3d(0.5, 0.5, 1.0));
751 mapping.clear();
752 outData->values().setZero();
753 mapping.computeMapping();
754 mapping.map(inDataID, outDataID);
755 value = outData->values()(0);
756 BOOST_TEST(mapping.hasComputedMapping() == true);
757 BOOST_TEST(value == 4.5);
758
759 vertex.setCoords(Vector3d(0.5, 1.0, 0.0));
760 mapping.clear();
761 outData->values().setZero();
762 mapping.computeMapping();
763 mapping.map(inDataID, outDataID);
764 value = outData->values()(0);
765 BOOST_TEST(mapping.hasComputedMapping() == true);
766 BOOST_TEST(value == 1.5);
767}
768
769void perform3DTestConsistentMappingVector(Mapping &mapping)
770{
771 int dimensions = 3;
772 int dataDimensions = dimensions;
773
774 using Eigen::Vector3d;
775
776 // Create mesh to map from
777 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
778 mesh::PtrData inData = inMesh->createData("InData", dataDimensions, 0_dataID);
779 int inDataID = inData->getID();
780 inMesh->createVertex(Vector3d(0.0, 0.0, 0.0));
781 inMesh->createVertex(Vector3d(1.0, 0.0, 0.0));
782 inMesh->createVertex(Vector3d(1.0, 1.0, 0.0));
783 inMesh->createVertex(Vector3d(0.0, 1.0, 0.0));
784
785 inMesh->createVertex(Vector3d(0.0, 0.0, 1.0));
786 inMesh->createVertex(Vector3d(1.0, 0.0, 1.0));
787 inMesh->createVertex(Vector3d(1.0, 1.0, 1.0));
788 inMesh->createVertex(Vector3d(0.0, 1.0, 1.0));
789
790 inMesh->createVertex(Vector3d(2.0, 0.0, 0.0));
791 inMesh->createVertex(Vector3d(3.0, 0.0, 0.0));
792 inMesh->createVertex(Vector3d(3.0, 1.0, 0.0));
793 inMesh->createVertex(Vector3d(2.0, 1.0, 0.0));
794
795 inMesh->createVertex(Vector3d(2.0, 0.0, 1.0));
796 inMesh->createVertex(Vector3d(3.0, 0.0, 1.0));
797 inMesh->createVertex(Vector3d(3.0, 1.0, 1.0));
798 inMesh->createVertex(Vector3d(2.0, 1.0, 1.0));
799
800 inMesh->createVertex(Vector3d(4.0, 0.0, 0.0));
801 inMesh->createVertex(Vector3d(5.0, 0.0, 0.0));
802 inMesh->createVertex(Vector3d(5.0, 1.0, 0.0));
803 inMesh->createVertex(Vector3d(4.0, 1.0, 0.0));
804
805 inMesh->createVertex(Vector3d(4.0, 0.0, 1.0));
806 inMesh->createVertex(Vector3d(5.0, 0.0, 1.0));
807 inMesh->createVertex(Vector3d(5.0, 1.0, 1.0));
808 inMesh->createVertex(Vector3d(4.0, 1.0, 1.0));
809
810 inMesh->allocateDataValues();
811 addGlobalIndex(inMesh);
812
813 // Create mesh to map to
814 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
815 mesh::PtrData outData = outMesh->createData("OutData", dataDimensions, 1_dataID);
816 int outDataID = outData->getID();
817 mesh::Vertex &vertex = outMesh->createVertex(Vector3d(0, 0, 0));
818 mesh::Vertex &vertex1 = outMesh->createVertex(Vector3d(3.5, 0.5, 0.5));
819 mesh::Vertex &vertex2 = outMesh->createVertex(Vector3d(2.5, 0.5, 1.0));
820 outMesh->createVertex(Vector3d(0, 0, 0.5));
821 outMesh->createVertex(Vector3d(5, 1, 0.5));
822 outMesh->allocateDataValues();
823 addGlobalIndex(outMesh);
824
825 // Setup mapping with mapping coordinates and geometry used
826 mapping.setMeshes(inMesh, outMesh);
827 BOOST_TEST(mapping.hasComputedMapping() == false);
828
829 for (unsigned int i = 0; i < inData->values().size() / dataDimensions; ++i) {
830 inData->values()(i * dataDimensions) = 7;
831 inData->values()(i * dataDimensions + 1) = 38;
832 inData->values()(i * dataDimensions + 2) = 19;
833 }
834
835 // Check for the independent interpolation of each component when applying a constant field
836 mapping.computeMapping();
837 mapping.map(inDataID, outDataID);
838 BOOST_TEST(mapping.hasComputedMapping() == true);
839 BOOST_TEST(outData->values()(0) == 7.0);
840 BOOST_TEST(outData->values()(1) == 38.0);
841 BOOST_TEST(outData->values()(2) == 19.0);
842 BOOST_TEST(outData->values()(9) == 7.0);
843 BOOST_TEST(outData->values()(10) == 38.0);
844 BOOST_TEST(outData->values()(11) == 19.0);
845
846 // Check for the consistency between individual components
847 for (unsigned int i = 0; i < inData->values().size() / dataDimensions; ++i) {
848 inData->values()(i * dataDimensions) = std::pow(i * dataDimensions, 2);
849 inData->values()(i * dataDimensions + 1) = 2 * std::pow(i * dataDimensions, 2);
850 inData->values()(i * dataDimensions + 2) = 3 * std::pow(i * dataDimensions, 2);
851 }
852
853 mapping.clear();
854 outData->values().setZero();
855 mapping.computeMapping();
856 mapping.map(inDataID, outDataID);
857 BOOST_TEST(mapping.hasComputedMapping() == true);
858 BOOST_TEST(outData->values()(12) == 3673.7308684337013);
859 BOOST_TEST(outData->values()(13) == 2 * outData->values()(12));
860 BOOST_TEST(outData->values()(14) == 3 * outData->values()(12));
861
862 // The remaining parts should already be covered by the other 3D/2D tests
863}
864
865void perform3DTestConservativeMapping(Mapping &mapping)
866{
867 int dimensions = 3;
868 using Eigen::Vector3d;
869
870 // Create mesh to map from
871 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
872 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
873 const DataID inDataID = inData->getID();
874 mesh::Vertex &vertex = inMesh->createVertex(Vector3d(0, 0, 0));
875 mesh::Vertex &vertex1 = inMesh->createVertex(Vector3d(3.5, 0.5, 0.5));
876 mesh::Vertex &vertex2 = inMesh->createVertex(Vector3d(2.5, 0.5, 1.0));
877 inMesh->createVertex(Vector3d(1, 1, 1));
878 inMesh->createVertex(Vector3d(5, 1, 0.5));
879 inMesh->allocateDataValues();
880 addGlobalIndex(inMesh);
881
882 // Create mesh to map from
883 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
884 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
885 const DataID outDataID = outData->getID();
886 outMesh->createVertex(Vector3d(0.0, 0.0, 0.0));
887 outMesh->createVertex(Vector3d(1.0, 0.0, 0.0));
888 outMesh->createVertex(Vector3d(1.0, 1.0, 0.0));
889 outMesh->createVertex(Vector3d(0.0, 1.0, 0.0));
890
891 outMesh->createVertex(Vector3d(0.0, 0.0, 1.0));
892 outMesh->createVertex(Vector3d(1.0, 0.0, 1.0));
893 outMesh->createVertex(Vector3d(1.0, 1.0, 1.0));
894 outMesh->createVertex(Vector3d(0.0, 1.0, 1.0));
895
896 outMesh->createVertex(Vector3d(2.0, 0.0, 0.0));
897 outMesh->createVertex(Vector3d(3.0, 0.0, 0.0));
898 outMesh->createVertex(Vector3d(3.0, 1.0, 0.0));
899 outMesh->createVertex(Vector3d(2.0, 1.0, 0.0));
900
901 outMesh->createVertex(Vector3d(2.0, 0.0, 1.0));
902 outMesh->createVertex(Vector3d(3.0, 0.0, 1.0));
903 outMesh->createVertex(Vector3d(3.0, 1.0, 1.0));
904 outMesh->createVertex(Vector3d(2.0, 1.0, 1.0));
905
906 outMesh->createVertex(Vector3d(4.0, 0.0, 0.0));
907 outMesh->createVertex(Vector3d(5.0, 0.0, 0.0));
908 outMesh->createVertex(Vector3d(5.0, 1.0, 0.0));
909 outMesh->createVertex(Vector3d(4.0, 1.0, 0.0));
910
911 outMesh->createVertex(Vector3d(4.0, 0.0, 1.0));
912 outMesh->createVertex(Vector3d(5.0, 0.0, 1.0));
913 outMesh->createVertex(Vector3d(5.0, 1.0, 1.0));
914 outMesh->createVertex(Vector3d(4.0, 1.0, 1.0));
915
916 outMesh->allocateDataValues();
917 addGlobalIndex(outMesh);
918
919 // Setup mapping with mapping coordinates and geometry used
920 mapping.setMeshes(inMesh, outMesh);
921 BOOST_TEST(mapping.hasComputedMapping() == false);
922
923 auto &values = inData->values();
924 // Some arbitrary input values in order to check conservation
925 values << 1.0, 2.0, 2.0, 1.0, 3.0;
926
927 mapping.computeMapping();
928 mapping.map(inDataID, outDataID);
929 BOOST_TEST(mapping.hasComputedMapping() == true);
930 BOOST_TEST(outData->values().sum() == inData->values().sum());
931
932 // Check for the exact reproduction of individual values
933 values << 12.0, 5.0, 7.0, 8.0, 9.0;
934
935 mapping.clear();
936 outData->values().setZero();
937 mapping.computeMapping();
938 mapping.map(inDataID, outDataID);
939 BOOST_TEST(mapping.hasComputedMapping() == true);
940 const double expectedSum = 41;
941 BOOST_TEST(outData->values().sum() == expectedSum);
942 BOOST_TEST(outData->values()(6) == 8);
943 BOOST_TEST(outData->values()(0) == 12);
944
945 // Check for consistency (applying a heavy load to the front layer z = 1)
946 values << 0.0, 50.0, 107.0, 108.0, 48.0;
947
948 mapping.clear();
949 outData->values().setZero();
950 mapping.computeMapping();
951 mapping.map(inDataID, outDataID);
952 BOOST_TEST(mapping.hasComputedMapping() == true);
953 BOOST_TEST(outData->values().sum() == inData->values().sum());
954 BOOST_TEST(outData->values()(0) == 0.0);
955 BOOST_TEST(outData->values()(6) > outData->values()(1));
956 BOOST_TEST(outData->values()(13) > outData->values()(9));
957
958 // Check for symmetry when applying a central load (not guaranteed ?)
959 values << 0.0, 100.0, 0.0, 0.0, 0.0;
960
961 mapping.clear();
962 outData->values().setZero();
963 mapping.computeMapping();
964 mapping.map(inDataID, outDataID);
965 BOOST_TEST(mapping.hasComputedMapping() == true);
966 BOOST_TEST(outData->values().sum() == inData->values().sum());
967 BOOST_TEST(outData->values()(9) == outData->values()(13));
968 BOOST_TEST(outData->values()(9) == outData->values()(10));
969 BOOST_TEST(outData->values()(16) == outData->values()(19));
970 // TODO: There is no symmetry between the (3, x, x) and (4, x , x) layer
971 // as the domain length is 5 and we have a partition center at x = 2.9 and x = 3.75
972 // BOOST_TEST(outData->values()(16) == outData->values()(9));
973 BOOST_TEST(outData->values()(20) == outData->values()(23));
974 BOOST_TEST(outData->values()(16) == outData->values()(20));
975 BOOST_TEST(outData->values()(9) == 15.470584170385226);
976}
977
978void perform3DTestConservativeMappingVector(Mapping &mapping)
979{
980 const int dimensions = 3;
981 const int dataDimension = dimensions;
982 using Eigen::Vector3d;
983
984 // Create mesh to map from
985 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
986 mesh::PtrData inData = inMesh->createData("InData", dataDimension, 0_dataID);
987 const DataID inDataID = inData->getID();
988 mesh::Vertex &vertex = inMesh->createVertex(Vector3d(0, 0, 0));
989 mesh::Vertex &vertex1 = inMesh->createVertex(Vector3d(3.5, 0.5, 0.5));
990 mesh::Vertex &vertex2 = inMesh->createVertex(Vector3d(2.5, 0.5, 1.0));
991 inMesh->createVertex(Vector3d(1, 1, 1));
992 inMesh->createVertex(Vector3d(5, 1, 0.5));
993 inMesh->allocateDataValues();
994 addGlobalIndex(inMesh);
995
996 // Create mesh to map from
997 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
998 mesh::PtrData outData = outMesh->createData("OutData", dataDimension, 1_dataID);
999 const DataID outDataID = outData->getID();
1000 outMesh->createVertex(Vector3d(0.0, 0.0, 0.0));
1001 outMesh->createVertex(Vector3d(1.0, 0.0, 0.0));
1002 outMesh->createVertex(Vector3d(1.0, 1.0, 0.0));
1003 outMesh->createVertex(Vector3d(0.0, 1.0, 0.0));
1004
1005 outMesh->createVertex(Vector3d(0.0, 0.0, 1.0));
1006 outMesh->createVertex(Vector3d(1.0, 0.0, 1.0));
1007 outMesh->createVertex(Vector3d(1.0, 1.0, 1.0));
1008 outMesh->createVertex(Vector3d(0.0, 1.0, 1.0));
1009
1010 outMesh->createVertex(Vector3d(2.0, 0.0, 0.0));
1011 outMesh->createVertex(Vector3d(3.0, 0.0, 0.0));
1012 outMesh->createVertex(Vector3d(3.0, 1.0, 0.0));
1013 outMesh->createVertex(Vector3d(2.0, 1.0, 0.0));
1014
1015 outMesh->createVertex(Vector3d(2.0, 0.0, 1.0));
1016 outMesh->createVertex(Vector3d(3.0, 0.0, 1.0));
1017 outMesh->createVertex(Vector3d(3.0, 1.0, 1.0));
1018 outMesh->createVertex(Vector3d(2.0, 1.0, 1.0));
1019
1020 outMesh->createVertex(Vector3d(4.0, 0.0, 0.0));
1021 outMesh->createVertex(Vector3d(5.0, 0.0, 0.0));
1022 outMesh->createVertex(Vector3d(5.0, 1.0, 0.0));
1023 outMesh->createVertex(Vector3d(4.0, 1.0, 0.0));
1024
1025 outMesh->createVertex(Vector3d(4.0, 0.0, 1.0));
1026 outMesh->createVertex(Vector3d(5.0, 0.0, 1.0));
1027 outMesh->createVertex(Vector3d(5.0, 1.0, 1.0));
1028 outMesh->createVertex(Vector3d(4.0, 1.0, 1.0));
1029
1030 outMesh->allocateDataValues();
1031 addGlobalIndex(outMesh);
1032
1033 // Setup mapping with mapping coordinates and geometry used
1034 mapping.setMeshes(inMesh, outMesh);
1035 BOOST_TEST(mapping.hasComputedMapping() == false);
1036
1037 // Check for the correct sum in the first component
1038 for (unsigned int i = 0; i < inData->values().size() / dataDimension; ++i) {
1039 inData->values()(i * dataDimension) = 7;
1040 inData->values()(i * dataDimension + 1) = 0;
1041 inData->values()(i * dataDimension + 2) = 0;
1042 }
1043
1044 mapping.computeMapping();
1045 mapping.map(inDataID, outDataID);
1046 BOOST_TEST(mapping.hasComputedMapping() == true);
1047 BOOST_TEST(outData->values().sum() == inData->values().sum());
1048 BOOST_TEST(outData->values().sum() == 35);
1049 for (unsigned int i = 0; i < outData->values().size() / dataDimension; ++i) {
1050 BOOST_TEST(outData->values()(i * dataDimension + 1) == 0);
1051 BOOST_TEST(outData->values()(i * dataDimension + 2) == 0);
1052 }
1053
1054 // Check for the correct sum in the second component
1055 for (unsigned int i = 0; i < inData->values().size() / dataDimension; ++i) {
1056 inData->values()(i * dataDimension) = 0;
1057 inData->values()(i * dataDimension + 1) = 27;
1058 inData->values()(i * dataDimension + 2) = 0;
1059 }
1060 mapping.clear();
1061 outData->values().setZero();
1062 mapping.computeMapping();
1063 mapping.map(inDataID, outDataID);
1064 BOOST_TEST(mapping.hasComputedMapping() == true);
1065 BOOST_TEST(outData->values().sum() == inData->values().sum());
1066 BOOST_TEST(outData->values().sum() == 135);
1067
1068 for (unsigned int i = 0; i < outData->values().size() / dataDimension; ++i) {
1069 BOOST_TEST(outData->values()(i * dataDimension) == 0);
1070 BOOST_TEST(outData->values()(i * dataDimension + 2) == 0);
1071 }
1072
1073 // Check for the correct sum in the third component
1074 for (unsigned int i = 0; i < inData->values().size() / dataDimension; ++i) {
1075 inData->values()(i * dataDimension) = 0;
1076 inData->values()(i * dataDimension + 1) = 0;
1077 inData->values()(i * dataDimension + 2) = 3;
1078 }
1079 mapping.clear();
1080 outData->values().setZero();
1081 mapping.computeMapping();
1082 mapping.map(inDataID, outDataID);
1083 BOOST_TEST(mapping.hasComputedMapping() == true);
1084 BOOST_TEST(outData->values().sum() == inData->values().sum());
1085 BOOST_TEST(outData->values().sum() == 15);
1086
1087 for (unsigned int i = 0; i < outData->values().size() / dataDimension; ++i) {
1088 BOOST_TEST(outData->values()(i * dataDimension) == 0);
1089 BOOST_TEST(outData->values()(i * dataDimension + 1) == 0);
1090 }
1091
1092 // Check for the correct relation between copmonents
1093 for (unsigned int i = 0; i < inData->values().size() / dataDimension; ++i) {
1094 inData->values()(i * dataDimension) = std::pow(i * dataDimension, 3);
1095 inData->values()(i * dataDimension + 1) = 5 * std::pow(i * dataDimension, 3);
1096 inData->values()(i * dataDimension + 2) = 10 * std::pow(i * dataDimension, 3);
1097 }
1098 mapping.clear();
1099 outData->values().setZero();
1100 mapping.computeMapping();
1101 mapping.map(inDataID, outDataID);
1102 BOOST_TEST(mapping.hasComputedMapping() == true);
1103 BOOST_TEST(outData->values().sum() == inData->values().sum());
1104
1105 for (unsigned int i = 0; i < outData->values().size() / dataDimension; ++i) {
1106 // Some result values are close to zero and we cannot compare the relation between these values
1107 if (outData->values()(i * dataDimension) > 1e-10) {
1108 BOOST_TEST(outData->values()(i * dataDimension + 1) == 5 * outData->values()(i * dataDimension));
1109 BOOST_TEST(outData->values()(i * dataDimension + 2) == 10 * outData->values()(i * dataDimension));
1110 }
1111 }
1112 BOOST_TEST(outData->values()(55) == 4087.8100404079933);
1113
1114 // The remaining parts should already be covered by the other 3D/2D tests
1115}
1116
1117BOOST_AUTO_TEST_CASE(PartitionOfUnityMappingTests)
1118{
1119 PRECICE_TEST(1_rank);
1120 mapping::CompactPolynomialC0 function(3);
1121 mapping::PartitionOfUnityMapping<CompactPolynomialC0> consistentMap2D(Mapping::CONSISTENT, 2, function, Polynomial::SEPARATE, 5, 0.4, false);
1122 perform2DTestConsistentMapping(consistentMap2D);
1123 mapping::PartitionOfUnityMapping<CompactPolynomialC0> consistentMap2DVector(Mapping::CONSISTENT, 2, function, Polynomial::SEPARATE, 5, 0.4, false);
1124 perform2DTestConsistentMappingVector(consistentMap2DVector);
1125 mapping::PartitionOfUnityMapping<CompactPolynomialC6> consistentMap2DDeadAxis(Mapping::CONSISTENT, 2, mapping::CompactPolynomialC6(6), Polynomial::SEPARATE, 5, 0.4, false);
1126 performTestConsistentMapDeadAxis(consistentMap2DDeadAxis, 2);
1127 mapping::PartitionOfUnityMapping<CompactPolynomialC6> consistentMap3DDeadAxis(Mapping::CONSISTENT, 3, mapping::CompactPolynomialC6(8), Polynomial::SEPARATE, 5, 0.265, false);
1128 performTestConsistentMapDeadAxis(consistentMap3DDeadAxis, 3);
1129 mapping::PartitionOfUnityMapping<CompactPolynomialC0> conservativeMap2D(Mapping::CONSERVATIVE, 2, function, Polynomial::SEPARATE, 5, 0.4, false);
1130 perform2DTestConservativeMapping(conservativeMap2D);
1131 mapping::PartitionOfUnityMapping<CompactPolynomialC0> conservativeMap2DVector(Mapping::CONSERVATIVE, 2, function, Polynomial::SEPARATE, 5, 0.4, false);
1132 perform2DTestConservativeMappingVector(conservativeMap2DVector);
1133 mapping::PartitionOfUnityMapping<CompactPolynomialC0> consistentMap3D(Mapping::CONSISTENT, 3, function, Polynomial::SEPARATE, 5, 0.265, false);
1134 perform3DTestConsistentMapping(consistentMap3D);
1135 mapping::PartitionOfUnityMapping<CompactPolynomialC0> consistentMap3DVector(Mapping::CONSISTENT, 3, function, Polynomial::SEPARATE, 5, 0.265, false);
1136 perform3DTestConsistentMappingVector(consistentMap3DVector);
1137 mapping::PartitionOfUnityMapping<CompactPolynomialC0> conservativeMap3D(Mapping::CONSERVATIVE, 3, function, Polynomial::SEPARATE, 5, 0.265, false);
1138 perform3DTestConservativeMapping(conservativeMap3D);
1139 mapping::PartitionOfUnityMapping<CompactPolynomialC0> conservativeMap3DVector(Mapping::CONSERVATIVE, 3, function, Polynomial::SEPARATE, 5, 0.265, false);
1140 perform3DTestConservativeMappingVector(conservativeMap3DVector);
1141}
1142
1143// Test for small meshes, where the number of requested vertices per cluster is bigger than the global
1144BOOST_AUTO_TEST_CASE(TestSingleClusterPartitionOfUnity)
1145{
1146 PRECICE_TEST(1_rank);
1147 mapping::CompactPolynomialC0 function(3);
1148 mapping::PartitionOfUnityMapping<CompactPolynomialC0> mapping(Mapping::CONSISTENT, 3, function, Polynomial::SEPARATE, 50, 0.4, false);
1149
1150 int dimensions = 3;
1151 using Eigen::Vector3d;
1152
1153 // Create mesh to map from
1154 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
1155 mesh::PtrData inData = inMesh->createData("InData", 1, 0_dataID);
1156 int inDataID = inData->getID();
1157 inMesh->createVertex(Vector3d(0.0, 0.0, 0.0));
1158 inMesh->createVertex(Vector3d(1.0, 0.0, 0.0));
1159 inMesh->createVertex(Vector3d(1.0, 1.0, 0.0));
1160 inMesh->createVertex(Vector3d(0.0, 1.0, 0.0));
1161
1162 inMesh->createVertex(Vector3d(0.0, 0.0, 1.0));
1163 inMesh->createVertex(Vector3d(1.0, 0.0, 1.0));
1164 inMesh->createVertex(Vector3d(1.0, 1.0, 1.0));
1165 inMesh->createVertex(Vector3d(0.0, 1.0, 1.0));
1166
1167 inMesh->createVertex(Vector3d(2.0, 0.0, 0.0));
1168 inMesh->createVertex(Vector3d(3.0, 0.0, 0.0));
1169 inMesh->createVertex(Vector3d(3.0, 1.0, 0.0));
1170 inMesh->createVertex(Vector3d(2.0, 1.0, 0.0));
1171
1172 inMesh->createVertex(Vector3d(2.0, 0.0, 1.0));
1173 inMesh->createVertex(Vector3d(3.0, 0.0, 1.0));
1174 inMesh->createVertex(Vector3d(3.0, 1.0, 1.0));
1175 inMesh->createVertex(Vector3d(2.0, 1.0, 1.0));
1176
1177 inMesh->createVertex(Vector3d(4.0, 0.0, 0.0));
1178 inMesh->createVertex(Vector3d(5.0, 0.0, 0.0));
1179 inMesh->createVertex(Vector3d(5.0, 1.0, 0.0));
1180 inMesh->createVertex(Vector3d(4.0, 1.0, 0.0));
1181
1182 inMesh->createVertex(Vector3d(4.0, 0.0, 1.0));
1183 inMesh->createVertex(Vector3d(5.0, 0.0, 1.0));
1184 inMesh->createVertex(Vector3d(5.0, 1.0, 1.0));
1185 inMesh->createVertex(Vector3d(4.0, 1.0, 1.0));
1186
1187 inMesh->allocateDataValues();
1188 addGlobalIndex(inMesh);
1189
1190 auto &values = inData->values();
1191 // Set the values in the parallel (z) plane 3*values
1192 values << 1.0, 2.0, 2.0, 1.0, 3.0, 6.0, 6.0, 3.0, 3.0, 4.0, 4.0, 3.0, 9.0, 12.0, 12.0, 9.0, 5.0, 6.0, 6.0, 5.0, 15.0, 18.0, 18.0, 15.0;
1193
1194 // Create mesh to map to
1195 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
1196 mesh::PtrData outData = outMesh->createData("OutData", 1, 1_dataID);
1197 int outDataID = outData->getID();
1198 mesh::Vertex &vertex = outMesh->createVertex(Vector3d(0, 0, 0));
1199 mesh::Vertex &vertex1 = outMesh->createVertex(Vector3d(3.5, 0.5, 0.5));
1200 mesh::Vertex &vertex2 = outMesh->createVertex(Vector3d(2.5, 0.5, 1.0));
1201 outMesh->createVertex(Vector3d(0, 0, 0.5));
1202 outMesh->createVertex(Vector3d(5, 1, 0.5));
1203 outMesh->allocateDataValues();
1204 addGlobalIndex(outMesh);
1205
1206 // Setup mapping with mapping coordinates and geometry used
1207 mapping.setMeshes(inMesh, outMesh);
1208 BOOST_TEST(mapping.hasComputedMapping() == false);
1209
1210 vertex.setCoords(Vector3d(0.0, 0.0, 0.0));
1211 mapping.computeMapping();
1212 mapping.map(inDataID, outDataID);
1213 double value = outData->values()(0);
1214 double value1 = outData->values()(1);
1215 double value2 = outData->values()(2);
1216 BOOST_TEST(mapping.hasComputedMapping() == true);
1217 BOOST_TEST(value == 1.0);
1218 BOOST_TEST(value1 == 9.0);
1219 BOOST_TEST(value2 == 10.5);
1220
1221 vertex.setCoords(Vector3d(0.0, 0.5, 0.5));
1222 mapping.clear();
1223 outData->values().setZero();
1224 mapping.computeMapping();
1225 mapping.map(inDataID, outDataID);
1226 value = outData->values()(0);
1227 BOOST_TEST(mapping.hasComputedMapping() == true);
1228 BOOST_TEST(value == 2.0);
1229
1230 vertex.setCoords(Vector3d(1.0, 0.0, 0.0));
1231 mapping.clear();
1232 outData->values().setZero();
1233 mapping.computeMapping();
1234 mapping.map(inDataID, outDataID);
1235 value = outData->values()(0);
1236 BOOST_TEST(mapping.hasComputedMapping() == true);
1237 BOOST_TEST(value == 2.0);
1238
1239 vertex.setCoords(Vector3d(1.0, 0.5, 0.5));
1240 mapping.clear();
1241 outData->values().setZero();
1242 mapping.computeMapping();
1243 mapping.map(inDataID, outDataID);
1244 value = outData->values()(0);
1245 BOOST_TEST(mapping.hasComputedMapping() == true);
1246 BOOST_TEST(value == 4.0);
1247
1248 vertex.setCoords(Vector3d(0.5, 0.5, 1.0));
1249 mapping.clear();
1250 outData->values().setZero();
1251 mapping.computeMapping();
1252 mapping.map(inDataID, outDataID);
1253 value = outData->values()(0);
1254 BOOST_TEST(mapping.hasComputedMapping() == true);
1255 BOOST_TEST(value > 4.6);
1256
1257 vertex.setCoords(Vector3d(0.5, 1.0, 0.0));
1258 mapping.clear();
1259 outData->values().setZero();
1260 mapping.computeMapping();
1261 mapping.map(inDataID, outDataID);
1262 value = outData->values()(0);
1263 BOOST_TEST(mapping.hasComputedMapping() == true);
1264 BOOST_TEST(value < 1.4);
1265}
1266
1267BOOST_AUTO_TEST_SUITE_END() // Serial
1268
1269BOOST_AUTO_TEST_SUITE(Parallel)
1270
1271void addGlobalIndex(mesh::PtrMesh &mesh, int offset = 0)
1272{
1273 for (mesh::Vertex &v : mesh->vertices()) {
1274 v.setGlobalIndex(v.getID() + offset);
1275 }
1276}
1277
1285
1286/*
1287MeshSpecification format:
1288{ {rank, owner rank, {x, y, z}, {v}}, ... }
1289
1290also see struct VertexSpecification.
1291
1292- -1 on rank means all ranks
1293- -1 on owner rank means no rank
1294- x, y, z is position of vertex, z is optional, 2D mesh will be created then
1295- v is the value of the respective vertex. Only 1D supported at this time.
1296
1297ReferenceSpecification format:
1298{ {rank, {v}, ... }
1299- -1 on rank means all ranks
1300- v is the expected value of n-th vertex on that particular rank
1301*/
1302using MeshSpecification = std::vector<VertexSpecification>;
1303
1305using ReferenceSpecification = std::vector<std::pair<int, std::vector<double>>>;
1306
1307void getDistributedMesh(const TestContext & context,
1308 MeshSpecification const &vertices,
1309 mesh::PtrMesh & mesh,
1310 mesh::PtrData & data,
1311 int globalIndexOffset = 0,
1312 bool meshIsSmaller = false)
1313{
1314 Eigen::VectorXd d;
1315
1316 int i = 0;
1317 for (auto &vertex : vertices) {
1318 if (vertex.rank == context.rank or vertex.rank == -1) {
1319 if (vertex.position.size() == 3) // 3-dimensional
1320 mesh->createVertex(Eigen::Vector3d(vertex.position.data()));
1321 else if (vertex.position.size() == 2) // 2-dimensional
1322 mesh->createVertex(Eigen::Vector2d(vertex.position.data()));
1323
1324 int valueDimension = vertex.value.size();
1325
1326 if (vertex.owner == context.rank)
1327 mesh->vertices().back().setOwner(true);
1328 else
1329 mesh->vertices().back().setOwner(false);
1330
1331 d.conservativeResize(i * valueDimension + valueDimension);
1332 // Get data in every dimension
1333 for (int dim = 0; dim < valueDimension; ++dim) {
1334 d(i * valueDimension + dim) = vertex.value.at(dim);
1335 }
1336 i++;
1337 }
1338 }
1339 addGlobalIndex(mesh, globalIndexOffset);
1340 mesh->allocateDataValues();
1341 // All tests use eight vertices
1342 // if (meshIsSmaller) {
1343 // mesh->setGlobalNumberOfVertices(7);
1344 // } else {
1345 // mesh->setGlobalNumberOfVertices(8);
1346 // }
1347 data->values() = d;
1348}
1349
1350void testDistributed(const TestContext & context,
1351 Mapping & mapping,
1352 MeshSpecification inMeshSpec,
1353 MeshSpecification outMeshSpec,
1354 ReferenceSpecification referenceSpec,
1355 int inGlobalIndexOffset = 0,
1356 bool meshIsSmaller = false)
1357{
1358 int meshDimension = inMeshSpec.at(0).position.size();
1359 int valueDimension = inMeshSpec.at(0).value.size();
1360
1361 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", meshDimension, testing::nextMeshID()));
1362 mesh::PtrData inData = inMesh->createData("InData", valueDimension, 0_dataID);
1363 int inDataID = inData->getID();
1364
1365 getDistributedMesh(context, inMeshSpec, inMesh, inData, inGlobalIndexOffset);
1366
1367 mesh::PtrMesh outMesh(new mesh::Mesh("outMesh", meshDimension, testing::nextMeshID()));
1368 mesh::PtrData outData = outMesh->createData("OutData", valueDimension, 1_dataID);
1369 int outDataID = outData->getID();
1370
1371 getDistributedMesh(context, outMeshSpec, outMesh, outData, 0, meshIsSmaller);
1372
1373 mapping.setMeshes(inMesh, outMesh);
1374 BOOST_TEST(mapping.hasComputedMapping() == false);
1375
1376 mapping.computeMapping();
1377 BOOST_TEST(mapping.hasComputedMapping() == true);
1378 mapping.map(inDataID, outDataID);
1379
1380 int index = 0;
1381 for (auto &referenceVertex : referenceSpec) {
1382 if (referenceVertex.first == context.rank or referenceVertex.first == -1) {
1383 for (int dim = 0; dim < valueDimension; ++dim) {
1384 BOOST_TEST_INFO("Index of vertex: " << index << " - Dimension: " << dim);
1385 BOOST_TEST(outData->values()(index * valueDimension + dim) == referenceVertex.second.at(dim));
1386 }
1387 ++index;
1388 }
1389 }
1390 BOOST_TEST(outData->values().size() == index * valueDimension);
1391}
1392
1393BOOST_AUTO_TEST_CASE(DistributedConsistent2D)
1394{
1395 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1396 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
1397
1398 MeshSpecification in{// Consistent mapping: The inMesh is communicated
1399 {-1, 0, {0, 0}, {1}},
1400 {-1, 0, {0, 1}, {2}},
1401 {-1, 1, {1, 0}, {3}},
1402 {-1, 1, {1, 1}, {4}},
1403 {-1, 2, {2, 0}, {5}},
1404 {-1, 2, {2, 1}, {6}},
1405 {-1, 3, {3, 0}, {7}},
1406 {-1, 3, {3, 1}, {8}}};
1407 MeshSpecification out{// The outMesh is local, distributed among all ranks
1408 {0, -1, {0, 0}, {0}},
1409 {0, -1, {0, 1}, {0}},
1410 {1, -1, {1, 0}, {0}},
1411 {1, -1, {1, 1}, {0}},
1412 {2, -1, {2, 0}, {0}},
1413 {2, -1, {2, 1}, {0}},
1414 {3, -1, {3, 0}, {0}},
1415 {3, -1, {3, 1}, {0}}};
1416
1417 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1418 {0, {1}},
1419 {0, {2}},
1420 {1, {3}},
1421 {1, {4}},
1422 {2, {5}},
1423 {2, {6}},
1424 {3, {7}},
1425 {3, {8}}};
1426
1427 mapping::PartitionOfUnityMapping<CompactPolynomialC6> consistentMap2D(Mapping::CONSISTENT, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1428 testDistributed(context, consistentMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1429}
1430
1431// Same as above, but including empty ranks not participating
1432BOOST_AUTO_TEST_CASE(DistributedConsistent2DEmptyOut)
1433{
1434 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1435 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
1436
1437 MeshSpecification in{// Consistent mapping: The inMesh is communicated
1438 {-1, 0, {0, 0}, {1}},
1439 {-1, 0, {0, 1}, {2}},
1440 {-1, 1, {1, 0}, {3}},
1441 {-1, 1, {1, 1}, {4}},
1442 {-1, 2, {2, 0}, {5}},
1443 {-1, 2, {2, 1}, {6}},
1444 {-1, 3, {3, 0}, {7}},
1445 {-1, 3, {3, 1}, {8}}};
1446 MeshSpecification out{// The outMesh is local, distributed among all ranks
1447 {0, 0, {0, 0}, {0}},
1448 {0, 0, {0, 1}, {0}},
1449 {1, 1, {1, 0}, {0}},
1450 {1, 1, {1, 1}, {0}},
1451 {2, 2, {2, 0}, {0}},
1452 {2, 2, {2, 1}, {0}},
1453 {2, 2, {3, 0}, {0}},
1454 {2, 2, {3, 1}, {0}}};
1455
1456 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1457 {0, {1}},
1458 {0, {2}},
1459 {1, {3}},
1460 {1, {4}},
1461 {2, {5}},
1462 {2, {6}},
1463 {2, {7}},
1464 {2, {8}}};
1465
1466 mapping::PartitionOfUnityMapping<CompactPolynomialC6> consistentMap2D(Mapping::CONSISTENT, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1467 testDistributed(context, consistentMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1468}
1469
1470// Same as above, but including empty ranks not participating
1471BOOST_AUTO_TEST_CASE(DistributedConsistent2DEmptyRank)
1472{
1473 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1474 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
1475
1476 MeshSpecification in{// Consistent mapping: The inMesh is communicated
1477 {-1, 0, {0, 0}, {1}},
1478 {-1, 0, {0, 1}, {2}},
1479 {-1, 1, {1, 0}, {3}},
1480 {-1, 1, {1, 1}, {4}},
1481 {-1, 2, {2, 0}, {5}},
1482 {-1, 2, {2, 1}, {6}}};
1483
1484 MeshSpecification out{// The outMesh is local, distributed among all ranks
1485 {0, 0, {0, 0}, {0}},
1486 {0, 0, {0, 1}, {0}},
1487 {1, 1, {1, 0}, {0}},
1488 {1, 1, {1, 1}, {0}},
1489 {2, 2, {2, 0}, {0}},
1490 {2, 2, {2, 1}, {0}}};
1491
1492 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1493 {0, {1}},
1494 {0, {2}},
1495 {1, {3}},
1496 {1, {4}},
1497 {2, {5}},
1498 {2, {6}}};
1499
1500 mapping::PartitionOfUnityMapping<CompactPolynomialC6> consistentMap2D(Mapping::CONSISTENT, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1501 testDistributed(context, consistentMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1502}
1503
1504BOOST_AUTO_TEST_CASE(DistributedConservative2D)
1505{
1506 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1507 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
1508
1509 MeshSpecification in{// Conservative mapping: The inMesh is local
1510 {0, -1, {0, 0}, {1}},
1511 {0, -1, {0, 1}, {2}},
1512 {1, -1, {1, 0}, {3}},
1513 {1, -1, {1, 1}, {4}},
1514 {2, -1, {2, 0}, {5}},
1515 {2, -1, {2, 1}, {6}},
1516 {3, -1, {3, 0}, {7}},
1517 {3, -1, {3, 1}, {8}}};
1518 MeshSpecification out{// The outMesh is remote, distributed among all ranks
1519 {0, 0, {0, 0}, {0}},
1520 {0, 0, {0, 1}, {0}},
1521 {1, 1, {1, 0}, {0}},
1522 {1, 1, {1, 1}, {0}},
1523 {2, 2, {2, 0}, {0}},
1524 {2, 2, {2, 1}, {0}},
1525 {3, 3, {3, 0}, {0}},
1526 {3, 3, {3, 1}, {0}}};
1527
1528 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1529 {0, {1}},
1530 {0, {2}},
1531 {1, {3}},
1532 {1, {4}},
1533 {2, {5}},
1534 {2, {6}},
1535 {3, {7}},
1536 {3, {8}}};
1537
1538 mapping::PartitionOfUnityMapping<CompactPolynomialC6> conservativeMap2D(Mapping::CONSERVATIVE, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1539 testDistributed(context, conservativeMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1540}
1541
1542// Same as above, but checking with empty ranks
1543BOOST_AUTO_TEST_CASE(DistributedConservative2DEmptyRank)
1544{
1545 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1546 std::vector<int> globalIndexOffsets = {0, 2, 4, 6};
1547
1548 MeshSpecification in{// Conservative mapping: The inMesh is local
1549 {0, -1, {0, 0}, {1}},
1550 {0, -1, {0, 1}, {2}},
1551 {1, -1, {1, 0}, {3}},
1552 {1, -1, {1, 1}, {4}},
1553 {2, -1, {2, 0}, {5}},
1554 {2, -1, {2, 1}, {6}},
1555 {2, -1, {3, 0}, {7}},
1556 {2, -1, {3, 1}, {8}}};
1557 MeshSpecification out{// The outMesh is remote, distributed among all ranks
1558 {-1, 0, {0, 0}, {0}},
1559 {-1, 0, {0, 1}, {0}},
1560 {-1, 1, {1, 0}, {0}},
1561 {-1, 1, {1, 1}, {0}},
1562 {-1, 2, {2, 0}, {0}},
1563 {-1, 2, {2, 1}, {0}},
1564 {-1, 3, {3, 0}, {0}},
1565 {-1, 3, {3, 1}, {0}}};
1566
1567 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1568 {0, {1}},
1569 {0, {2}},
1570 {0, {0}},
1571 {0, {0}},
1572 {0, {0}},
1573 {0, {0}},
1574 {0, {0}},
1575 {0, {0}},
1576 {1, {0}},
1577 {1, {0}},
1578 {1, {3}},
1579 {1, {4}},
1580 {1, {0}},
1581 {1, {0}},
1582 {1, {0}},
1583 {1, {0}},
1584 {2, {0}},
1585 {2, {0}},
1586 {2, {0}},
1587 {2, {0}},
1588 {2, {5}},
1589 {2, {6}},
1590 {2, {7}},
1591 {2, {8}},
1592 {3, {0}},
1593 {3, {0}},
1594 {3, {0}},
1595 {3, {0}},
1596 {3, {0}},
1597 {3, {0}},
1598 {3, {0}},
1599 {3, {0}}};
1600
1601 mapping::PartitionOfUnityMapping<CompactPolynomialC6> conservativeMap2D(Mapping::CONSERVATIVE, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1602 testDistributed(context, conservativeMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1603}
1604
1605// Same as above, but checking the primary rank for all output vertices
1606BOOST_AUTO_TEST_CASE(DistributedConservative2DTwoRanks)
1607{
1608 PRECICE_TEST(""_on(2_ranks).setupIntraComm());
1609 std::vector<int> globalIndexOffsets = {0, 2};
1610
1611 MeshSpecification in{// Conservative mapping: The inMesh is local
1612 {0, -1, {0, 0}, {1}},
1613 {0, -1, {0, 1}, {2}},
1614 {0, -1, {1, 0}, {3}},
1615 {0, -1, {1, 1}, {4}},
1616 {1, -1, {2, 0}, {5}},
1617 {1, -1, {2, 1}, {6}},
1618 {1, -1, {3, 0}, {7}},
1619 {1, -1, {3, 1}, {8}}};
1620 MeshSpecification out{// The outMesh is remote, distributed among all ranks
1621 {-1, 0, {0, 0}, {0}},
1622 {-1, 0, {0, 1}, {0}},
1623 {-1, 0, {1, 0}, {0}},
1624 {-1, 0, {1, 1}, {0}},
1625 {-1, 0, {2, 0}, {0}},
1626 {-1, 0, {2, 1}, {0}},
1627 {-1, 0, {3, 0}, {0}},
1628 {-1, 0, {3, 1}, {0}}};
1629
1630 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
1631 {0, {1}},
1632 {0, {2}},
1633 {0, {3}},
1634 {0, {4}},
1635 {0, {0}},
1636 {0, {0}},
1637 {0, {0}},
1638 {0, {0}},
1639 {1, {0}},
1640 {1, {0}},
1641 {1, {0}},
1642 {1, {0}},
1643 {1, {5}},
1644 {1, {6}},
1645 {1, {7}},
1646 {1, {8}}};
1647
1648 mapping::PartitionOfUnityMapping<CompactPolynomialC6> conservativeMap2D(Mapping::CONSERVATIVE, 2, CompactPolynomialC6(3.), Polynomial::SEPARATE, 5, 0.3, false);
1649 testDistributed(context, conservativeMap2D, in, out, ref, globalIndexOffsets.at(context.rank));
1650}
1651
1652void testTagging(const TestContext &context,
1653 MeshSpecification inMeshSpec,
1654 MeshSpecification outMeshSpec,
1655 MeshSpecification shouldTagFirstRound,
1656 MeshSpecification shouldTagSecondRound,
1657 bool consistent)
1658{
1659 int meshDimension = inMeshSpec.at(0).position.size();
1660 int valueDimension = inMeshSpec.at(0).value.size();
1661
1662 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", meshDimension, testing::nextMeshID()));
1663 mesh::PtrData inData = inMesh->createData("InData", valueDimension, 0_dataID);
1664 getDistributedMesh(context, inMeshSpec, inMesh, inData);
1665
1666 mesh::PtrMesh outMesh(new mesh::Mesh("outMesh", meshDimension, testing::nextMeshID()));
1667 mesh::PtrData outData = outMesh->createData("OutData", valueDimension, 1_dataID);
1668 getDistributedMesh(context, outMeshSpec, outMesh, outData);
1669
1670 Mapping::Constraint constr = consistent ? Mapping::CONSISTENT : Mapping::CONSERVATIVE;
1671 mapping::PartitionOfUnityMapping<CompactPolynomialC4> mapping(constr, 2, CompactPolynomialC4(2), Polynomial::SEPARATE, 2, 0.3, false);
1672 inMesh->computeBoundingBox();
1673 outMesh->computeBoundingBox();
1674
1675 mapping.setMeshes(inMesh, outMesh);
1676 mapping.tagMeshFirstRound();
1677
1678 for (const auto &v : inMesh->vertices()) {
1679 auto pos = std::find_if(shouldTagFirstRound.begin(), shouldTagFirstRound.end(),
1680 [meshDimension, &v](const VertexSpecification &spec) {
1681 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1682 });
1683 bool found = pos != shouldTagFirstRound.end();
1684 BOOST_TEST(found >= v.isTagged(),
1685 "FirstRound: Vertex " << v << " is tagged, but should not be.");
1686 BOOST_TEST(found <= v.isTagged(),
1687 "FirstRound: Vertex " << v << " is not tagged, but should be.");
1688 }
1689
1690 mapping.tagMeshSecondRound();
1691
1692 for (const auto &v : inMesh->vertices()) {
1693 auto posFirst = std::find_if(shouldTagFirstRound.begin(), shouldTagFirstRound.end(),
1694 [meshDimension, &v](const VertexSpecification &spec) {
1695 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1696 });
1697 bool foundFirst = posFirst != shouldTagFirstRound.end();
1698 auto posSecond = std::find_if(shouldTagSecondRound.begin(), shouldTagSecondRound.end(),
1699 [meshDimension, &v](const VertexSpecification &spec) {
1700 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1701 });
1702 bool foundSecond = posSecond != shouldTagSecondRound.end();
1703 BOOST_TEST(foundFirst <= v.isTagged(), "SecondRound: Vertex " << v
1704 << " is not tagged, but should be from the first round.");
1705 BOOST_TEST(foundSecond <= v.isTagged(), "SecondRound: Vertex " << v
1706 << " is not tagged, but should be.");
1707 BOOST_TEST((foundSecond or foundFirst) >= v.isTagged(), "SecondRound: Vertex " << v
1708 << " is tagged, but should not be.");
1709 }
1710}
1711
1712BOOST_AUTO_TEST_CASE(testTagFirstRound)
1713{
1714 PRECICE_TEST(""_on(4_ranks).setupIntraComm())
1715
1716 MeshSpecification outMeshSpec = {
1717 {0, -1, {0, 0}, {0}}};
1718 MeshSpecification inMeshSpec = {
1719 {0, -1, {-1, 0}, {1}}, //inside
1720 {0, -1, {-3, 0}, {1}}, //outside
1721 {0, 0, {1, 0}, {1}}, //inside, owner
1722 {0, -1, {3, 0}, {1}}, //outside
1723 {0, -1, {0, -1}, {1}}, //inside
1724 {0, -1, {0, -3}, {1}}, //outside
1725 {0, -1, {0, 1}, {1}}, //inside
1726 {0, -1, {0, 3}, {1}} //outside
1727 };
1728 MeshSpecification shouldTagFirstRound = {
1729 {0, -1, {-1, 0}, {1}},
1730 {0, -1, {1, 0}, {1}},
1731 {0, -1, {0, -1}, {1}},
1732 {0, -1, {0, 1}, {1}}};
1733 // No tagging happening in round two
1734 MeshSpecification shouldTagSecondRound;
1735
1736 testTagging(context, inMeshSpec, outMeshSpec, shouldTagFirstRound, shouldTagSecondRound, true);
1737 // For conservative just swap meshes.
1738 testTagging(context, outMeshSpec, inMeshSpec, shouldTagFirstRound, shouldTagSecondRound, false);
1739}
1740
1741BOOST_AUTO_TEST_SUITE_END() // Parallel
1742
1743BOOST_AUTO_TEST_SUITE_END()
1744BOOST_AUTO_TEST_SUITE_END()
out
std::ostream & out
Definition EventUtils.cpp:193
BOOST_AUTO_TEST_SUITE
BOOST_AUTO_TEST_SUITE(PreProcess)
BOOST_AUTO_TEST_SUITE_END
BOOST_AUTO_TEST_SUITE_END()
testTagging
void testTagging(const TestContext &context, MeshSpecification inMeshSpec, MeshSpecification outMeshSpec, MeshSpecification shouldTagFirstRound, MeshSpecification shouldTagSecondRound, bool consistent)
Definition PartitionOfUnityMappingTest.cpp:1652
addGlobalIndex
void addGlobalIndex(mesh::PtrMesh &mesh, int offset=0)
Definition PartitionOfUnityMappingTest.cpp:29
perform2DTestConservativeMappingVector
void perform2DTestConservativeMappingVector(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:476
perform3DTestConsistentMappingVector
void perform3DTestConsistentMappingVector(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:769
perform2DTestConsistentMappingVector
void perform2DTestConsistentMappingVector(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:180
perform2DTestConsistentMapping
void perform2DTestConsistentMapping(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:48
BOOST_AUTO_TEST_CASE
BOOST_AUTO_TEST_CASE(PartitionOfUnityMappingTests)
Definition PartitionOfUnityMappingTest.cpp:1117
perform3DTestConsistentMapping
void perform3DTestConsistentMapping(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:623
testDistributed
void testDistributed(const TestContext &context, Mapping &mapping, MeshSpecification inMeshSpec, MeshSpecification outMeshSpec, ReferenceSpecification referenceSpec, int inGlobalIndexOffset=0, bool meshIsSmaller=false)
Definition PartitionOfUnityMappingTest.cpp:1350
sumComponentWise
double sumComponentWise(const Eigen::VectorXd &vec, int component, int dataDimension)
Definition PartitionOfUnityMappingTest.cpp:36
perform3DTestConservativeMapping
void perform3DTestConservativeMapping(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:865
perform2DTestConservativeMapping
void perform2DTestConservativeMapping(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:379
perform3DTestConservativeMappingVector
void perform3DTestConservativeMappingVector(Mapping &mapping)
Definition PartitionOfUnityMappingTest.cpp:978
getDistributedMesh
void getDistributedMesh(const TestContext &context, MeshSpecification const &vertices, mesh::PtrMesh &mesh, mesh::PtrData &data, int globalIndexOffset=0, bool meshIsSmaller=false)
Definition PartitionOfUnityMappingTest.cpp:1307
performTestConsistentMapDeadAxis
void performTestConsistentMapDeadAxis(Mapping &mapping, int dim)
Definition PartitionOfUnityMappingTest.cpp:317
index
unsigned int index
Definition SocketCommunication.cpp:717
PRECICE_TEST
#define PRECICE_TEST(...)
Definition Testing.hpp:27
PRECICE_ASSERT
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
std::vector::at
T at(T... args)
std::deque::back
T back(T... args)
precice::mapping::CompactPolynomialC0
Wendland radial basis function with compact support.
Definition BasisFunctions.hpp:368
precice::mapping::CompactPolynomialC4
Wendland radial basis function with compact support.
Definition BasisFunctions.hpp:475
precice::mapping::CompactPolynomialC6
Wendland radial basis function with compact support.
Definition BasisFunctions.hpp:529
precice::mapping::Mapping
Abstract base class for mapping of data from one mesh to another.
Definition Mapping.hpp:15
precice::mapping::Mapping::Constraint
Constraint
Specifies additional constraints for a mapping.
Definition Mapping.hpp:29
precice::mapping::Mapping::clear
virtual void clear()=0
Removes a computed mapping.
precice::mapping::Mapping::setMeshes
void setMeshes(const mesh::PtrMesh &input, const mesh::PtrMesh &output)
Sets input and output meshes carrying data to be mapped.
Definition Mapping.cpp:27
precice::mapping::Mapping::computeMapping
virtual void computeMapping()=0
Computes the mapping coefficients from the in- and output mesh.
precice::mapping::Mapping::hasComputedMapping
bool hasComputedMapping() const
Returns true, if the mapping has been computed.
Definition Mapping.cpp:252
precice::mapping::Mapping::map
void map(int inputDataID, int outputDataID)
Definition Mapping.cpp:126
precice::mapping::PartitionOfUnityMapping::tagMeshSecondRound
void tagMeshSecondRound() final override
nothing to do here
Definition PartitionOfUnityMapping.hpp:360
precice::mapping::PartitionOfUnityMapping::tagMeshFirstRound
void tagMeshFirstRound() final override
tag the vertices required for the mapping
Definition PartitionOfUnityMapping.hpp:294
precice::mapping::PartitionOfUnityMapping::clear
void clear() final override
Clears a computed mapping by deleting the content of the _clusters vector.
Definition PartitionOfUnityMapping.hpp:415
precice::mesh::Data::getID
DataID getID() const
Returns the ID of the data set (supposed to be unique).
Definition Data.cpp:93
precice::mesh::Data::values
Eigen::VectorXd & values()
Returns a reference to the data values.
Definition Data.cpp:28
precice::mesh::Mesh
Container and creator for meshes.
Definition Mesh.hpp:39
precice::mesh::Mesh::setGlobalNumberOfVertices
void setGlobalNumberOfVertices(int num)
Definition Mesh.hpp:264
precice::mesh::Mesh::vertices
VertexContainer & vertices()
Returns modifieable container holding all vertices.
Definition Mesh.cpp:53
precice::mesh::Mesh::nVertices
std::size_t nVertices() const
Returns the number of vertices.
Definition Mesh.cpp:63
precice::mesh::Mesh::createData
PtrData & createData(const std::string &name, int dimension, DataID id, int waveformDegree=time::Time::DEFAULT_WAVEFORM_DEGREE)
Create only data for vertex.
Definition Mesh.cpp:151
precice::mesh::Mesh::computeBoundingBox
void computeBoundingBox()
Computes the boundingBox for the vertices.
Definition Mesh.cpp:242
precice::mesh::Mesh::createVertex
Vertex & createVertex(const Eigen::VectorXd &coords)
Creates and initializes a Vertex object.
Definition Mesh.cpp:103
precice::mesh::Mesh::allocateDataValues
void allocateDataValues()
Allocates memory for the vertex data values and corresponding gradient values.
Definition Mesh.cpp:233
precice::mesh::Vertex
Vertex of a mesh.
Definition Vertex.hpp:16
precice::mesh::Vertex::setCoords
void setCoords(const VECTOR_T &coordinates)
Sets the coordinates of the vertex.
Definition Vertex.hpp:102
precice::mesh::Vertex::setOwner
void setOwner(bool owner)
Definition Vertex.cpp:27
precice::testing::TestContext::rank
Rank rank
the rank of the current participant
Definition TestContext.hpp:160
std::find_if
T find_if(T... args)
precice::mapping
contains data mapping from points to meshes.
Definition AxialGeoMultiscaleMapping.cpp:5
precice::math::equals
constexpr bool equals(const Eigen::MatrixBase< DerivedA > &A, const Eigen::MatrixBase< DerivedB > &B, double tolerance=NUMERICAL_ZERO_DIFFERENCE)
Compares two Eigen::MatrixBase for equality up to tolerance.
Definition differences.hpp:13
precice::mesh
provides Mesh, Data and primitives.
Definition ScaleActionTest.cpp:20
precice::testing
contains the testing framework.
Definition helper.hpp:9
precice
Main namespace of the precice library.
Definition Acceleration.cpp:5
precice::DataID
int DataID
Definition Types.hpp:25
std::pow
T pow(T... args)
VertexSpecification
Holds rank, owner, position and value of a single vertex.
Definition PartitionOfUnityMappingTest.cpp:1279