preCICE v3.2.0
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
AccelerationSerialTest.cpp
Go to the documentation of this file.
1#include <Eigen/Core>
2#include <algorithm>
17#include "testing/Meshes.hpp"
19#include "testing/Testing.hpp"
21
22using namespace precice;
23using namespace precice::acceleration;
24
26
27BOOST_AUTO_TEST_SUITE(AccelerationTests)
28
31
32 // AccelerationSerialTestsFixture() {}
33};
34
35BOOST_FIXTURE_TEST_SUITE(AccelerationSerialTests, AccelerationSerialTestsFixture)
36
37#ifndef PRECICE_NO_MPI
38
39void testIQNIMVJPP(bool exchangeSubsteps)
40{
42 // use two vectors and see if underrelaxation works
43 double initialRelaxation = 0.01;
44 int maxIterationsUsed = 50;
45 int timeWindowsReused = 6;
46 int reusedTimeWindowsAtRestart = 0;
47 int chunkSize = 0;
48 int filter = Acceleration::QR1FILTER;
49 int restartType = IQNIMVJAcceleration::NO_RESTART;
50 double singularityLimit = 1e-10;
51 double svdTruncationEps = 0.0;
52 bool enforceInitialRelaxation = false;
53 bool alwaysBuildJacobian = false;
54 const double windowStart = 0;
55 const double windowEnd = 1;
56
57 std::vector<int> dataIDs;
58 dataIDs.push_back(0);
59 dataIDs.push_back(1);
60 std::vector<double> factors;
61 factors.resize(2, 1.0);
63 auto dummyMesh = testing::makeDummy2DMesh(4);
64
65 IQNIMVJAcceleration pp(initialRelaxation, enforceInitialRelaxation, maxIterationsUsed,
66 timeWindowsReused, filter, singularityLimit, dataIDs, prec, alwaysBuildJacobian,
67 restartType, chunkSize, reusedTimeWindowsAtRestart, svdTruncationEps, !exchangeSubsteps);
68
69 Eigen::VectorXd fcol1;
70
71 mesh::PtrData displacements(new mesh::Data("dvalues", -1, 1));
72 mesh::PtrData forces(new mesh::Data("fvalues", -1, 1));
73
74 // init displacements & forces
75 displacements->emplaceSampleAtTime(windowStart, {1.0, 1.0, 1.0, 1.0});
76 displacements->emplaceSampleAtTime(windowEnd, {1.0, 1.0, 1.0, 1.0});
77 forces->emplaceSampleAtTime(windowStart, {0.2, 0.2, 0.2, 0.2});
78 forces->emplaceSampleAtTime(windowEnd, {0.2, 0.2, 0.2, 0.2});
79
80 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
81 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
82
83 DataMap data;
86 dpcd->storeIteration();
87 fpcd->storeIteration();
88
89 pp.initialize(data);
90
91 displacements->emplaceSampleAtTime(windowEnd, {1.0, 2.0, 3.0, 4.0});
92 forces->emplaceSampleAtTime(windowEnd, {0.1, 0.1, 0.1, 0.1});
93
94 pp.performAcceleration(data, windowStart, windowEnd);
95
96 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(0), 1.00000000000000000000));
97 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(1), 1.01000000000000000888));
98 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(2), 1.02000000000000001776));
99 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(3), 1.03000000000000002665));
100 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(0), 0.199000000000000010214));
101 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(1), 0.199000000000000010214));
102 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(2), 0.199000000000000010214));
103 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(3), 0.199000000000000010214));
104
105 // Update the waveform as well
106 displacements->emplaceSampleAtTime(windowEnd, {10, 10, 10, 10});
107 forces->setSampleAtTime(windowEnd, forces->sample());
108
109 pp.performAcceleration(data, windowStart, windowEnd);
110
111 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(0), -5.63401340929695848558e-01));
112 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(1), 6.10309919173602111186e-01));
113 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(2), 1.78402117927690184729e+00));
114 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(3), 2.95773243938020247157e+00));
115 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(0), 8.28025852497733250157e-02));
116 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(1), 8.28025852497733250157e-02));
117 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(2), 8.28025852497733250157e-02));
118 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(3), 8.28025852497733250157e-02));
119}
120
121PRECICE_TEST_SETUP(1_rank)
122BOOST_AUTO_TEST_CASE(testIQNIMVJPPWithSubsteps)
123{
124 PRECICE_TEST();
125 testIQNIMVJPP(true);
126}
127
128PRECICE_TEST_SETUP(1_rank)
129BOOST_AUTO_TEST_CASE(testIQNIMVJPPWithoutSubsteps)
130{
131 PRECICE_TEST();
132 testIQNIMVJPP(false);
133}
134
135void testVIQNPP(bool exchangeSubsteps)
136{
138 // use two vectors and see if underrelaxation works
139
140 double initialRelaxation = 0.01;
141 int maxIterationsUsed = 50;
142 int timeWindowsReused = 6;
143 int filter = acceleration::BaseQNAcceleration::QR1FILTER;
144 double singularityLimit = 1e-10;
145 bool enforceInitialRelaxation = false;
146 std::vector<int> dataIDs;
147 const double windowStart = 0;
148 const double windowEnd = 1;
149
150 dataIDs.push_back(0);
151 dataIDs.push_back(1);
152 std::vector<double> factors;
153 factors.resize(2, 1.0);
155
156 std::map<int, double> scalings;
157 scalings.insert(std::make_pair(0, 1.0));
158 scalings.insert(std::make_pair(1, 1.0));
159 auto dummyMesh = testing::makeDummy2DMesh(4);
160
161 IQNILSAcceleration pp(initialRelaxation, enforceInitialRelaxation, maxIterationsUsed,
162 timeWindowsReused, filter, singularityLimit, dataIDs, prec, !exchangeSubsteps);
163
164 mesh::PtrData displacements(new mesh::Data("dvalues", -1, 1));
165 mesh::PtrData forces(new mesh::Data("fvalues", -1, 1));
166
167 // init displacements & forces
168 displacements->emplaceSampleAtTime(windowStart, {1.0, 1.0, 1.0, 1.0});
169 displacements->emplaceSampleAtTime(windowEnd, {1.0, 1.0, 1.0, 1.0});
170 forces->emplaceSampleAtTime(windowStart, {0.2, 0.2, 0.2, 0.2});
171 forces->emplaceSampleAtTime(windowEnd, {0.2, 0.2, 0.2, 0.2});
172
173 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
174 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
175
176 dpcd->storeIteration();
177 fpcd->storeIteration();
178
179 DataMap data;
180 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
181 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
182
183 pp.initialize(data);
184
185 displacements->emplaceSampleAtTime(windowEnd, {1.0, 2.0, 3.0, 4.0});
186 forces->emplaceSampleAtTime(windowEnd, {0.1, 0.1, 0.1, 0.1});
187
188 pp.performAcceleration(data, windowStart, windowEnd);
189
190 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(0), 1.00));
191 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(1), 1.01));
192 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(2), 1.02));
193 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(3), 1.03));
194 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(0), 0.199));
195 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(1), 0.199));
196 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(2), 0.199));
197 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(3), 0.199));
198
199 Eigen::VectorXd newdvalues;
200 utils::append(newdvalues, 10.0);
201 utils::append(newdvalues, 10.0);
202 utils::append(newdvalues, 10.0);
203 utils::append(newdvalues, 10.0);
204
205 displacements->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), newdvalues));
206 forces->setSampleAtTime(windowEnd, forces->sample());
207
208 pp.performAcceleration(data, windowStart, windowEnd);
209
210 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(0), -5.63401340929692295845e-01));
211 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(1), 6.10309919173607440257e-01));
212 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(2), 1.78402117927690717636e+00));
213 BOOST_TEST(testing::equals(data.at(0)->timeStepsStorage().sample(windowEnd)(3), 2.95773243938020513610e+00));
214 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(0), 8.28025852497733944046e-02));
215 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(1), 8.28025852497733944046e-02));
216 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(2), 8.28025852497733944046e-02));
217 BOOST_TEST(testing::equals(data.at(1)->timeStepsStorage().sample(windowEnd)(3), 8.28025852497733944046e-02));
218}
219
220PRECICE_TEST_SETUP(1_rank)
221BOOST_AUTO_TEST_CASE(testVIQNPPWithSubsteps)
222{
223 PRECICE_TEST();
224 testVIQNPP(true);
225}
226
227PRECICE_TEST_SETUP(1_rank)
228BOOST_AUTO_TEST_CASE(testVIQNPPWithoutSubsteps)
229{
230 PRECICE_TEST();
231 testVIQNPP(false);
232}
233
234PRECICE_TEST_SETUP(1_rank)
235BOOST_AUTO_TEST_CASE(testConstantUnderrelaxationWithSubsteps)
236{
237 PRECICE_TEST();
238 // use two vectors and see if underrelaxation works
239 double relaxation = 0.4;
240 std::vector<int> dataIDs{0, 1};
241 auto dummyMesh = testing::makeDummy3DMesh(4);
242 const double windowStart = 0;
243 const double windowEnd = 1;
244
245 ConstantRelaxationAcceleration acc(relaxation, dataIDs);
246
247 mesh::PtrData displacements = std::make_shared<mesh::Data>("dvalues", -1, 1);
248 mesh::PtrData forces = std::make_shared<mesh::Data>("fvalues", -1, 1);
249
250 // init displacements & forces
251 displacements->emplaceSampleAtTime(windowStart, {1.0, 2.0, 3.0, 4.0});
252 forces->emplaceSampleAtTime(windowStart, {0.2, 0.2, 0.2, 0.2});
253
254 bool exchangeSubsteps = false;
255
256 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
257 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
258
259 DataMap data;
260 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
261 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
262 dpcd->storeIteration();
263 fpcd->storeIteration();
264
265 acc.initialize(data);
266
267 displacements->emplaceSampleAtTime(windowEnd, {3.5, 2.0, 2.0, 1.0});
268 forces->emplaceSampleAtTime(windowEnd, {0.1, 0.1, 0.1, 0.1});
269
270 acc.performAcceleration(data, windowStart, windowEnd);
271
272 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 2);
273 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2);
274 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 2.6);
275 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 2.8);
276 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.16);
277 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.16);
278 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.16);
279 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.16);
280
281 displacements->emplaceSampleAtTime(windowEnd, {10, 10, 10, 10});
282 forces->setSampleAtTime(windowEnd, forces->sample());
283
284 acc.performAcceleration(data, windowStart, windowEnd);
285
286 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 4.6);
287 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 5.2);
288 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 5.8);
289 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 6.4);
290 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.184);
291 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.184);
292 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.184);
293 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.184);
294}
295
296PRECICE_TEST_SETUP(1_rank)
297BOOST_AUTO_TEST_CASE(testAitkenUnderrelaxationWithoutSubsteps)
298{
299 PRECICE_TEST();
300
301 double relaxation = 0.4;
302 std::vector<int> dataIDs{0, 1};
303 std::vector<double> factors{1, 1};
304 auto dummyMesh = testing::makeDummy3DMesh(4);
305 const double windowStart = 0;
306 const double windowEnd = 1;
307
309 AitkenAcceleration acc(relaxation, dataIDs, prec);
310
311 mesh::PtrData displacements = std::make_shared<mesh::Data>("dvalues", -1, 1);
312 mesh::PtrData forces = std::make_shared<mesh::Data>("fvalues", -1, 1);
313
314 // //init displacements & forces
315 displacements->emplaceSampleAtTime(windowStart, {1.0, 2.0, 3.0, 4.0});
316 forces->emplaceSampleAtTime(windowStart, {0.2, 0.2, 0.2, 0.2});
317
318 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, false);
319 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, false);
320
321 DataMap data;
322 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
323 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
324 dpcd->storeIteration();
325 fpcd->storeIteration();
326
327 acc.initialize(data);
328
329 displacements->emplaceSampleAtTime(windowEnd, {3.5, 2.0, 2.0, 1.0});
330 forces->emplaceSampleAtTime(windowEnd, {0.1, 0.1, 0.1, 0.1});
331
332 acc.performAcceleration(data, windowStart, windowEnd);
333
334 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 2);
335 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2);
336 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 2.6);
337 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 2.8);
338 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.16);
339 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.16);
340 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.16);
341 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.16);
342
343 displacements->emplaceSampleAtTime(windowEnd, {10, 10, 10, 10});
344 forces->setSampleAtTime(windowEnd, forces->sample());
345
346 acc.performAcceleration(data, windowStart, windowEnd);
347
348 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 1.2689851805508461);
349 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2.2390979382674185);
350 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 3.2092106959839914);
351 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 4.1793234537005644);
352 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.19880451030866292);
353 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.19880451030866292);
354 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.19880451030866292);
355 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.19880451030866292);
356}
357
358PRECICE_TEST_SETUP(1_rank)
359BOOST_AUTO_TEST_CASE(testAitkenUnderrelaxationWithPreconditioner)
360{
361 PRECICE_TEST();
362
363 double relaxation = 0.8;
364 std::vector<int> dataIDs{0, 1, 2, 3};
365 auto dummyMesh = testing::makeDummy2DMesh(2);
366
367 double windowStart = 0;
368 double windowEnd = 1;
369 const double dt = 1;
370
372 AitkenAcceleration acc(relaxation, dataIDs, prec);
373
374 mesh::PtrData data1 = std::make_shared<mesh::Data>("dvalues", -1, 1);
375 mesh::PtrData data2 = std::make_shared<mesh::Data>("fvalues", -1, 1);
376 mesh::PtrData data3 = std::make_shared<mesh::Data>("gvalues", -1, 2);
377 mesh::PtrData data4 = std::make_shared<mesh::Data>("hvalues", -1, 2);
378
379 // init data
380 data1->emplaceSampleAtTime(windowStart, {40, 80});
381 data2->emplaceSampleAtTime(windowStart, {5, 5});
382 data3->emplaceSampleAtTime(windowStart, {1, 2, 3, 4});
383 data4->emplaceSampleAtTime(windowStart, {20, 40, 60, 80});
384
385 cplscheme::PtrCouplingData dpcd = makeCouplingData(data1, dummyMesh, false);
386 cplscheme::PtrCouplingData fpcd = makeCouplingData(data2, dummyMesh, false);
387 cplscheme::PtrCouplingData gpcd = makeCouplingData(data3, dummyMesh, false);
388 cplscheme::PtrCouplingData hpcd = makeCouplingData(data4, dummyMesh, false);
389
390 DataMap data;
391 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
392 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
393 data.insert(std::pair<int, cplscheme::PtrCouplingData>(2, gpcd));
394 data.insert(std::pair<int, cplscheme::PtrCouplingData>(3, hpcd));
395 dpcd->storeIteration();
396 fpcd->storeIteration();
397 gpcd->storeIteration();
398 hpcd->storeIteration();
399
400 acc.initialize(data);
401
402 data1->emplaceSampleAtTime(windowEnd, {1, 7});
403 data2->emplaceSampleAtTime(windowEnd, {10, 10});
404 data3->emplaceSampleAtTime(windowEnd, {10, 11, 12, 13});
405 data4->emplaceSampleAtTime(windowEnd, {40, 60, 80, 100});
406
407 acc.performAcceleration(data, windowStart, windowEnd);
408
409 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 8.8);
410 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 21.6);
411 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 9);
412 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 9);
413 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(0) == 8.2);
414 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(1) == 9.2);
415 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(2) == 10.2);
416 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(0) == 36);
417 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(1) == 56);
418 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(2) == 76);
419 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(3) == 96);
420
421 data1->emplaceSampleAtTime(windowEnd, {2, 14});
422 data2->emplaceSampleAtTime(windowEnd, {8, 8});
423 data3->emplaceSampleAtTime(windowEnd, {13, 14, 15, 16});
424 data4->emplaceSampleAtTime(windowEnd, {41, 61, 81, 90});
425
426 acc.performAcceleration(data, windowStart, windowEnd);
427
428 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == -17.745640722103754);
429 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == -20.295060201548626);
430 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 9.5588663727976648);
431 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 9.5588663727976648);
432 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(0) == 19.235465491190659);
433 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(1) == 20.235465491190659);
434 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(2) == 21.235465491190659);
435 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(0) == 51.912064609583652);
436 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(1) == 71.912064609583652);
437 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(2) == 91.912064609583652);
438 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(3) == 95.196221242658879);
439
440 data1->emplaceSampleAtTime(windowEnd, {2.1, 14.1});
441 data2->emplaceSampleAtTime(windowEnd, {8, 8});
442 data3->emplaceSampleAtTime(windowEnd, {13.05, 14.07, 15.1, 16.1});
443 data4->emplaceSampleAtTime(windowEnd, {42, 60, 81.3, 91});
444
445 acc.iterationsConverged(data, windowStart);
446
447 // move to next window
448 windowStart += dt;
449 windowEnd += dt;
450
451 // move to next window
452 windowStart += dt;
453 windowEnd += dt;
454
455 data1->emplaceSampleAtTime(windowEnd, {3, 16});
456 data2->emplaceSampleAtTime(windowEnd, {7, 7});
457 data3->emplaceSampleAtTime(windowEnd, {18, 19, 20, 21});
458 data4->emplaceSampleAtTime(windowEnd, {50, 70, 90, 110});
459
460 acc.performAcceleration(data, windowStart, windowEnd);
461
462 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 10.4);
463 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 28.8);
464 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 6.6);
465 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 6.6);
466 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(0) == 14.6);
467 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(1) == 15.6);
468 BOOST_TEST(data.at(2)->timeStepsStorage().sample(windowEnd)(2) == 16.6);
469 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(0) == 44);
470 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(1) == 64);
471 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(2) == 84);
472 BOOST_TEST(data.at(3)->timeStepsStorage().sample(windowEnd)(3) == 104);
473}
474
475PRECICE_TEST_SETUP(1_rank)
476BOOST_AUTO_TEST_CASE(testConstantUnderrelaxationWithGradientWithSubsteps)
477{
478 PRECICE_TEST();
479 // use two vectors and see if underrelaxation works
480 double relaxation = 0.4;
481 std::vector<int> dataIDs{0, 1};
482 const int dim = 3;
483 auto dummyMesh = testing::makeDummy3DMesh(4);
484 const double windowStart = 0;
485 const double windowEnd = 1;
486
487 ConstantRelaxationAcceleration acc(relaxation, dataIDs);
488
489 mesh::PtrData displacements = std::make_shared<mesh::Data>("dvalues", -1, 1);
490 mesh::PtrData forces = std::make_shared<mesh::Data>("fvalues", -1, 1);
491
492 // init displacements
493 displacements->requireDataGradient();
494 Eigen::MatrixXd displacementGradient(displacements->gradients());
495 displacementGradient.resize(dim, 4);
496 for (unsigned int r = 0; r < dim; ++r) {
497 for (unsigned int c = 0; c < 4; ++c)
498 displacementGradient(r, c) = r + r * c;
499 }
500 displacements->setSampleAtTime(windowStart, time::Sample(displacements->getDimensions(), Eigen::Vector4d{1.0, 2.0, 3.0, 4.0}, displacementGradient));
501 // init forces
502 forces->requireDataGradient();
503 Eigen::MatrixXd forcesGradient(forces->gradients());
504 forcesGradient.resize(dim, 4);
505 forcesGradient.setConstant(-2);
506 forces->setSampleAtTime(windowStart, time::Sample(forces->getDimensions(), Eigen::Vector4d{0.2, 0.2, 0.2, 0.2}, forcesGradient));
507
508 bool exchangeSubsteps = true;
509
510 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
511 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
512
513 DataMap data;
514 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
515 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
516 dpcd->storeIteration();
517 fpcd->storeIteration();
518
519 acc.initialize(data);
520
521 displacements->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), Eigen::Vector4d{3.5, 2.0, 2.0, 1.0}, Eigen::MatrixXd(displacements->gradients()).setConstant(2.5)));
522 forces->setSampleAtTime(windowEnd, time::Sample(forces->getDimensions(), Eigen::Vector4d{0.1, 0.1, 0.1, 0.1}, Eigen::MatrixXd(forces->gradients()).setConstant(3)));
523
524 acc.performAcceleration(data, windowStart, windowEnd);
525
526 // Test value data
527 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 2);
528 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2);
529 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 2.6);
530 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 2.8);
531 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.16);
532 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.16);
533 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.16);
534 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.16);
535
536 // Test gradient data
537 BOOST_TEST(data.at(0)->gradients()(0, 0) == 1);
538 BOOST_TEST(data.at(0)->gradients()(0, 1) == 1);
539 BOOST_TEST(data.at(0)->gradients()(0, 2) == 1);
540 BOOST_TEST(data.at(0)->gradients()(1, 0) == 1.6);
541 BOOST_TEST(data.at(0)->gradients()(1, 1) == 2.2);
542 BOOST_TEST(data.at(0)->gradients()(1, 2) == 2.8);
543 BOOST_TEST(data.at(1)->gradients()(0, 0) == 0);
544 BOOST_TEST(data.at(1)->gradients()(0, 1) == 0);
545 BOOST_TEST(data.at(1)->gradients()(0, 2) == 0);
546 BOOST_TEST(data.at(1)->gradients()(1, 0) == 0);
547 BOOST_TEST(data.at(1)->gradients()(1, 1) == 0);
548 BOOST_TEST(data.at(1)->gradients()(1, 2) == 0);
549
550 displacements->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), Eigen::Vector4d{10, 10, 10, 10}, Eigen::MatrixXd(displacements->gradients()).setConstant(4)));
551 forces->setSampleAtTime(windowEnd, forces->sample());
552
553 acc.performAcceleration(data, windowStart, windowEnd);
554
555 // Check that store iteration works properly
556 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 4.6);
557 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 5.2);
558 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 5.8);
559 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 6.4);
560 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.184);
561 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.184);
562 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.184);
563 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.184);
564
565 BOOST_TEST(data.at(0)->gradients()(0, 0) == 1.6);
566 BOOST_TEST(data.at(0)->gradients()(0, 1) == 1.6);
567 BOOST_TEST(data.at(0)->gradients()(0, 2) == 1.6);
568 BOOST_TEST(data.at(0)->gradients()(1, 0) == 2.2);
569 BOOST_TEST(data.at(0)->gradients()(1, 1) == 2.8);
570 BOOST_TEST(data.at(0)->gradients()(1, 2) == 3.4);
571}
572
573PRECICE_TEST_SETUP(1_rank)
574BOOST_AUTO_TEST_CASE(testConstantUnderrelaxationWithoutSubsteps)
575{
576 PRECICE_TEST();
577 // use two vectors and see if underrelaxation works
578 double relaxation = 0.4;
579 std::vector<int> dataIDs{0, 1};
580 auto dummyMesh = testing::makeDummy3DMesh(4);
581 const double windowStart = 0;
582 const double windowEnd = 1;
583
584 ConstantRelaxationAcceleration acc(relaxation, dataIDs);
585
586 mesh::PtrData displacements = std::make_shared<mesh::Data>("dvalues", -1, 1);
587 mesh::PtrData forces = std::make_shared<mesh::Data>("fvalues", -1, 1);
588
589 displacements->emplaceSampleAtTime(windowStart, {1.0, 2.0, 3.0, 4.0});
590 forces->emplaceSampleAtTime(windowStart, {0.2, 0.2, 0.2, 0.2});
591
592 bool exchangeSubsteps = false;
593
594 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
595 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
596
597 DataMap data;
598 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
599 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
600 dpcd->storeIteration();
601 fpcd->storeIteration();
602
603 acc.initialize(data);
604
605 displacements->emplaceSampleAtTime(windowEnd, {3.5, 2.0, 2.0, 1.0});
606 forces->emplaceSampleAtTime(windowEnd, {0.1, 0.1, 0.1, 0.1});
607
608 acc.performAcceleration(data, windowStart, windowEnd);
609
610 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 2);
611 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2);
612 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 2.6);
613 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 2.8);
614 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.16);
615 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.16);
616 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.16);
617 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.16);
618
619 displacements->emplaceSampleAtTime(windowEnd, {10, 10, 10, 10});
620
621 forces->setSampleAtTime(windowEnd, forces->sample());
622
623 acc.performAcceleration(data, windowStart, windowEnd);
624
625 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 4.6);
626 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 5.2);
627 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 5.8);
628 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 6.4);
629 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.184);
630 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.184);
631 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.184);
632 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.184);
633}
634
635PRECICE_TEST_SETUP(1_rank)
636BOOST_AUTO_TEST_CASE(testConstantUnderrelaxationWithGradientWithoutSubsteps)
637{
638 PRECICE_TEST();
639 // use two vectors and see if underrelaxation works
640 double relaxation = 0.4;
641 std::vector<int> dataIDs{0, 1};
642 const int dim = 3;
643 auto dummyMesh = testing::makeDummy3DMesh(4);
644 const double windowStart = 0;
645 const double windowEnd = 1;
646
647 ConstantRelaxationAcceleration acc(relaxation, dataIDs);
648
649 mesh::PtrData displacements = std::make_shared<mesh::Data>("dvalues", -1, 1);
650 mesh::PtrData forces = std::make_shared<mesh::Data>("fvalues", -1, 1);
651
652 // init displacements
653 displacements->requireDataGradient();
654 Eigen::MatrixXd displacementGradient(displacements->gradients());
655 displacementGradient.resize(dim, 4);
656 for (unsigned int r = 0; r < dim; ++r) {
657 for (unsigned int c = 0; c < 4; ++c)
658 displacementGradient(r, c) = r + r * c;
659 }
660 displacements->setSampleAtTime(windowStart, time::Sample(displacements->getDimensions(), Eigen::Vector4d{1.0, 2.0, 3.0, 4.0}, displacementGradient));
661 // init forces
662 forces->requireDataGradient();
663 Eigen::MatrixXd forcesGradient(forces->gradients());
664 forcesGradient.resize(dim, 4);
665 forcesGradient.setConstant(-2);
666 forces->setSampleAtTime(windowStart, time::Sample(forces->getDimensions(), Eigen::Vector4d{0.2, 0.2, 0.2, 0.2}, forcesGradient));
667
668 bool exchangeSubsteps = false;
669
670 cplscheme::PtrCouplingData dpcd = makeCouplingData(displacements, dummyMesh, exchangeSubsteps);
671 cplscheme::PtrCouplingData fpcd = makeCouplingData(forces, dummyMesh, exchangeSubsteps);
672
673 DataMap data;
674 data.insert(std::pair<int, cplscheme::PtrCouplingData>(0, dpcd));
675 data.insert(std::pair<int, cplscheme::PtrCouplingData>(1, fpcd));
676 dpcd->storeIteration();
677 fpcd->storeIteration();
678
679 acc.initialize(data);
680
681 displacements->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), Eigen::Vector4d{3.5, 2.0, 2.0, 1.0}, Eigen::MatrixXd(displacements->gradients()).setConstant(2.5)));
682 forces->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), Eigen::Vector4d{0.1, 0.1, 0.1, 0.1}, Eigen::MatrixXd(displacements->gradients()).setConstant(3)));
683
684 acc.performAcceleration(data, windowStart, windowEnd);
685
686 // Test value data
687 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 2);
688 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 2);
689 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 2.6);
690 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 2.8);
691 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.16);
692 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.16);
693 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.16);
694 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.16);
695
696 // Test gradient data
697 BOOST_TEST(data.at(0)->gradients()(0, 0) == 1);
698 BOOST_TEST(data.at(0)->gradients()(0, 1) == 1);
699 BOOST_TEST(data.at(0)->gradients()(0, 2) == 1);
700 BOOST_TEST(data.at(0)->gradients()(1, 0) == 1.6);
701 BOOST_TEST(data.at(0)->gradients()(1, 1) == 2.2);
702 BOOST_TEST(data.at(0)->gradients()(1, 2) == 2.8);
703 BOOST_TEST(data.at(1)->gradients()(0, 0) == 0);
704 BOOST_TEST(data.at(1)->gradients()(0, 1) == 0);
705 BOOST_TEST(data.at(1)->gradients()(0, 2) == 0);
706 BOOST_TEST(data.at(1)->gradients()(1, 0) == 0);
707 BOOST_TEST(data.at(1)->gradients()(1, 1) == 0);
708 BOOST_TEST(data.at(1)->gradients()(1, 2) == 0);
709
710 displacements->setSampleAtTime(windowEnd, time::Sample(displacements->getDimensions(), Eigen::Vector4d{10, 10, 10, 10}, Eigen::MatrixXd(displacements->gradients()).setConstant(4)));
711 forces->setSampleAtTime(windowEnd, forces->sample());
712
713 acc.performAcceleration(data, windowStart, windowEnd);
714
715 // Check that store iteration works properly
716 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(0) == 4.6);
717 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(1) == 5.2);
718 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(2) == 5.8);
719 BOOST_TEST(data.at(0)->timeStepsStorage().sample(windowEnd)(3) == 6.4);
720 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(0) == 0.184);
721 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(1) == 0.184);
722 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(2) == 0.184);
723 BOOST_TEST(data.at(1)->timeStepsStorage().sample(windowEnd)(3) == 0.184);
724
725 BOOST_TEST(data.at(0)->gradients()(0, 0) == 1.6);
726 BOOST_TEST(data.at(0)->gradients()(0, 1) == 1.6);
727 BOOST_TEST(data.at(0)->gradients()(0, 2) == 1.6);
728 BOOST_TEST(data.at(0)->gradients()(1, 0) == 2.2);
729 BOOST_TEST(data.at(0)->gradients()(1, 1) == 2.8);
730 BOOST_TEST(data.at(0)->gradients()(1, 2) == 3.4);
731}
732
733#endif // not PRECICE_NO_MPI
734
BOOST_AUTO_TEST_CASE(testIQNIMVJPPWithSubsteps)
void testIQNIMVJPP(bool exchangeSubsteps)
void testVIQNPP(bool exchangeSubsteps)
BOOST_AUTO_TEST_SUITE(PreProcess)
BOOST_AUTO_TEST_SUITE_END()
#define PRECICE_TEST()
Definition Testing.hpp:39
#define PRECICE_TEST_SETUP(...)
Creates and attaches a TestSetup to a Boost test case.
Definition Testing.hpp:29
void performAcceleration(DataMap &cpldata, double windowStart, double windowEnd) final override
void initialize(const DataMap &cpldata) final override
void iterationsConverged(const DataMap &cpldata, double windowStart) final override
void initialize(const DataMap &cplData) final override
Initializes the acceleration.
void performAcceleration(DataMap &cplData, double windowStart, double windowEnd) final override
Performs one acceleration step.
void performAcceleration(DataMap &cplData, double windowStart, double windowEnd) override
Interface quasi-Newton with interface least-squares approximation.
Multi vector quasi-Newton update scheme.
Preconditioner that uses the constant user-defined factors to scale the quasi-Newton system.
Preconditioner that uses the recent residual to scale the quasi-Newton system.
Describes a set of data values belonging to the vertices of a mesh.
Definition Data.hpp:27
T insert(T... args)
T make_pair(T... args)
T make_unique(T... args)
contains implementations of acceleration schemes.
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:93
cplscheme::PtrCouplingData makeCouplingData(mesh::PtrData data, mesh::PtrMesh mesh, bool exchangeSubsteps)
Definition helper.hpp:10
auto makeDummy3DMesh(size_t nVertices)
Definition Meshes.hpp:28
auto makeDummy2DMesh(size_t nVertices)
Definition Meshes.hpp:21
void append(Eigen::VectorXd &v, double value)
Main namespace of the precice library.
T push_back(T... args)
T resize(T... args)
std::map< int, cplscheme::PtrCouplingData > DataMap