preCICE v3.1.2
Loading...
Searching...
No Matches
Data.cpp
Go to the documentation of this file.
1#include "Data.hpp"
2#include <algorithm>
3#include <utility>
4
8#include "utils/assertion.hpp"
9
10namespace precice::mesh {
11
14 DataID id,
15 int dimensions,
16 int spatialDimensions,
17 int waveformDegree)
18 : _waveform(waveformDegree),
19 _name(std::move(name)),
20 _id(id),
21 _dimensions(dimensions),
22 _spatialDimensions(spatialDimensions),
23 _sample(_dimensions)
24{
25 PRECICE_ASSERT(dimensions > 0, dimensions);
26}
27
28Eigen::VectorXd &Data::values()
29{
30 return _sample.values;
31}
32
33const Eigen::VectorXd &Data::values() const
34{
35 return _sample.values;
36}
37
38Eigen::MatrixXd &Data::gradients()
39{
40 return _sample.gradients;
41}
42
43const Eigen::MatrixXd &Data::gradients() const
44{
45 return _sample.gradients;
46}
47
49{
50 return _sample;
51}
52
54{
55 return _sample;
56}
57
58Eigen::VectorXd Data::sampleAtTime(double time) const
59{
60 return _waveform.sample(time);
61}
62
67
72
74{
75 if (stamples().size() > 1) { // Needed to avoid CompositionalCouplingScheme callong moveToNextWindow on same Data multiple times. Could be simplified by replacing Storage::move() with clearBefore(double time). See https://github.com/precice/precice/issues/1821.
77 PRECICE_ASSERT(stamples().size() == 1);
78 sample() = stamples().back().sample;
79 }
80}
81
82void Data::setSampleAtTime(double time, const time::Sample &sample)
83{
84 _sample = sample; // @todo at some point we should not need this anymore, when mapping, acceleration ... directly work on _timeStepsStorage
86}
87
89{
90 return _name;
91}
92
94{
95 return _id;
96}
97
99{
100 _sample.values.setZero();
101 if (_hasGradient) {
102 _sample.gradients.setZero();
103 }
104}
105
107{
108 return _hasGradient;
109}
110
112{
113 _hasGradient = true;
114};
115
117{
118 return _dimensions;
119}
120
121void Data::allocateValues(int expectedCount)
122{
123 using SizeType = std::remove_cv<decltype(expectedCount)>::type;
124 // Allocate data values
125 const SizeType expectedSize = expectedCount * _dimensions;
126 const auto actualSize = static_cast<SizeType>(_sample.values.size());
127 // Shrink Buffer
128 if (expectedSize < actualSize) {
129 _sample.values.resize(expectedSize);
130 }
131 // Enlarge Buffer
132 if (expectedSize > actualSize) {
133 const auto leftToAllocate = expectedSize - actualSize;
134 utils::append(_sample.values, Eigen::VectorXd(Eigen::VectorXd::Zero(leftToAllocate)));
135 }
136 PRECICE_DEBUG("Data {} now has {} values", _name, _sample.values.size());
137
138 // Allocate gradient data values
139 if (_hasGradient) {
140 const SizeType spaceDimensions = _spatialDimensions;
141
142 const SizeType expectedColumnSize = expectedCount * _dimensions;
143 const auto actualColumnSize = static_cast<SizeType>(_sample.gradients.cols());
144
145 // Shrink Buffer
146 if (expectedColumnSize < actualColumnSize) {
147 _sample.gradients.resize(spaceDimensions, expectedColumnSize);
148 }
149
150 // Enlarge Buffer
151 if (expectedColumnSize > actualColumnSize) {
152 const auto columnLeftToAllocate = expectedColumnSize - actualColumnSize;
153 utils::append(_sample.gradients, Eigen::MatrixXd(Eigen::MatrixXd::Zero(spaceDimensions, columnLeftToAllocate)));
154 }
155 PRECICE_DEBUG("Gradient Data {} now has {} x {} values", _name, _sample.gradients.rows(), _sample.gradients.cols());
156 }
157}
158
160{
161 return _spatialDimensions;
162}
163
164} // namespace precice::mesh
#define PRECICE_DEBUG(...)
Definition LogMacros.hpp:64
std::string name
#define PRECICE_ASSERT(...)
Definition assertion.hpp:87
time::Sample & sample()
Returns a reference to the _sample.
Definition Data.cpp:48
const std::string & getName() const
Returns the name of the data set, as set in the config file.
Definition Data.cpp:88
std::string _name
Name of the data set.
Definition Data.hpp:143
time::Storage & timeStepsStorage()
Returns a reference to the _timeStepsStorage of _waveform.
Definition Data.cpp:68
Data(std::string name, DataID id, int dimension, int spatialDimensions=-1, int waveformDegree=time::Time::DEFAULT_WAVEFORM_DEGREE)
Constructor.
Definition Data.cpp:12
int getDimensions() const
Returns the dimension (i.e., number of components) of one data value (i.e number of columns of one gr...
Definition Data.cpp:116
void requireDataGradient()
Set the additional requirement of gradient data.
Definition Data.cpp:111
bool hasGradient() const
Returns if the data contains gradient data.
Definition Data.cpp:106
int _spatialDimensions
Spatial Dimension of one element -> number of rows (only 2, 3 allowed for 2D, 3D).
Definition Data.hpp:152
auto stamples() const
Returns a the stamples from _timeStepsStorage.
Definition Data.hpp:100
void allocateValues(int expectedCount)
Allocates memory for the data values and corresponding gradient values.
Definition Data.cpp:121
int getSpatialDimensions() const
Returns the mesh dimension (i.e., number of rows) of one gradient data value .
Definition Data.cpp:159
int _dimensions
Dimensionality of one data value.
Definition Data.hpp:149
int getWaveformDegree() const
get degree of _waveform.
Definition Data.cpp:63
Eigen::VectorXd sampleAtTime(double time) const
Samples _waveform at given time.
Definition Data.cpp:58
void toZero()
Sets all values to zero.
Definition Data.cpp:98
Eigen::MatrixXd & gradients()
Returns a reference to the gradient data values.
Definition Data.cpp:38
DataID getID() const
Returns the ID of the data set (supposed to be unique).
Definition Data.cpp:93
DataID _id
ID of the data set (supposed to be unique).
Definition Data.hpp:146
void moveToNextWindow()
Definition Data.cpp:73
bool _hasGradient
Whether gradient data is available or not.
Definition Data.hpp:155
Eigen::VectorXd & values()
Returns a reference to the data values.
Definition Data.cpp:28
time::Waveform _waveform
Waveform wrapping this Data.
Definition Data.hpp:140
time::Sample _sample
Definition Data.hpp:157
void setSampleAtTime(double time, const time::Sample &sample)
Add sample at given time to _timeStepsStorage.
Definition Data.cpp:82
void setSampleAtTime(double time, const Sample &sample)
Store Sample at a specific time.
Definition Storage.cpp:27
void move()
Move this Storage by deleting all stamples except the one at the end of the window.
Definition Storage.cpp:86
int getInterpolationDegree() const
Definition Storage.cpp:61
Eigen::VectorXd sample(const double time) const
Evaluate waveform at specific point in time. Uses interpolation if necessary.
Definition Waveform.cpp:26
time::Storage & timeStepsStorage()
Returns a reference to the _timeStepsStorage.
Definition Waveform.cpp:16
provides Mesh, Data and primitives.
void append(Eigen::VectorXd &v, double value)
int DataID
Definition Types.hpp:25
STL namespace.
Eigen::MatrixXd gradients
The gradients of the data. Use gradients.col(i) to get the gradient at vertex i.
Definition Sample.hpp:52
Eigen::VectorXd values
Definition Sample.hpp:49