5#include <boost/container/flat_set.hpp>
34template <
typename RADIAL_BASIS_FUNCTION_T>
57 RADIAL_BASIS_FUNCTION_T function,
127template <
typename RADIAL_BASIS_FUNCTION_T>
131 RADIAL_BASIS_FUNCTION_T function,
135 : _center(center), _radius(radius), _polynomial(polynomial), _weightingFunction(radius)
146 auto inIDs = inputMesh->index().getVerticesInsideBox(center, radius);
150 _inputIDs.insert(inIDs.begin(), inIDs.end());
151 _outputIDs.insert(outIDs.begin(), outIDs.end());
158 PRECICE_DEBUG(
"SphericalVertexCluster input size: {}", inIDs.size());
159 PRECICE_DEBUG(
"SphericalVertexCluster output size: {}", outIDs.size());
173template <
typename RADIAL_BASIS_FUNCTION_T>
180 if (_normalizedWeights.size() == 0)
181 _normalizedWeights.resize(_outputIDs.size());
184 auto localID = _outputIDs.index_of(_outputIDs.find(
id));
186 PRECICE_ASSERT(
static_cast<Eigen::Index
>(localID) < _normalizedWeights.size(), localID, _normalizedWeights.size());
187 _normalizedWeights[localID] = normalizedWeight;
190template <
typename RADIAL_BASIS_FUNCTION_T>
196 PRECICE_ASSERT(_normalizedWeights.size() ==
static_cast<Eigen::Index
>(_outputIDs.size()));
199 const unsigned int nComponents = inData.
dataDims;
200 const auto & localInData = inData.
values;
203 Eigen::VectorXd in(_rbfSolver.getOutputSize());
206 for (
unsigned int c = 0; c < nComponents; ++c) {
209 for (
unsigned int i = 0; i < _outputIDs.size(); ++i) {
210 const auto dataIndex = *(_outputIDs.nth(i));
211 PRECICE_ASSERT(dataIndex * nComponents + c < localInData.size(), dataIndex * nComponents + c, localInData.size());
212 PRECICE_ASSERT(_normalizedWeights[i] > 0, _normalizedWeights[i], i);
214 in[i] = localInData[dataIndex * nComponents + c] * _normalizedWeights[i];
218 auto result = _rbfSolver.solveConservative(in, _polynomial);
219 PRECICE_ASSERT(result.size() ==
static_cast<Eigen::Index
>(_inputIDs.size()));
222 for (
unsigned int i = 0; i < _inputIDs.size(); ++i) {
223 const auto dataIndex = *(_inputIDs.nth(i));
224 PRECICE_ASSERT(dataIndex * nComponents + c < outData.size(), dataIndex * nComponents + c, outData.size());
225 outData[dataIndex * nComponents + c] += result(i);
230template <
typename RADIAL_BASIS_FUNCTION_T>
236 PRECICE_ASSERT(_normalizedWeights.size() ==
static_cast<Eigen::Index
>(_outputIDs.size()));
239 const unsigned int nComponents = inData.
dataDims;
240 const auto & localInData = inData.
values;
242 Eigen::VectorXd in(_rbfSolver.getInputSize());
245 for (
unsigned int c = 0; c < nComponents; ++c) {
248 for (
unsigned int i = 0; i < _inputIDs.size(); i++) {
249 const auto dataIndex = *(_inputIDs.nth(i));
250 PRECICE_ASSERT(dataIndex * nComponents + c < localInData.size(), dataIndex * nComponents + c, localInData.size());
251 in[i] = localInData[dataIndex * nComponents + c];
255 auto result = _rbfSolver.solveConsistent(in, _polynomial);
256 PRECICE_ASSERT(
static_cast<Eigen::Index
>(_outputIDs.size()) == result.size());
259 for (
unsigned int i = 0; i < _outputIDs.size(); ++i) {
260 const auto dataIndex = *(_outputIDs.nth(i));
261 PRECICE_ASSERT(dataIndex * nComponents + c < outData.size(), dataIndex * nComponents + c, outData.size());
264 outData[dataIndex * nComponents + c] += result(i) * _normalizedWeights[i];
269template <
typename RADIAL_BASIS_FUNCTION_T>
275 return _weightingFunction.evaluate(
std::sqrt(res));
278template <
typename RADIAL_BASIS_FUNCTION_T>
281 return _inputIDs.size();
284template <
typename RADIAL_BASIS_FUNCTION_T>
287 return _center.rawCoords();
290template <
typename RADIAL_BASIS_FUNCTION_T>
293 return _outputIDs.size() == 0 || _inputIDs.size() == 0;
296template <
typename RADIAL_BASIS_FUNCTION_T>
303 _hasComputedMapping =
false;
#define PRECICE_DEBUG(...)
#define PRECICE_TRACE(...)
#define PRECICE_ASSERT(...)
This class provides a lightweight logger.
Wendland radial basis function with compact support.
boost::container::flat_set< VertexID > _outputIDs
Global VertexIDs associated to the output mesh (see constructor)
unsigned int getNumberOfInputVertices() const
Number of input vertices this partition operates on.
Eigen::VectorXd _normalizedWeights
mesh::Vertex _center
center vertex of the cluster
const double _radius
radius of the vertex cluster
std::array< double, 3 > getCenterCoords() const
The center coordinate of this cluster.
double computeWeight(const mesh::Vertex &v) const
Compute the weight for a given vertex.
bool _hasComputedMapping
Boolean switch in order to indicate that a mapping was computed.
CompactPolynomialC2 _weightingFunction
The weighting function.
void setNormalizedWeight(double normalizedWeight, VertexID vertexID)
Set the normalized weight for the given vertexID in the outputMesh.
Polynomial _polynomial
Polynomial treatment in the RBF solver.
RadialBasisFctSolver< RADIAL_BASIS_FUNCTION_T > _rbfSolver
The RBF solver.
void mapConsistent(const time::Sample &inData, Eigen::VectorXd &outData) const
Evaluates a consistent mapping and agglomerates the result in the given output data.
void clear()
Invalidates and erases data structures the cluster holds.
SphericalVertexCluster(mesh::Vertex center, double radius, RADIAL_BASIS_FUNCTION_T function, Polynomial polynomial, mesh::PtrMesh inputMesh, mesh::PtrMesh outputMesh)
boost::container::flat_set< VertexID > _inputIDs
Global VertexIDs associated to the input mesh (see constructor)
void mapConservative(const time::Sample &inData, Eigen::VectorXd &outData) const
Evaluates a conservative mapping and agglomerates the result in the given output data.
precice::logging::Logger _log
logger, as usual
const RawCoords & rawCoords() const
Direct access to the coordinates.
Eigen::VectorXd getCoords() const
Returns the coordinates of the vertex.
void stop()
Stops a running event.
double computeSquaredDifference(const std::array< double, 3 > &u, std::array< double, 3 > v, const std::array< bool, 3 > &activeAxis={{true, true, true}})
Deletes all dead directions from fullVector and returns a vector of reduced dimensionality.
Polynomial
How to handle the polynomial?
constexpr double NUMERICAL_ZERO_DIFFERENCE
Main namespace of the precice library.
int dataDims
The dimensionality of the data.