preCICE v3.1.2
Loading...
Searching...
No Matches
RadialBasisFctMappingTest.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"
15#include "mesh/Data.hpp"
16#include "mesh/Mesh.hpp"
18#include "mesh/Utils.hpp"
19#include "mesh/Vertex.hpp"
21#include "testing/Testing.hpp"
22
23using namespace precice;
24using namespace precice::mesh;
25using namespace precice::mapping;
26using namespace precice::testing;
28
29BOOST_AUTO_TEST_SUITE(MappingTests)
30BOOST_AUTO_TEST_SUITE(RadialBasisFunctionMapping)
31
33
34
36 int rank;
37 int owner;
40};
41
42/*
43MeshSpecification format:
44{ {rank, owner rank, {x, y, z}, {v}}, ... }
45
46also see struct VertexSpecification.
47
48- -1 on rank means all ranks
49- -1 on owner rank means no rank
50- x, y, z is position of vertex, z is optional, 2D mesh will be created then
51- v is the value of the respective vertex. Only 1D supported at this time.
52
53ReferenceSpecification format:
54{ {rank, {v}, ... }
55- -1 on rank means all ranks
56- v is the expected value of n-th vertex on that particular rank
57*/
63
66
68 MeshSpecification &meshSpec,
69 int globalIndexOffset = 0,
70 bool meshIsSmaller = false)
71{
72 mesh::PtrMesh distributedMesh(new mesh::Mesh(meshSpec.meshName, meshSpec.meshDimension, testing::nextMeshID()));
73 int i = 0;
74 for (auto &vertex : meshSpec.vertices) {
75 if (vertex.rank == context.rank or vertex.rank == -1) {
76 if (vertex.position.size() == 3) // 3-dimensional
77 distributedMesh->createVertex(Eigen::Vector3d(vertex.position.data()));
78 else if (vertex.position.size() == 2) // 2-dimensional
79 distributedMesh->createVertex(Eigen::Vector2d(vertex.position.data()));
80
81 if (vertex.owner == context.rank)
82 distributedMesh->vertices().back().setOwner(true);
83 else
84 distributedMesh->vertices().back().setOwner(false);
85
86 i++;
87 }
88 }
89 addGlobalIndex(distributedMesh, globalIndexOffset);
90 // All tests use eight vertices
91 if (meshIsSmaller) {
92 distributedMesh->setGlobalNumberOfVertices(7);
93 } else {
94 distributedMesh->setGlobalNumberOfVertices(8);
95 }
96 return distributedMesh;
97}
98
99Eigen::VectorXd getDistributedData(const TestContext & context,
100 MeshSpecification const &meshSpec)
101{
102 Eigen::VectorXd d;
103
104 int i = 0;
105 for (auto &vertex : meshSpec.vertices) {
106 if (vertex.rank == context.rank or vertex.rank == -1) {
107 int valueDimension = vertex.value.size();
108 d.conservativeResize(i * valueDimension + valueDimension);
109 // Get data in every dimension
110 for (int dim = 0; dim < valueDimension; ++dim) {
111 d(i * valueDimension + dim) = vertex.value.at(dim);
112 }
113 i++;
114 }
115 }
116 return d;
117}
118
119void testDistributed(const TestContext & context,
120 Mapping & mapping,
121 MeshSpecification inMeshSpec,
122 MeshSpecification outMeshSpec,
123 ReferenceSpecification referenceSpec,
124 int inGlobalIndexOffset = 0,
125 bool meshIsSmaller = false)
126{
127 int valueDimension = inMeshSpec.vertices.at(0).value.size();
128
129 mesh::PtrMesh inMesh = getDistributedMesh(context, inMeshSpec, inGlobalIndexOffset);
130 Eigen::VectorXd inValues = getDistributedData(context, inMeshSpec);
131
132 mesh::PtrMesh outMesh = getDistributedMesh(context, outMeshSpec, 0, meshIsSmaller);
133 Eigen::VectorXd outValues = getDistributedData(context, outMeshSpec);
134 mapping.setMeshes(inMesh, outMesh);
135 BOOST_TEST(mapping.hasComputedMapping() == false);
136
137 mapping.computeMapping();
138 BOOST_TEST(mapping.hasComputedMapping() == true);
139 time::Sample inSample(valueDimension, inValues);
140 mapping.map(inSample, outValues);
141
142 int index = 0;
143 for (auto &referenceVertex : referenceSpec) {
144 if (referenceVertex.first == context.rank or referenceVertex.first == -1) {
145 for (int dim = 0; dim < valueDimension; ++dim) {
146 BOOST_TEST_INFO("Index of vertex: " << index << " - Dimension: " << dim);
147 BOOST_TEST(outValues(index * valueDimension + dim) == referenceVertex.second.at(dim));
148 }
149 ++index;
150 }
151 }
152 BOOST_TEST(outValues.size() == index * valueDimension);
153}
154
155constexpr int meshDims2D{2};
156
158BOOST_AUTO_TEST_CASE(DistributedConsistent2DV1)
159{
160 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
161 Multiquadrics fct(5.0);
162
164 {-1, 0, {0, 0}, {1}},
165 {-1, 0, {0, 1}, {2}},
166 {-1, 1, {1, 0}, {3}},
167 {-1, 1, {1, 1}, {4}},
168 {-1, 2, {2, 0}, {5}},
169 {-1, 2, {2, 1}, {6}},
170 {-1, 3, {3, 0}, {7}},
171 {-1, 3, {3, 1}, {8}}};
172
174 std::move(inVertexList),
176 "inMesh"};
177
179 {0, -1, {0, 0}, {0}},
180 {0, -1, {0, 1}, {0}},
181 {1, -1, {1, 0}, {0}},
182 {1, -1, {1, 1}, {0}},
183 {2, -1, {2, 0}, {0}},
184 {2, -1, {2, 1}, {0}},
185 {3, -1, {3, 0}, {0}},
186 {3, -1, {3, 1}, {0}}};
187
189 std::move(outVertexList),
191 "outMesh"};
192 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
193 {0, {1}},
194 {0, {2}},
195 {1, {3}},
196 {1, {4}},
197 {2, {5}},
198 {2, {6}},
199 {3, {7}},
200 {3, {8}}};
201
202 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
203 testDistributed(context, mapping_on, in, out, ref);
204 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
205 testDistributed(context, mapping_sep, in, out, ref);
206 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
207 testDistributed(context, mapping_off, in, out, ref);
208}
209
210BOOST_AUTO_TEST_CASE(DistributedConsistent2DV1Vector)
211{
212 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
213 Multiquadrics fct(5.0);
214
215 std::vector<VertexSpecification> inVertexList{// Consistent mapping: The inMesh is communicated
216 {-1, 0, {0, 0}, {1, 4}},
217 {-1, 0, {0, 1}, {2, 5}},
218 {-1, 1, {1, 0}, {3, 6}},
219 {-1, 1, {1, 1}, {4, 7}},
220 {-1, 2, {2, 0}, {5, 8}},
221 {-1, 2, {2, 1}, {6, 9}},
222 {-1, 3, {3, 0}, {7, 10}},
223 {-1, 3, {3, 1}, {8, 11}}};
224 MeshSpecification in{// The outMesh is local, distributed among all ranks
225 std::move(inVertexList),
227 "inMesh"};
228
229 std::vector<VertexSpecification> outVertexList{// The outMesh is local, distributed among all ranks
230 {0, -1, {0, 0}, {0, 0}},
231 {0, -1, {0, 1}, {0, 0}},
232 {1, -1, {1, 0}, {0, 0}},
233 {1, -1, {1, 1}, {0, 0}},
234 {2, -1, {2, 0}, {0, 0}},
235 {2, -1, {2, 1}, {0, 0}},
236 {3, -1, {3, 0}, {0, 0}},
237 {3, -1, {3, 1}, {0, 0}}};
239 std::move(outVertexList),
241 "outMesh"};
242 ReferenceSpecification ref{// Tests for {0, 1} on the first rank, {1, 2} on the second, ...
243 {0, {1, 4}},
244 {0, {2, 5}},
245 {1, {3, 6}},
246 {1, {4, 7}},
247 {2, {5, 8}},
248 {2, {6, 9}},
249 {3, {7, 10}},
250 {3, {8, 11}}};
251
252 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
253 testDistributed(context, mapping_on, in, out, ref);
254 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
255 testDistributed(context, mapping_sep, in, out, ref);
256 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
257 testDistributed(context, mapping_off, in, out, ref);
258}
259
261BOOST_AUTO_TEST_CASE(DistributedConsistent2DV2)
262{
263 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
264 Multiquadrics fct(5.0);
265
266 std::vector<VertexSpecification> inVertexList{// Consistent mapping: The inMesh is communicated, rank 2 owns no vertices
267 {-1, 0, {0, 0}, {1}},
268 {-1, 0, {0, 1}, {2}},
269 {-1, 1, {1, 0}, {3}},
270 {-1, 1, {1, 1}, {4}},
271 {-1, 1, {2, 0}, {5}},
272 {-1, 3, {2, 1}, {6}},
273 {-1, 3, {3, 0}, {7}},
274 {-1, 3, {3, 1}, {8}}};
275
276 MeshSpecification in = {
277 std::move(inVertexList),
279 "inMesh"};
280 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 1 is empty
281 {0, -1, {0, 0}, {0}},
282 {0, -1, {0, 1}, {0}},
283 {0, -1, {1, 0}, {0}},
284 {2, -1, {1, 1}, {0}},
285 {2, -1, {2, 0}, {0}},
286 {2, -1, {2, 1}, {0}},
287 {3, -1, {3, 0}, {0}},
288 {3, -1, {3, 1}, {0}}};
290 std::move(outVertexList),
292 "outMesh"};
293 ReferenceSpecification ref{// Tests for {0, 1, 2} on the first rank,
294 // second rank (consistent with the outMesh) is empty, ...
295 {0, {1}},
296 {0, {2}},
297 {0, {3}},
298 {2, {4}},
299 {2, {5}},
300 {2, {6}},
301 {3, {7}},
302 {3, {8}}};
303 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
304 testDistributed(context, mapping_on, in, out, ref);
305 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
306 testDistributed(context, mapping_sep, in, out, ref);
307 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
308 testDistributed(context, mapping_off, in, out, ref);
309}
310
312BOOST_AUTO_TEST_CASE(DistributedConsistent2DV3)
313{
314 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
315 Multiquadrics fct(5.0);
316
317 std::vector<int> globalIndexOffsets = {0, 0, 0, 4};
318
320 // Rank 0 has part of the mesh, owns a subpart
321 {0, 0, {0, 0}, {1}},
322 {0, 0, {0, 1}, {2}},
323 {0, 0, {1, 0}, {3}},
324 {0, -1, {1, 1}, {4}},
325 {0, -1, {2, 0}, {5}},
326 {0, -1, {2, 1}, {6}},
327 // Rank 1 has no vertices
328 // Rank 2 has the entire mesh, but owns just 3 and 5.
329 {2, -1, {0, 0}, {1}},
330 {2, -1, {0, 1}, {2}},
331 {2, -1, {1, 0}, {3}},
332 {2, 2, {1, 1}, {4}},
333 {2, -1, {2, 0}, {5}},
334 {2, 2, {2, 1}, {6}},
335 {2, -1, {3, 0}, {7}},
336 {2, -1, {3, 1}, {8}},
337 // Rank 3 has the last 4 vertices, owns 4, 6 and 7
338 {3, 3, {2, 0}, {5}},
339 {3, -1, {2, 1}, {6}},
340 {3, 3, {3, 0}, {7}},
341 {3, 3, {3, 1}, {8}},
342 };
344 std::move(inVertexList),
346 "inMesh"};
347 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 1 is empty
348 {0, -1, {0, 0}, {0}},
349 {0, -1, {0, 1}, {0}},
350 {0, -1, {1, 0}, {0}},
351 {2, -1, {1, 1}, {0}},
352 {2, -1, {2, 0}, {0}},
353 {2, -1, {2, 1}, {0}},
354 {3, -1, {3, 0}, {0}},
355 {3, -1, {3, 1}, {0}}};
357 std::move(outVertexList),
359 "outMesh"};
360 ReferenceSpecification ref{// Tests for {0, 1, 2} on the first rank,
361 // second rank (consistent with the outMesh) is empty, ...
362 {0, {1}},
363 {0, {2}},
364 {0, {3}},
365 {2, {4}},
366 {2, {5}},
367 {2, {6}},
368 {3, {7}},
369 {3, {8}}};
370
371 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
372 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
373 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
374 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
375 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
376 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
377}
378
380BOOST_AUTO_TEST_CASE(DistributedConsistent2DV3Vector)
381{
382 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
383 Multiquadrics fct(5.0);
384
385 std::vector<int> globalIndexOffsets = {0, 0, 0, 4};
386
388 // Rank 0 has part of the mesh, owns a subpart
389 {0, 0, {0, 0}, {1, 4}},
390 {0, 0, {0, 1}, {2, 5}},
391 {0, 0, {1, 0}, {3, 6}},
392 {0, -1, {1, 1}, {4, 7}},
393 {0, -1, {2, 0}, {5, 8}},
394 {0, -1, {2, 1}, {6, 9}},
395 // Rank 1 has no vertices
396 // Rank 2 has the entire mesh, but owns just 3 and 5.
397 {2, -1, {0, 0}, {1, 4}},
398 {2, -1, {0, 1}, {2, 5}},
399 {2, -1, {1, 0}, {3, 6}},
400 {2, 2, {1, 1}, {4, 7}},
401 {2, -1, {2, 0}, {5, 8}},
402 {2, 2, {2, 1}, {6, 9}},
403 {2, -1, {3, 0}, {7, 10}},
404 {2, -1, {3, 1}, {8, 11}},
405 // Rank 3 has the last 4 vertices, owns 4, 6 and 7
406 {3, 3, {2, 0}, {5, 8}},
407 {3, -1, {2, 1}, {6, 9}},
408 {3, 3, {3, 0}, {7, 10}},
409 {3, 3, {3, 1}, {8, 11}},
410 };
412 std::move(inVertexList),
414 "inMesh"};
415 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 1 is empty
416 {0, -1, {0, 0}, {0, 0}},
417 {0, -1, {0, 1}, {0, 0}},
418 {0, -1, {1, 0}, {0, 0}},
419 {2, -1, {1, 1}, {0, 0}},
420 {2, -1, {2, 0}, {0, 0}},
421 {2, -1, {2, 1}, {0, 0}},
422 {3, -1, {3, 0}, {0, 0}},
423 {3, -1, {3, 1}, {0, 0}}};
425 std::move(outVertexList),
427 "outMesh"};
428
429 ReferenceSpecification ref{// Tests for {0, 1, 2} on the first rank,
430 // second rank (consistent with the outMesh) is empty, ...
431 {0, {1, 4}},
432 {0, {2, 5}},
433 {0, {3, 6}},
434 {2, {4, 7}},
435 {2, {5, 8}},
436 {2, {6, 9}},
437 {3, {7, 10}},
438 {3, {8, 11}}};
439 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
440 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
441 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
442 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
443 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
444 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
445}
446
448BOOST_AUTO_TEST_CASE(DistributedConsistent2DV4)
449{
450 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
452
453 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
454
456 // Rank 0 has no vertices
457 // Rank 1 has the entire mesh, owns a subpart
458 {1, 1, {0, 0}, {1.1}},
459 {1, 1, {0, 1}, {2.5}},
460 {1, 1, {1, 0}, {3}},
461 {1, 1, {1, 1}, {4}},
462 {1, -1, {2, 0}, {5}},
463 {1, -1, {2, 1}, {6}},
464 {1, -1, {3, 0}, {7}},
465 {1, -1, {3, 1}, {8}},
466 // Rank 2 has the entire mesh, owns a subpart
467 {2, -1, {0, 0}, {1.1}},
468 {2, -1, {0, 1}, {2.5}},
469 {2, -1, {1, 0}, {3}},
470 {2, -1, {1, 1}, {4}},
471 {2, 2, {2, 0}, {5}},
472 {2, 2, {2, 1}, {6}},
473 {2, 2, {3, 0}, {7}},
474 {2, 2, {3, 1}, {8}},
475 // Rank 3 has no vertices
476 };
478 std::move(inVertexList),
480 "inMesh"};
481 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 0 and 3 are empty
482 // not same order as input mesh and vertex (2,0) appears twice
483 {1, -1, {2, 0}, {0}},
484 {1, -1, {1, 0}, {0}},
485 {1, -1, {0, 1}, {0}},
486 {1, -1, {1, 1}, {0}},
487 {1, -1, {0, 0}, {0}},
488 {2, -1, {2, 0}, {0}},
489 {2, -1, {2, 1}, {0}},
490 {2, -1, {3, 0}, {0}},
491 {2, -1, {3, 1}, {0}}};
493 std::move(outVertexList),
495 "outMesh"};
496 ReferenceSpecification ref{{1, {5}},
497 {1, {3}},
498 {1, {2.5}},
499 {1, {4}},
500 {1, {1.1}},
501 {2, {5}},
502 {2, {6}},
503 {2, {7}},
504 {2, {8}}};
505 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
506 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
507 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
508 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
509 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
510 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
511}
512
513// same as 2DV4, but all ranks have vertices
514BOOST_AUTO_TEST_CASE(DistributedConsistent2DV5)
515{
516 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
518
519 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
520 std::vector<VertexSpecification> inVertexList = {
521 // Every rank has the entire mesh and owns a subpart
522 {0, 0, {0, 0}, {1.1}},
523 {0, 0, {0, 1}, {2.5}},
524 {0, -1, {1, 0}, {3}},
525 {0, -1, {1, 1}, {4}},
526 {0, -1, {2, 0}, {5}},
527 {0, -1, {2, 1}, {6}},
528 {0, -1, {3, 0}, {7}},
529 {0, -1, {3, 1}, {8}},
530 {1, -1, {0, 0}, {1.1}},
531 {1, -1, {0, 1}, {2.5}},
532 {1, 1, {1, 0}, {3}},
533 {1, 1, {1, 1}, {4}},
534 {1, -1, {2, 0}, {5}},
535 {1, -1, {2, 1}, {6}},
536 {1, -1, {3, 0}, {7}},
537 {1, -1, {3, 1}, {8}},
538 {2, -1, {0, 0}, {1.1}},
539 {2, -1, {0, 1}, {2.5}},
540 {2, -1, {1, 0}, {3}},
541 {2, -1, {1, 1}, {4}},
542 {2, 2, {2, 0}, {5}},
543 {2, 2, {2, 1}, {6}},
544 {2, -1, {3, 0}, {7}},
545 {2, -1, {3, 1}, {8}},
546 {3, -1, {0, 0}, {1.1}},
547 {3, -1, {0, 1}, {2.5}},
548 {3, -1, {1, 0}, {3}},
549 {3, -1, {1, 1}, {4}},
550 {3, -1, {2, 0}, {5}},
551 {3, -1, {2, 1}, {6}},
552 {3, 3, {3, 0}, {7}},
553 {3, 3, {3, 1}, {8}},
554 };
556 std::move(inVertexList),
558 "inMesh"};
559 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 0 and 3 are empty
560 // not same order as input mesh and vertex (2,0) appears twice
561 {1, -1, {2, 0}, {0}},
562 {1, -1, {1, 0}, {0}},
563 {1, -1, {0, 1}, {0}},
564 {1, -1, {1, 1}, {0}},
565 {1, -1, {0, 0}, {0}},
566 {2, -1, {2, 0}, {0}},
567 {2, -1, {2, 1}, {0}},
568 {2, -1, {3, 0}, {0}},
569 {2, -1, {3, 1}, {0}}};
571 std::move(outVertexList),
573 "outMesh"};
574 ReferenceSpecification ref{{1, {5}},
575 {1, {3}},
576 {1, {2.5}},
577 {1, {4}},
578 {1, {1.1}},
579 {2, {5}},
580 {2, {6}},
581 {2, {7}},
582 {2, {8}}};
583 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
584 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
585 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
586 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
587 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
588 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
589}
590
592BOOST_AUTO_TEST_CASE(DistributedConsistent2DV6,
593 *boost::unit_test::tolerance(1e-7))
594{
595 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
597 std::vector<int> globalIndexOffsets = {0, 0, 0, 0};
598
600 // Rank 0 has no vertices
601 // Rank 1 has the entire mesh, owns a subpart
602 {1, 1, {0, 0}, {1}},
603 {1, 1, {0, 1}, {2}},
604 {1, 1, {1, 0}, {3}},
605 {1, 1, {1, 1}, {4}},
606 {1, -1, {2, 0}, {5}},
607 {1, -1, {2, 1}, {6}},
608 {1, -1, {3, 0}, {7}},
609 {1, -1, {3, 1}, {8}},
610 // Rank 2 has the entire mesh, owns a subpart
611 {2, -1, {0, 0}, {1}},
612 {2, -1, {0, 1}, {2}},
613 {2, -1, {1, 0}, {3}},
614 {2, -1, {1, 1}, {4}},
615 {2, 2, {2, 0}, {5}},
616 {2, 2, {2, 1}, {6}},
617 {2, 2, {3, 0}, {7}},
618 {2, 2, {3, 1}, {8}},
619 // Rank 3 has no vertices
620 };
622 std::move(inVertexList),
624 "inMesh"};
625 std::vector<VertexSpecification> outVertexList{// The outMesh is local, rank 0 and 3 are empty
626 // not same order as input mesh and vertex (2,0) appears twice
627 {1, -1, {2, 0}, {0}},
628 {1, -1, {1, 0}, {0}},
629 {1, -1, {0, 1}, {0}},
630 {1, -1, {1, 1}, {0}},
631 {1, -1, {0, 0}, {0}},
632 {2, -1, {2, 0}, {0}},
633 {2, -1, {2, 1}, {0}},
634 {2, -1, {3, 0}, {0}},
635 {2, -1, {3, 1}, {0}}};
637 std::move(outVertexList),
639 "outMesh"};
640 ReferenceSpecification ref{{1, {5}},
641 {1, {3}},
642 {1, {2}},
643 {1, {4}},
644 {1, {1}},
645 {2, {5}},
646 {2, {6}},
647 {2, {7}},
648 {2, {8}}};
649 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_on(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::ON);
650 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
651 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_sep(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
652 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
653 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_off(Mapping::CONSISTENT, 2, fct, {{false, false, false}}, Polynomial::OFF);
654 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
655}
656
658BOOST_AUTO_TEST_CASE(DistributedConservative2DV1)
659{
660 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
661 Multiquadrics fct(5.0);
662
663 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local
664 {0, -1, {0, 0}, {1}},
665 {0, -1, {0, 1}, {2}},
666 {1, -1, {1, 0}, {3}},
667 {1, -1, {1, 1}, {4}},
668 {2, -1, {2, 0}, {5}},
669 {2, -1, {2, 1}, {6}},
670 {3, -1, {3, 0}, {7}},
671 {3, -1, {3, 1}, {8}}};
673 std::move(inVertexList),
675 "inMesh"};
676 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed
677 {-1, 0, {0, 0}, {0}},
678 {-1, 0, {0, 1}, {0}},
679 {-1, 1, {1, 0}, {0}},
680 {-1, 1, {1, 1}, {0}},
681 {-1, 2, {2, 0}, {0}},
682 {-1, 2, {2, 1}, {0}},
683 {-1, 3, {3, 0}, {0}},
684 {-1, 3, {3, 1}, {0}}};
686 std::move(outVertexList),
688 "outMesh"};
689 ReferenceSpecification ref{// Tests for {0, 1, 0, 0, 0, 0, 0, 0} on the first rank,
690 // {0, 0, 2, 3, 0, 0, 0, 0} on the second, ...
691 {0, {1}},
692 {0, {2}},
693 {0, {0}},
694 {0, {0}},
695 {0, {0}},
696 {0, {0}},
697 {0, {0}},
698 {0, {0}},
699 {1, {0}},
700 {1, {0}},
701 {1, {3}},
702 {1, {4}},
703 {1, {0}},
704 {1, {0}},
705 {1, {0}},
706 {1, {0}},
707 {2, {0}},
708 {2, {0}},
709 {2, {0}},
710 {2, {0}},
711 {2, {5}},
712 {2, {6}},
713 {2, {0}},
714 {2, {0}},
715 {3, {0}},
716 {3, {0}},
717 {3, {0}},
718 {3, {0}},
719 {3, {0}},
720 {3, {0}},
721 {3, {7}},
722 {3, {8}}};
723
724 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
725 testDistributed(context, mapping_on, in, out, ref, context.rank * 2);
726 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
727 testDistributed(context, mapping_sep, in, out, ref, context.rank * 2);
728 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
729 testDistributed(context, mapping_off, in, out, ref, context.rank * 2);
730}
731
733BOOST_AUTO_TEST_CASE(DistributedConservative2DV1Vector)
734{
735 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
736 Multiquadrics fct(5.0);
737 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local
738 {0, -1, {0, 0}, {1, 4}},
739 {0, -1, {0, 1}, {2, 5}},
740 {1, -1, {1, 0}, {3, 6}},
741 {1, -1, {1, 1}, {4, 7}},
742 {2, -1, {2, 0}, {5, 8}},
743 {2, -1, {2, 1}, {6, 9}},
744 {3, -1, {3, 0}, {7, 10}},
745 {3, -1, {3, 1}, {8, 11}}};
747 std::move(inVertexList),
749 "inMesh"};
750 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed
751 {-1, 0, {0, 0}, {0, 0}},
752 {-1, 0, {0, 1}, {0, 0}},
753 {-1, 1, {1, 0}, {0, 0}},
754 {-1, 1, {1, 1}, {0, 0}},
755 {-1, 2, {2, 0}, {0, 0}},
756 {-1, 2, {2, 1}, {0, 0}},
757 {-1, 3, {3, 0}, {0, 0}},
758 {-1, 3, {3, 1}, {0, 0}}};
760 std::move(outVertexList),
762 "outMesh"};
763 ReferenceSpecification ref{// Tests for {0, 1, 0, 0, 0, 0, 0, 0} on the first rank,
764 // {0, 0, 2, 3, 0, 0, 0, 0} on the second, ...
765 {0, {1, 4}},
766 {0, {2, 5}},
767 {0, {0, 0}},
768 {0, {0, 0}},
769 {0, {0, 0}},
770 {0, {0, 0}},
771 {0, {0, 0}},
772 {0, {0, 0}},
773 {1, {0, 0}},
774 {1, {0, 0}},
775 {1, {3, 6}},
776 {1, {4, 7}},
777 {1, {0, 0}},
778 {1, {0, 0}},
779 {1, {0, 0}},
780 {1, {0, 0}},
781 {2, {0, 0}},
782 {2, {0, 0}},
783 {2, {0, 0}},
784 {2, {0, 0}},
785 {2, {5, 8}},
786 {2, {6, 9}},
787 {2, {0, 0}},
788 {2, {0, 0}},
789 {3, {0, 0}},
790 {3, {0, 0}},
791 {3, {0, 0}},
792 {3, {0, 0}},
793 {3, {0, 0}},
794 {3, {0, 0}},
795 {3, {7, 10}},
796 {3, {8, 11}}};
797 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
798 testDistributed(context, mapping_on, in, out, ref, context.rank * 2);
799 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
800 testDistributed(context, mapping_sep, in, out, ref, context.rank * 2);
801 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
802 testDistributed(context, mapping_off, in, out, ref, context.rank * 2);
803}
804
806BOOST_AUTO_TEST_CASE(DistributedConservative2DV2)
807{
808 PRECICE_TEST(""_on(4_ranks).setupIntraComm())
809 Multiquadrics fct(5.0);
810
811 std::vector<int> globalIndexOffsets = {0, 0, 4, 6};
812
813 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local but rank 0 has no vertices
814 {1, -1, {0, 0}, {1}},
815 {1, -1, {0, 1}, {2}},
816 {1, -1, {1, 0}, {3}},
817 {1, -1, {1, 1}, {4}},
818 {2, -1, {2, 0}, {5}},
819 {2, -1, {2, 1}, {6}},
820 {3, -1, {3, 0}, {7}},
821 {3, -1, {3, 1}, {8}}};
823 std::move(inVertexList),
825 "inMesh"};
826 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed, rank 0 owns no vertex
827 {-1, 1, {0, 0}, {0}},
828 {-1, 1, {0, 1}, {0}},
829 {-1, 1, {1, 0}, {0}},
830 {-1, 1, {1, 1}, {0}},
831 {-1, 2, {2, 0}, {0}},
832 {-1, 2, {2, 1}, {0}},
833 {-1, 3, {3, 0}, {0}},
834 {-1, 3, {3, 1}, {0}}};
836 std::move(outVertexList),
838 "outMesh"};
839 ReferenceSpecification ref{// Tests for {0, 0, 0, 0, 0, 0, 0, 0} on the first rank,
840 // {1, 2, 2, 3, 0, 0, 0, 0} on the second, ...
841 {0, {0}},
842 {0, {0}},
843 {0, {0}},
844 {0, {0}},
845 {0, {0}},
846 {0, {0}},
847 {0, {0}},
848 {0, {0}},
849 {1, {1}},
850 {1, {2}},
851 {1, {3}},
852 {1, {4}},
853 {1, {0}},
854 {1, {0}},
855 {1, {0}},
856 {1, {0}},
857 {2, {0}},
858 {2, {0}},
859 {2, {0}},
860 {2, {0}},
861 {2, {5}},
862 {2, {6}},
863 {2, {0}},
864 {2, {0}},
865 {3, {0}},
866 {3, {0}},
867 {3, {0}},
868 {3, {0}},
869 {3, {0}},
870 {3, {0}},
871 {3, {7}},
872 {3, {8}}};
873
874 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
875 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
876 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
877 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
878 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
879 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
880}
881
883BOOST_AUTO_TEST_CASE(DistributedConservative2DV3)
884{
885 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
886 Multiquadrics fct(2.0);
887 std::vector<int> globalIndexOffsets = {0, 0, 3, 5};
888
889 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local but rank 0 has no vertices
890 {1, -1, {0, 0}, {1}},
891 {1, -1, {1, 0}, {3}},
892 {1, -1, {1, 1}, {4}},
893 {2, -1, {2, 0}, {5}},
894 {2, -1, {2, 1}, {6}},
895 {3, -1, {3, 0}, {7}},
896 {3, -1, {3, 1}, {8}}}; // Sum of all vertices is 34
898 std::move(inVertexList),
900 "inMesh"};
901 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed, rank 0 owns no vertex
902 {-1, 1, {0, 0}, {0}},
903 {-1, 1, {0, 1}, {0}},
904 {-1, 1, {1, 0}, {0}},
905 {-1, 1, {1, 1}, {0}},
906 {-1, 2, {2, 0}, {0}},
907 {-1, 2, {2, 1}, {0}},
908 {-1, 3, {3, 0}, {0}},
909 {-1, 3, {3, 1}, {0}}};
911 std::move(outVertexList),
913 "outMesh"};
914 ReferenceSpecification ref{// Tests for {0, 0, 0, 0, 0, 0, 0, 0} on the first rank,
915 // {1, 2, 2, 3, 0, 0, 0, 0} on the second, ...
916 {0, {0}},
917 {0, {0}},
918 {0, {0}},
919 {0, {0}},
920 {0, {0}},
921 {0, {0}},
922 {0, {0}},
923 {0, {0}},
924 {1, {1}},
925 {1, {0}},
926 {1, {3}},
927 {1, {4}},
928 {1, {0}},
929 {1, {0}},
930 {1, {0}},
931 {1, {0}},
932 {2, {0}},
933 {2, {0}},
934 {2, {0}},
935 {2, {0}},
936 {2, {5}},
937 {2, {6}},
938 {2, {0}},
939 {2, {0}},
940 {3, {0}},
941 {3, {0}},
942 {3, {0}},
943 {3, {0}},
944 {3, {0}},
945 {3, {0}},
946 {3, {7}},
947 {3, {8}}};
948 // Sum of reference is also 34
949 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
950 testDistributed(context, mapping_on, in, out, ref, globalIndexOffsets.at(context.rank));
951 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
952 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank));
953 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
954 testDistributed(context, mapping_off, in, out, ref, globalIndexOffsets.at(context.rank));
955}
956
958BOOST_AUTO_TEST_CASE(DistributedConservative2DV4,
959 *boost::unit_test::tolerance(1e-6))
960{
961 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
963 std::vector<int> globalIndexOffsets = {0, 2, 4, 6};
964
965 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local
966 {0, -1, {0, 0}, {1}},
967 {0, -1, {0, 1}, {2}},
968 {1, -1, {1, 0}, {3}},
969 {1, -1, {1, 1}, {4}},
970 {2, -1, {2, 0}, {5}},
971 {2, -1, {2, 1}, {6}},
972 {3, -1, {3, 0}, {7}},
973 {3, -1, {3, 1}, {8}}}; // Sum is 36
975 std::move(inVertexList),
977 "inMesh"};
978 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed, rank 0 has no vertex at all
979 {-1, 1, {0, 1}, {0}},
980 {-1, 1, {1, 0}, {0}},
981 {-1, 1, {1, 1}, {0}},
982 {-1, 2, {2, 0}, {0}},
983 {-1, 2, {2, 1}, {0}},
984 {-1, 3, {3, 0}, {0}},
985 {-1, 3, {3, 1}, {0}}};
987 std::move(outVertexList),
989 "outMesh"};
990 ReferenceSpecification ref{// Tests for {0, 0, 0, 0, 0, 0, 0, 0} on the first rank,
991 // {2, 3, 4, 3, 0, 0, 0, 0} on the second, ...
992 {0, {0}},
993 {0, {0}},
994 {0, {0}},
995 {0, {0}},
996 {0, {0}},
997 {0, {0}},
998 {0, {0}},
999 {1, {3.1307987056295525}},
1000 {1, {4.0734031184906971}},
1001 {1, {3.0620533966233006}},
1002 {1, {0}},
1003 {1, {0}},
1004 {1, {0}},
1005 {1, {0}},
1006 {2, {0}},
1007 {2, {0}},
1008 {2, {0}},
1009 {2, {4.4582564956060873}},
1010 {2, {5.8784343572772633}},
1011 {2, {0}},
1012 {2, {0}},
1013 {3, {0}},
1014 {3, {0}},
1015 {3, {0}},
1016 {3, {0}},
1017 {3, {0}},
1018 {3, {7.4683403859032156}},
1019 {3, {7.9287135404698859}}}; // Sum is ~36
1020 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
1021 testDistributed(context, mapping_sep, in, out, ref, globalIndexOffsets.at(context.rank), true);
1022 // Polynomial == OFF won't reach the desired accuracy
1023}
1024
1026BOOST_AUTO_TEST_CASE(testDistributedConservative2DV5)
1027{
1028 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1029 Multiquadrics fct(5.0);
1030 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local
1031 {0, -1, {0, 0}, {1}},
1032 {0, -1, {0, 1}, {2}},
1033 {1, -1, {1, 0}, {3}},
1034 {1, -1, {1, 1}, {4}},
1035 {2, -1, {2, 0}, {5}},
1036 {2, -1, {2, 1}, {6}},
1037 {3, -1, {3, 0}, {7}},
1038 {3, -1, {3, 1}, {8}}};
1040 std::move(inVertexList),
1041 meshDims2D,
1042 "inMesh"};
1043 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed and non-contigous
1044 {-1, 0, {0, 0}, {0}},
1045 {-1, 1, {0, 1}, {0}},
1046 {-1, 1, {1, 0}, {0}},
1047 {-1, 0, {1, 1}, {0}},
1048 {-1, 2, {2, 0}, {0}},
1049 {-1, 2, {2, 1}, {0}},
1050 {-1, 3, {3, 0}, {0}},
1051 {-1, 3, {3, 1}, {0}}};
1053 std::move(outVertexList),
1054 meshDims2D,
1055 "outMesh"};
1056 ReferenceSpecification ref{// Tests for {0, 1, 0, 0, 0, 0, 0, 0} on the first rank,
1057 // {0, 0, 2, 3, 0, 0, 0, 0} on the second, ...
1058 {0, {1}},
1059 {0, {0}},
1060 {0, {0}},
1061 {0, {4}},
1062 {0, {0}},
1063 {0, {0}},
1064 {0, {0}},
1065 {0, {0}},
1066 {1, {0}},
1067 {1, {2}},
1068 {1, {3}},
1069 {1, {0}},
1070 {1, {0}},
1071 {1, {0}},
1072 {1, {0}},
1073 {1, {0}},
1074 {2, {0}},
1075 {2, {0}},
1076 {2, {0}},
1077 {2, {0}},
1078 {2, {5}},
1079 {2, {6}},
1080 {2, {0}},
1081 {2, {0}},
1082 {3, {0}},
1083 {3, {0}},
1084 {3, {0}},
1085 {3, {0}},
1086 {3, {0}},
1087 {3, {0}},
1088 {3, {7}},
1089 {3, {8}}};
1090
1091 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
1092 testDistributed(context, mapping_on, in, out, ref, context.rank * 2);
1093 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
1094 testDistributed(context, mapping_sep, in, out, ref, context.rank * 2);
1095 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
1096 testDistributed(context, mapping_off, in, out, ref, context.rank * 2);
1097}
1098
1100BOOST_AUTO_TEST_CASE(testDistributedConservative2DV5Vector)
1101{
1102 PRECICE_TEST(""_on(4_ranks).setupIntraComm());
1103 Multiquadrics fct(5.0);
1104
1105 std::vector<VertexSpecification> inVertexList{// Conservative mapping: The inMesh is local
1106 {0, -1, {0, 0}, {1, 4}},
1107 {0, -1, {0, 1}, {2, 5}},
1108 {1, -1, {1, 0}, {3, 6}},
1109 {1, -1, {1, 1}, {4, 7}},
1110 {2, -1, {2, 0}, {5, 8}},
1111 {2, -1, {2, 1}, {6, 9}},
1112 {3, -1, {3, 0}, {7, 10}},
1113 {3, -1, {3, 1}, {8, 11}}};
1115 std::move(inVertexList),
1116 meshDims2D,
1117 "inMesh"};
1118 std::vector<VertexSpecification> outVertexList{// The outMesh is distributed and non-contigous
1119 {-1, 0, {0, 0}, {0, 0}},
1120 {-1, 1, {0, 1}, {0, 0}},
1121 {-1, 1, {1, 0}, {0, 0}},
1122 {-1, 0, {1, 1}, {0, 0}},
1123 {-1, 2, {2, 0}, {0, 0}},
1124 {-1, 2, {2, 1}, {0, 0}},
1125 {-1, 3, {3, 0}, {0, 0}},
1126 {-1, 3, {3, 1}, {0, 0}}};
1128 std::move(outVertexList),
1129 meshDims2D,
1130 "outMesh"};
1131 ReferenceSpecification ref{// Tests for {0, 1, 0, 0, 0, 0, 0, 0} on the first rank,
1132 // {0, 0, 2, 3, 0, 0, 0, 0} on the second, ...
1133 {0, {1, 4}},
1134 {0, {0, 0}},
1135 {0, {0, 0}},
1136 {0, {4, 7}},
1137 {0, {0, 0}},
1138 {0, {0, 0}},
1139 {0, {0, 0}},
1140 {0, {0, 0}},
1141 {1, {0, 0}},
1142 {1, {2, 5}},
1143 {1, {3, 6}},
1144 {1, {0, 0}},
1145 {1, {0, 0}},
1146 {1, {0, 0}},
1147 {1, {0, 0}},
1148 {1, {0, 0}},
1149 {2, {0, 0}},
1150 {2, {0, 0}},
1151 {2, {0, 0}},
1152 {2, {0, 0}},
1153 {2, {5, 8}},
1154 {2, {6, 9}},
1155 {2, {0, 0}},
1156 {2, {0, 0}},
1157 {3, {0, 0}},
1158 {3, {0, 0}},
1159 {3, {0, 0}},
1160 {3, {0, 0}},
1161 {3, {0, 0}},
1162 {3, {0, 0}},
1163 {3, {7, 10}},
1164 {3, {8, 11}}};
1165 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_on(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::ON);
1166 testDistributed(context, mapping_on, in, out, ref, context.rank * 2);
1167 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_sep(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
1168 testDistributed(context, mapping_sep, in, out, ref, context.rank * 2);
1169 RadialBasisFctMapping<RadialBasisFctSolver<Multiquadrics>> mapping_off(Mapping::CONSERVATIVE, 2, fct, {{false, false, false}}, Polynomial::OFF);
1170 testDistributed(context, mapping_off, in, out, ref, context.rank * 2);
1171}
1172
1173void testTagging(const TestContext &context,
1174 MeshSpecification inMeshSpec,
1175 MeshSpecification outMeshSpec,
1176 MeshSpecification shouldTagFirstRound,
1177 MeshSpecification shouldTagSecondRound,
1178 bool consistent)
1179{
1180 int meshDimension = inMeshSpec.meshDimension;
1181
1182 mesh::PtrMesh inMesh = getDistributedMesh(context, inMeshSpec);
1183 Eigen::VectorXd inValues = getDistributedData(context, inMeshSpec);
1184
1185 mesh::PtrMesh outMesh = getDistributedMesh(context, outMeshSpec);
1186 Eigen::VectorXd outValues = getDistributedData(context, outMeshSpec);
1187
1188 Gaussian fct(4.5); // Support radius approx. 1
1190 RadialBasisFctMapping<RadialBasisFctSolver<Gaussian>> mapping(constr, 2, fct, {{false, false, false}}, Polynomial::SEPARATE);
1191 inMesh->computeBoundingBox();
1192 outMesh->computeBoundingBox();
1193
1194 mapping.setMeshes(inMesh, outMesh);
1195 mapping.tagMeshFirstRound();
1196
1197 for (const auto &v : inMesh->vertices()) {
1198 auto pos = std::find_if(shouldTagFirstRound.vertices.begin(), shouldTagFirstRound.vertices.end(),
1199 [meshDimension, &v](const VertexSpecification &spec) {
1200 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1201 });
1202 bool found = pos != shouldTagFirstRound.vertices.end();
1203 BOOST_TEST(found >= v.isTagged(),
1204 "FirstRound: Vertex " << v << " is tagged, but should not be.");
1205 BOOST_TEST(found <= v.isTagged(),
1206 "FirstRound: Vertex " << v << " is not tagged, but should be.");
1207 }
1208
1209 mapping.tagMeshSecondRound();
1210
1211 for (const auto &v : inMesh->vertices()) {
1212 auto posFirst = std::find_if(shouldTagFirstRound.vertices.begin(), shouldTagFirstRound.vertices.end(),
1213 [meshDimension, &v](const VertexSpecification &spec) {
1214 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1215 });
1216 bool foundFirst = posFirst != shouldTagFirstRound.vertices.end();
1217 auto posSecond = std::find_if(shouldTagSecondRound.vertices.begin(), shouldTagSecondRound.vertices.end(),
1218 [meshDimension, &v](const VertexSpecification &spec) {
1219 return std::equal(spec.position.data(), spec.position.data() + meshDimension, v.getCoords().data());
1220 });
1221 bool foundSecond = posSecond != shouldTagSecondRound.vertices.end();
1222 BOOST_TEST(foundFirst <= v.isTagged(), "SecondRound: Vertex " << v
1223 << " is not tagged, but should be from the first round.");
1224 BOOST_TEST(foundSecond <= v.isTagged(), "SecondRound: Vertex " << v
1225 << " is not tagged, but should be.");
1226 BOOST_TEST((foundSecond or foundFirst) >= v.isTagged(), "SecondRound: Vertex " << v
1227 << " is tagged, but should not be.");
1228 }
1229}
1230
1231BOOST_AUTO_TEST_CASE(testTagFirstRound)
1232{
1233 PRECICE_TEST(""_on(4_ranks).setupIntraComm())
1234 // *
1235 // + <-- owned
1236 //* * x * *
1237 // *
1238 // *
1240 {0, -1, {0, 0}, {0}}};
1241 MeshSpecification outMeshSpec{
1242 std::move(outVertexList),
1243 meshDims2D,
1244 "outMesh"};
1246 {0, -1, {-1, 0}, {1}}, // inside
1247 {0, -1, {-2, 0}, {1}}, // outside
1248 {0, 0, {1, 0}, {1}}, // inside, owner
1249 {0, -1, {2, 0}, {1}}, // outside
1250 {0, -1, {0, -1}, {1}}, // inside
1251 {0, -1, {0, -2}, {1}}, // outside
1252 {0, -1, {0, 1}, {1}}, // inside
1253 {0, -1, {0, 2}, {1}} // outside
1254 };
1255 MeshSpecification inMeshSpec{
1256 std::move(inVertexList),
1257 meshDims2D,
1258 "inMesh"};
1259 std::vector<VertexSpecification> firstRoundVertices = {
1260 {0, -1, {-1, 0}, {1}},
1261 {0, -1, {1, 0}, {1}},
1262 {0, -1, {0, -1}, {1}},
1263 {0, -1, {0, 1}, {1}}};
1264
1265 MeshSpecification shouldTagFirstRound = {
1266 firstRoundVertices,
1267 static_cast<int>(firstRoundVertices.at(0).position.size()),
1268 ""};
1269 std::vector<VertexSpecification> secondRoundVertices = {
1270 {0, -1, {2, 0}, {1}}};
1271 MeshSpecification shouldTagSecondRound = {
1272 secondRoundVertices,
1273 static_cast<int>(secondRoundVertices.at(0).position.size()),
1274 ""};
1275 testTagging(context, inMeshSpec, outMeshSpec, shouldTagFirstRound, shouldTagSecondRound, true);
1276 // For conservative just swap meshes.
1277 testTagging(context, outMeshSpec, inMeshSpec, shouldTagFirstRound, shouldTagSecondRound, false);
1278}
1279
1280BOOST_AUTO_TEST_SUITE_END() // Parallel
1281
1283
1284#undef doLocalCode
1285#define doLocalCode(Type, function, polynomial) \
1286 { \
1287 RadialBasisFctMapping<RadialBasisFctSolver<Type>> consistentMap2D(Mapping::CONSISTENT, 2, function, {{false, false, false}}, polynomial); \
1288 perform2DTestConsistentMapping(consistentMap2D); \
1289 RadialBasisFctMapping<RadialBasisFctSolver<Type>> consistentMap2DVector(Mapping::CONSISTENT, 2, function, {{false, false, false}}, polynomial); \
1290 perform2DTestConsistentMappingVector(consistentMap2DVector); \
1291 RadialBasisFctMapping<RadialBasisFctSolver<Type>> consistentMap3D(Mapping::CONSISTENT, 3, function, {{false, false, false}}, polynomial); \
1292 perform3DTestConsistentMapping(consistentMap3D); \
1293 RadialBasisFctMapping<RadialBasisFctSolver<Type>> scaledConsistentMap2D(Mapping::SCALED_CONSISTENT_SURFACE, 2, function, {{false, false, false}}, polynomial); \
1294 perform2DTestScaledConsistentMapping(scaledConsistentMap2D); \
1295 RadialBasisFctMapping<RadialBasisFctSolver<Type>> scaledConsistentMap3D(Mapping::SCALED_CONSISTENT_SURFACE, 3, function, {{false, false, false}}, polynomial); \
1296 perform3DTestScaledConsistentMapping(scaledConsistentMap3D); \
1297 RadialBasisFctMapping<RadialBasisFctSolver<Type>> conservativeMap2D(Mapping::CONSERVATIVE, 2, function, {{false, false, false}}, polynomial); \
1298 perform2DTestConservativeMapping(conservativeMap2D); \
1299 RadialBasisFctMapping<RadialBasisFctSolver<Type>> conservativeMap2DVector(Mapping::CONSERVATIVE, 2, function, {{false, false, false}}, polynomial); \
1300 perform2DTestConservativeMappingVector(conservativeMap2DVector); \
1301 RadialBasisFctMapping<RadialBasisFctSolver<Type>> conservativeMap3D(Mapping::CONSERVATIVE, 3, function, {{false, false, false}}, polynomial); \
1302 perform3DTestConservativeMapping(conservativeMap3D); \
1303 }
1304
1305BOOST_AUTO_TEST_CASE(MapThinPlateSplines)
1306{
1307 PRECICE_TEST(1_rank);
1308 ThinPlateSplines fct;
1309 doLocalCode(ThinPlateSplines, fct, Polynomial::ON);
1310 doLocalCode(ThinPlateSplines, fct, Polynomial::SEPARATE);
1311}
1312
1313BOOST_AUTO_TEST_CASE(MapMultiquadrics)
1314{
1315 PRECICE_TEST(1_rank);
1316 Multiquadrics fct(1e-3);
1317 doLocalCode(Multiquadrics, fct, Polynomial::ON);
1318 doLocalCode(Multiquadrics, fct, Polynomial::SEPARATE);
1319}
1320
1321BOOST_AUTO_TEST_CASE(MapInverseMultiquadrics)
1322{
1323 PRECICE_TEST(1_rank);
1324 InverseMultiquadrics fct(1e-3);
1325 doLocalCode(InverseMultiquadrics, fct, Polynomial::SEPARATE);
1326}
1327
1328BOOST_AUTO_TEST_CASE(MapVolumeSplines)
1329{
1330 PRECICE_TEST(1_rank);
1331 VolumeSplines fct;
1332 doLocalCode(VolumeSplines, fct, Polynomial::ON);
1333 doLocalCode(VolumeSplines, fct, Polynomial::SEPARATE);
1334}
1335
1337{
1338 PRECICE_TEST(1_rank);
1339 Gaussian fct(1.0);
1340 doLocalCode(Gaussian, fct, Polynomial::SEPARATE);
1341}
1342
1343BOOST_AUTO_TEST_CASE(MapCompactThinPlateSplinesC2)
1344{
1345 PRECICE_TEST(1_rank);
1346 double supportRadius = 1.2;
1347 CompactThinPlateSplinesC2 fct(supportRadius);
1348 doLocalCode(CompactThinPlateSplinesC2, fct, Polynomial::SEPARATE);
1349}
1350
1351BOOST_AUTO_TEST_CASE(MapCompactPolynomialC0)
1352{
1353 PRECICE_TEST(1_rank);
1354 double supportRadius = 1.2;
1355 CompactPolynomialC0 fct(supportRadius);
1356 doLocalCode(CompactPolynomialC0, fct, Polynomial::SEPARATE);
1357}
1358
1359BOOST_AUTO_TEST_CASE(MapCompactPolynomialC2)
1360{
1361 PRECICE_TEST(1_rank);
1362 double supportRadius = 1.2;
1363 CompactPolynomialC2 fct(supportRadius);
1364 doLocalCode(CompactPolynomialC2, fct, Polynomial::SEPARATE);
1365}
1366
1367BOOST_AUTO_TEST_CASE(MapCompactPolynomialC4)
1368{
1369 PRECICE_TEST(1_rank);
1370 double supportRadius = 1.2;
1371 CompactPolynomialC4 fct(supportRadius);
1372 doLocalCode(CompactPolynomialC4, fct, Polynomial::SEPARATE);
1373}
1374
1375BOOST_AUTO_TEST_CASE(MapCompactPolynomialC6)
1376{
1377 PRECICE_TEST(1_rank);
1378 double supportRadius = 1.2;
1379 CompactPolynomialC6 fct(supportRadius);
1380 doLocalCode(CompactPolynomialC6, fct, Polynomial::SEPARATE);
1381}
1382
1383BOOST_AUTO_TEST_CASE(MapCompactPolynomialC8)
1384{
1385 PRECICE_TEST(1_rank);
1386 double supportRadius = 1.2;
1387 CompactPolynomialC8 fct(supportRadius);
1388 doLocalCode(CompactPolynomialC8, fct, Polynomial::SEPARATE);
1389}
1390#undef doLocalCode
1391
1393{
1394 using Eigen::Vector2d;
1395 int dimensions = 2;
1396
1397 bool xDead = false;
1398 bool yDead = true;
1399 bool zDead = false;
1400
1401 ThinPlateSplines fct;
1402 RadialBasisFctMapping<RadialBasisFctSolver<ThinPlateSplines>> mapping(constraint, dimensions, fct,
1403 {{xDead, yDead, zDead}}, polynomial);
1404
1405 // Create mesh to map from
1406 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
1407 inMesh->createVertex(Vector2d(0.0, 1.0));
1408 inMesh->createVertex(Vector2d(1.0, 1.0));
1409 inMesh->createVertex(Vector2d(2.0, 1.0));
1410 inMesh->createVertex(Vector2d(3.0, 1.0));
1411 addGlobalIndex(inMesh);
1412 inMesh->setGlobalNumberOfVertices(inMesh->nVertices());
1413
1414 Eigen::VectorXd values(4);
1415 values << 1.0, 2.0, 2.0, 1.0;
1416
1417 // Create mesh to map to
1418 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
1419 outMesh->createVertex(Vector2d(0, 1.));
1420 outMesh->createVertex(Vector2d(3, 1.));
1421 outMesh->createVertex(Vector2d(1.3, 1.));
1422 outMesh->createVertex(Vector2d(5, 1.));
1423 addGlobalIndex(outMesh);
1424 outMesh->setGlobalNumberOfVertices(outMesh->nVertices());
1425
1426 // Setup mapping with mapping coordinates and geometry used
1427 mapping.setMeshes(inMesh, outMesh);
1428 BOOST_TEST(mapping.hasComputedMapping() == false);
1429
1430 mapping.computeMapping();
1431 Eigen::VectorXd outValues(4);
1432 time::Sample inSample(1, values);
1433 mapping.map(inSample, outValues);
1434 BOOST_TEST(mapping.hasComputedMapping() == true);
1435 if (constraint == Mapping::CONSISTENT) {
1436 if (polynomial == Polynomial::OFF) {
1437 BOOST_TEST(testing::equals(outValues(2), 2.0522549299731567, 1e-7));
1438 } else if (polynomial == Polynomial::SEPARATE) {
1439 BOOST_TEST(testing::equals(outValues(2), 2.0896514371485777, 1e-7));
1440 } else {
1441 BOOST_TEST(testing::equals(outValues(2), 2.1180354377884774, 1e-7));
1442 }
1443 } else {
1444 if (polynomial == Polynomial::OFF) {
1445 BOOST_TEST(testing::equals(outValues(1), 1.8471144693068295, 1e-7));
1446 } else if (polynomial == Polynomial::SEPARATE) {
1447 BOOST_TEST(testing::equals(outValues(1), 1.8236736422730249, 1e-7));
1448 } else {
1449 BOOST_TEST(testing::equals(outValues(1), 1.7587181970483183, 1e-7));
1450 }
1451 }
1452}
1453
1455{
1456 using Eigen::Vector3d;
1457 int dimensions = 3;
1458
1459 ThinPlateSplines fct;
1460 bool xDead = false;
1461 bool yDead = true;
1462 bool zDead = false;
1464 Mapping mapping(constraint, dimensions, fct, {{xDead, yDead, zDead}}, polynomial);
1465
1466 // Create mesh to map from
1467 mesh::PtrMesh inMesh(new mesh::Mesh("InMesh", dimensions, testing::nextMeshID()));
1468 inMesh->createVertex(Vector3d(0.0, 3.0, 0.0));
1469 inMesh->createVertex(Vector3d(1.0, 3.0, 0.0));
1470 inMesh->createVertex(Vector3d(0.0, 3.0, 1.0));
1471 inMesh->createVertex(Vector3d(1.0, 3.0, 1.0));
1472 addGlobalIndex(inMesh);
1473 inMesh->setGlobalNumberOfVertices(inMesh->nVertices());
1474
1475 Eigen::VectorXd values(4);
1476 values << 1.0, 2.0, 3.0, 4.0;
1477
1478 // Create mesh to map to
1479 mesh::PtrMesh outMesh(new mesh::Mesh("OutMesh", dimensions, testing::nextMeshID()));
1480 outMesh->createVertex(Vector3d(0.0, 2.9, 0.0));
1481 outMesh->createVertex(Vector3d(0.8, 2.9, 0.1));
1482 outMesh->createVertex(Vector3d(0.1, 2.9, 0.9));
1483 outMesh->createVertex(Vector3d(1.1, 2.9, 1.1));
1484 addGlobalIndex(outMesh);
1485 outMesh->setGlobalNumberOfVertices(outMesh->nVertices());
1486
1487 // Setup mapping with mapping coordinates and geometry used
1488 mapping.setMeshes(inMesh, outMesh);
1489 BOOST_TEST(mapping.hasComputedMapping() == false);
1490
1491 mapping.computeMapping();
1492 Eigen::VectorXd outValues(4);
1493 time::Sample inSample(1, values);
1494 mapping.map(inSample, outValues);
1495 BOOST_TEST(mapping.hasComputedMapping() == true);
1496
1497 if (constraint == Mapping::CONSISTENT) {
1498 if (polynomial == Polynomial::OFF) {
1499 const double tolerance = 1e-7;
1500 Eigen::VectorXd expected(4);
1501 expected << 1.0, -0.454450524334, 0.99146426249, 6.98958304876;
1502 BOOST_TEST(testing::equals(outValues, expected, tolerance));
1503 } else {
1504 Eigen::VectorXd expected(4);
1505 expected << 1.0, 2.0, 2.9, 4.3;
1506 BOOST_TEST(testing::equals(outValues, expected));
1507 }
1508 } else {
1509 if (polynomial == Polynomial::OFF) {
1510 const double tolerance = 1e-6;
1511 Eigen::VectorXd expected(4);
1512 expected << 1.17251596926, 4.10368825944, 3.56931954192, 3.40160932341;
1513 BOOST_TEST(testing::equals(outValues, expected, tolerance));
1514 } else if (polynomial == Polynomial::ON) {
1515 const double tolerance = 1e-6;
1516 Eigen::VectorXd expected(4);
1517 expected << 0.856701171969, 2.38947124326, 3.34078733786, 3.41304024691;
1518 BOOST_TEST(testing::equals(outValues, expected, tolerance));
1519 } else {
1520 const double tolerance = 1e-6;
1521 Eigen::VectorXd expected(4);
1522 expected << 0.380480856704, 2.83529451713, 3.73088270249, 3.05334192368;
1523 BOOST_TEST(testing::equals(outValues, expected, tolerance));
1524 }
1525 }
1526}
1527
1528BOOST_AUTO_TEST_CASE(DeadAxis2Consistent)
1529{
1530 PRECICE_TEST(1_rank);
1531 testDeadAxis2d(Polynomial::ON, Mapping::CONSISTENT);
1532 testDeadAxis2d(Polynomial::OFF, Mapping::CONSISTENT);
1533 testDeadAxis2d(Polynomial::SEPARATE, Mapping::CONSISTENT);
1534}
1535
1536BOOST_AUTO_TEST_CASE(DeadAxis3DConsistent)
1537{
1538 PRECICE_TEST(1_rank);
1539 testDeadAxis3d(Polynomial::ON, Mapping::CONSISTENT);
1540 testDeadAxis3d(Polynomial::OFF, Mapping::CONSISTENT);
1541 testDeadAxis3d(Polynomial::SEPARATE, Mapping::CONSISTENT);
1542}
1543
1544BOOST_AUTO_TEST_CASE(DeadAxis2Conservative)
1545{
1546 PRECICE_TEST(1_rank);
1547 testDeadAxis2d(Polynomial::ON, Mapping::CONSERVATIVE);
1548 testDeadAxis2d(Polynomial::OFF, Mapping::CONSERVATIVE);
1549 testDeadAxis2d(Polynomial::SEPARATE, Mapping::CONSERVATIVE);
1550}
1551
1552BOOST_AUTO_TEST_CASE(DeadAxis3DConervative)
1553{
1554 PRECICE_TEST(1_rank);
1555 testDeadAxis3d(Polynomial::ON, Mapping::CONSERVATIVE);
1556 testDeadAxis3d(Polynomial::OFF, Mapping::CONSERVATIVE);
1557 testDeadAxis3d(Polynomial::SEPARATE, Mapping::CONSERVATIVE);
1558}
1559
1560BOOST_AUTO_TEST_SUITE_END() // Serial
1561
1562BOOST_AUTO_TEST_SUITE_END() // RadialBasisFunctionMapping
std::ostream & out
BOOST_AUTO_TEST_SUITE(PreProcess)
BOOST_AUTO_TEST_SUITE_END()
void addGlobalIndex(mesh::PtrMesh &mesh, int offset=0)
void testTagging(const TestContext &context, MeshSpecification inMeshSpec, MeshSpecification outMeshSpec, MeshSpecification shouldTagFirstRound, MeshSpecification shouldTagSecondRound, bool consistent)
void testDeadAxis3d(Polynomial polynomial, Mapping::Constraint constraint)
void testDistributed(const TestContext &context, Mapping &mapping, MeshSpecification inMeshSpec, MeshSpecification outMeshSpec, ReferenceSpecification referenceSpec, int inGlobalIndexOffset=0, bool meshIsSmaller=false)
void testDeadAxis2d(Polynomial polynomial, Mapping::Constraint constraint)
constexpr int meshDims2D
BOOST_AUTO_TEST_CASE(DistributedConsistent2DV1)
Test with a homogeneous distribution of mesh among ranks.
Eigen::VectorXd getDistributedData(const TestContext &context, MeshSpecification const &meshSpec)
mesh::PtrMesh getDistributedMesh(const TestContext &context, MeshSpecification &meshSpec, int globalIndexOffset=0, bool meshIsSmaller=false)
#define doLocalCode(Type, function, polynomial)
unsigned int index
#define PRECICE_TEST(...)
Definition Testing.hpp:27
T at(T... args)
T back(T... args)
T begin(T... args)
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Wendland radial basis function with compact support.
Radial basis function with compact support.
Radial basis function with global and compact support.
Radial basis function with global support.
Abstract base class for mapping of data from one mesh to another.
Definition Mapping.hpp:15
Constraint
Specifies additional constraints for a mapping.
Definition Mapping.hpp:29
void setMeshes(const mesh::PtrMesh &input, const mesh::PtrMesh &output)
Sets input and output meshes carrying data to be mapped.
Definition Mapping.cpp:27
virtual void computeMapping()=0
Computes the mapping coefficients from the in- and output mesh.
bool hasComputedMapping() const
Returns true, if the mapping has been computed.
Definition Mapping.cpp:252
void map(int inputDataID, int outputDataID)
Definition Mapping.cpp:126
Radial basis function with global support.
Mapping with radial basis functions.
Radial basis function with global support.
Radial basis function with global support.
Container and creator for meshes.
Definition Mesh.hpp:39
void setGlobalNumberOfVertices(int num)
Definition Mesh.hpp:264
VertexContainer & vertices()
Returns modifieable container holding all vertices.
Definition Mesh.cpp:53
std::size_t nVertices() const
Returns the number of vertices.
Definition Mesh.cpp:63
void computeBoundingBox()
Computes the boundingBox for the vertices.
Definition Mesh.cpp:242
Vertex & createVertex(const Eigen::VectorXd &coords)
Creates and initializes a Vertex object.
Definition Mesh.cpp:103
void setOwner(bool owner)
Definition Vertex.cpp:27
Rank rank
the rank of the current participant
T end(T... args)
T find_if(T... args)
contains data mapping from points to meshes.
Polynomial
How to handle the polynomial?
provides Mesh, Data and primitives.
contains the testing framework.
Definition helper.hpp:9
boost::test_tools::predicate_result equals(const std::vector< float > &VectorA, const std::vector< float > &VectorB, float tolerance)
equals to be used in tests. Compares two std::vectors using a given tolerance. Prints both operands o...
Definition Testing.cpp:65
Main namespace of the precice library.
T size(T... args)
std::vector< VertexSpecification > vertices
Holds rank, owner, position and value of a single vertex.