Skip to content
Snippets Groups Projects
Commit bc654ee5 authored by Stéphane Del Pino's avatar Stéphane Del Pino
Browse files

Add a mesh interpolation functionality

parent 9b310bb4
No related branches found
No related tags found
1 merge request!119Feature/mesh interpoler
...@@ -12,6 +12,7 @@ ...@@ -12,6 +12,7 @@
#include <mesh/DiamondDualMeshManager.hpp> #include <mesh/DiamondDualMeshManager.hpp>
#include <mesh/GmshReader.hpp> #include <mesh/GmshReader.hpp>
#include <mesh/Mesh.hpp> #include <mesh/Mesh.hpp>
#include <mesh/MeshInterpoler.hpp>
#include <utils/Exceptions.hpp> #include <utils/Exceptions.hpp>
#include <Kokkos_Core.hpp> #include <Kokkos_Core.hpp>
...@@ -99,6 +100,19 @@ MeshModule::MeshModule() ...@@ -99,6 +100,19 @@ MeshModule::MeshModule()
)); ));
this->_addBuiltinFunction("interpolate",
std::make_shared<BuiltinFunctionEmbedder<
std::shared_ptr<const IMesh>(const std::shared_ptr<const IMesh>&,
const std::shared_ptr<const IMesh>&, const double&)>>(
[](const std::shared_ptr<const IMesh>& source_mesh,
const std::shared_ptr<const IMesh>& destination_mesh,
const double& theta) -> std::shared_ptr<const IMesh> {
return MeshInterpoler{}.interpolate(source_mesh, destination_mesh, theta);
}
));
this->_addBuiltinFunction("cartesian1dMesh", this->_addBuiltinFunction("cartesian1dMesh",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr< std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IMesh>(const TinyVector<1>, const TinyVector<1>, const std::vector<uint64_t>&)>>( const IMesh>(const TinyVector<1>, const TinyVector<1>, const std::vector<uint64_t>&)>>(
......
...@@ -20,6 +20,7 @@ add_library( ...@@ -20,6 +20,7 @@ add_library(
MeshFaceBoundary.cpp MeshFaceBoundary.cpp
MeshFlatFaceBoundary.cpp MeshFlatFaceBoundary.cpp
MeshFlatNodeBoundary.cpp MeshFlatNodeBoundary.cpp
MeshInterpoler.cpp
MeshLineNodeBoundary.cpp MeshLineNodeBoundary.cpp
MeshNodeBoundary.cpp MeshNodeBoundary.cpp
MeshRandomizer.cpp MeshRandomizer.cpp
......
#include <mesh/MeshInterpoler.hpp>
#include <mesh/Connectivity.hpp>
#include <mesh/Mesh.hpp>
template <typename ConnectivityType>
std::shared_ptr<const Mesh<ConnectivityType>>
MeshInterpoler::_interpolate(const Mesh<ConnectivityType>& source_mesh,
const Mesh<ConnectivityType>& destination_mesh,
const double& theta) const
{
if (source_mesh.shared_connectivity() == destination_mesh.shared_connectivity()) {
const ConnectivityType& connectivity = source_mesh.connectivity();
NodeValue<TinyVector<ConnectivityType::Dimension>> theta_xr{connectivity};
const NodeValue<const TinyVector<ConnectivityType::Dimension>> source_xr = source_mesh.xr();
const NodeValue<const TinyVector<ConnectivityType::Dimension>> destination_xr = destination_mesh.xr();
const double one_minus_theta = 1 - theta;
parallel_for(
connectivity.numberOfNodes(), PUGS_LAMBDA(const NodeId node_id) {
theta_xr[node_id] = one_minus_theta * source_xr[node_id] + theta * destination_xr[node_id];
});
return std::make_shared<Mesh<ConnectivityType>>(source_mesh.shared_connectivity(), theta_xr);
} else {
throw NormalError("interpolated meshes must share the same connectivity");
}
}
std::shared_ptr<const IMesh>
MeshInterpoler::interpolate(const std::shared_ptr<const IMesh>& p_source_mesh,
const std::shared_ptr<const IMesh>& p_destination_mesh,
const double& theta) const
{
if (p_source_mesh->dimension() != p_destination_mesh->dimension()) {
throw NormalError("incompatible mesh dimensions");
} else {
switch (p_source_mesh->dimension()) {
case 1: {
using MeshType = Mesh<Connectivity<1>>;
return this->_interpolate(dynamic_cast<const MeshType&>(*p_source_mesh),
dynamic_cast<const MeshType&>(*p_destination_mesh), theta);
}
case 2: {
using MeshType = Mesh<Connectivity<2>>;
return this->_interpolate(dynamic_cast<const MeshType&>(*p_source_mesh),
dynamic_cast<const MeshType&>(*p_destination_mesh), theta);
}
case 3: {
using MeshType = Mesh<Connectivity<3>>;
return this->_interpolate(dynamic_cast<const MeshType&>(*p_source_mesh),
dynamic_cast<const MeshType&>(*p_destination_mesh), theta);
}
default: {
throw UnexpectedError("invalid mesh dimension");
}
}
}
}
#ifndef MESH_INTERPOLER_HPP
#define MESH_INTERPOLER_HPP
class IMesh;
template <typename ConnectivityType>
class Mesh;
#include <memory>
class MeshInterpoler
{
private:
template <typename ConnectivityType>
std::shared_ptr<const Mesh<ConnectivityType>> _interpolate(const Mesh<ConnectivityType>& source_mesh,
const Mesh<ConnectivityType>& destination_mesh,
const double& theta) const;
public:
std::shared_ptr<const IMesh> interpolate(const std::shared_ptr<const IMesh>& p_source_mesh,
const std::shared_ptr<const IMesh>& p_destination_mesh,
const double& theta) const;
MeshInterpoler() = default;
~MeshInterpoler() = default;
};
#endif // MESH_INTERPOLER_HPP
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment