Select Git revision
SchemeModule.cpp
-
Axelle Drouard authoredAxelle Drouard authored
SchemeModule.cpp 23.24 KiB
#include <language/modules/SchemeModule.hpp>
#include <analysis/GaussLegendreQuadratureDescriptor.hpp>
#include <analysis/GaussLobattoQuadratureDescriptor.hpp>
#include <analysis/GaussQuadratureDescriptor.hpp>
#include <language/modules/BinaryOperatorRegisterForVh.hpp>
#include <language/modules/MathFunctionRegisterForVh.hpp>
#include <language/modules/UnaryOperatorRegisterForVh.hpp>
#include <language/utils/BinaryOperatorProcessorBuilder.hpp>
#include <language/utils/BuiltinFunctionEmbedder.hpp>
#include <language/utils/TypeDescriptor.hpp>
#include <mesh/Connectivity.hpp>
#include <mesh/IBoundaryDescriptor.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshData.hpp>
#include <mesh/MeshDataManager.hpp>
#include <mesh/MeshRandomizer.hpp>
#include <scheme/AcousticSolver.hpp>
#include <scheme/AxisBoundaryConditionDescriptor.hpp>
#include <scheme/DirichletBoundaryConditionDescriptor.hpp>
#include <scheme/DiscreteFunctionDescriptorP0.hpp>
#include <scheme/DiscreteFunctionDescriptorP0Vector.hpp>
#include <scheme/DiscreteFunctionIntegrator.hpp>
#include <scheme/DiscreteFunctionInterpoler.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <scheme/DiscreteFunctionUtils.hpp>
#include <scheme/DiscreteFunctionVectorIntegrator.hpp>
#include <scheme/DiscreteFunctionVectorInterpoler.hpp>
#include <scheme/ExternalBoundaryConditionDescriptor.hpp>
#include <scheme/FixedBoundaryConditionDescriptor.hpp>
#include <scheme/FourierBoundaryConditionDescriptor.hpp>
#include <scheme/FreeBoundaryConditionDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
#include <scheme/IDiscreteFunction.hpp>
#include <scheme/IDiscreteFunctionDescriptor.hpp>
#include <scheme/NeumannBoundaryConditionDescriptor.hpp>
#include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
#include <utils/Socket.hpp>
#include <memory>
SchemeModule::SchemeModule()
{
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IDiscreteFunction>>);
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IDiscreteFunctionDescriptor>>);
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IQuadratureDescriptor>>);
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IBoundaryConditionDescriptor>>);
this->_addBuiltinFunction("P0", std::make_shared<
BuiltinFunctionEmbedder<std::shared_ptr<const IDiscreteFunctionDescriptor>()>>(
[]() -> std::shared_ptr<const IDiscreteFunctionDescriptor> {
return std::make_shared<DiscreteFunctionDescriptorP0>();
}
));
this->_addBuiltinFunction("P0Vector",
std::make_shared<
BuiltinFunctionEmbedder<std::shared_ptr<const IDiscreteFunctionDescriptor>()>>(
[]() -> std::shared_ptr<const IDiscreteFunctionDescriptor> {
return std::make_shared<DiscreteFunctionDescriptorP0Vector>();
}
));
this->_addBuiltinFunction("Gauss", std::make_shared<
BuiltinFunctionEmbedder<std::shared_ptr<const IQuadratureDescriptor>(uint64_t)>>(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussQuadratureDescriptor>(degree);
}
));
this->_addBuiltinFunction("GaussLobatto",
std::make_shared<
BuiltinFunctionEmbedder<std::shared_ptr<const IQuadratureDescriptor>(uint64_t)>>(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussLobattoQuadratureDescriptor>(degree);
}
));
this->_addBuiltinFunction("GaussLegendre",
std::make_shared<
BuiltinFunctionEmbedder<std::shared_ptr<const IQuadratureDescriptor>(uint64_t)>>(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussLegendreQuadratureDescriptor>(degree);
}
));
this
->_addBuiltinFunction("integrate",
std::make_shared<BuiltinFunctionEmbedder<
std::shared_ptr<const IDiscreteFunction>(std::shared_ptr<const IMesh>,
std::shared_ptr<const IQuadratureDescriptor>,
std::shared_ptr<const IDiscreteFunctionDescriptor>,
const std::vector<FunctionSymbolId>&)>>(
[](std::shared_ptr<const IMesh> mesh,
std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
const std::vector<FunctionSymbolId>& function_id_list)
-> std::shared_ptr<const IDiscreteFunction> {
return DiscreteFunctionVectorIntegrator{mesh, quadrature_descriptor,
discrete_function_descriptor, function_id_list}
.integrate();
}
));
this->_addBuiltinFunction("integrate",
std::make_shared<BuiltinFunctionEmbedder<
std::shared_ptr<const IDiscreteFunction>(std::shared_ptr<const IMesh>,
std::shared_ptr<const IQuadratureDescriptor>,
const FunctionSymbolId&)>>(
[](std::shared_ptr<const IMesh> mesh,
std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
const FunctionSymbolId& function_id) -> std::shared_ptr<const IDiscreteFunction> {
return DiscreteFunctionIntegrator{mesh, quadrature_descriptor, function_id}.integrate();
}
));
this->_addBuiltinFunction(
"interpolate",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IDiscreteFunction>(std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunctionDescriptor>,
const std::vector<FunctionSymbolId>&)>>(
[](std::shared_ptr<const IMesh> mesh,
std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
const std::vector<FunctionSymbolId>& function_id_list) -> std::shared_ptr<const IDiscreteFunction> {
return DiscreteFunctionVectorInterpoler{mesh, discrete_function_descriptor, function_id_list}.interpolate();
}
));
this->_addBuiltinFunction(
"interpolate",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IDiscreteFunction>(std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunctionDescriptor>,
const FunctionSymbolId&)>>(
[](std::shared_ptr<const IMesh> mesh,
std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
const FunctionSymbolId& function_id) -> std::shared_ptr<const IDiscreteFunction> {
switch (discrete_function_descriptor->type()) {
case DiscreteFunctionType::P0: {
return DiscreteFunctionInterpoler{mesh, discrete_function_descriptor, function_id}.interpolate();
}
case DiscreteFunctionType::P0Vector: {
return DiscreteFunctionVectorInterpoler{mesh, discrete_function_descriptor, {function_id}}.interpolate();
}
default: {
throw NormalError("invalid function descriptor type");
}
}
}
));
this->_addBuiltinFunction("randomizeMesh",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IMesh>(std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&)>>(
[](std::shared_ptr<const IMesh> p_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list) -> std::shared_ptr<const IMesh> {
MeshRandomizerHandler handler;
return handler.getRandomizedMesh(*p_mesh, bc_descriptor_list);
}
));
this->_addBuiltinFunction("randomizeMesh",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IMesh>(std::shared_ptr<const IMesh>,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&,
const FunctionSymbolId&)>>(
[](std::shared_ptr<const IMesh> p_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list,
const FunctionSymbolId& function_symbol_id) -> std::shared_ptr<const IMesh> {
MeshRandomizerHandler handler;
return handler.getRandomizedMesh(*p_mesh, bc_descriptor_list, function_symbol_id);
}
));
this
->_addBuiltinFunction("fixed",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<const IBoundaryConditionDescriptor>(
std::shared_ptr<const IBoundaryDescriptor>)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<FixedBoundaryConditionDescriptor>(boundary);
}
));
this
->_addBuiltinFunction("axis",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<const IBoundaryConditionDescriptor>(
std::shared_ptr<const IBoundaryDescriptor>)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<AxisBoundaryConditionDescriptor>(boundary);
}
));
this
->_addBuiltinFunction("symmetry",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<const IBoundaryConditionDescriptor>(
std::shared_ptr<const IBoundaryDescriptor>)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<SymmetryBoundaryConditionDescriptor>(boundary);
}
));
this->_addBuiltinFunction("pressure",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IBoundaryConditionDescriptor>(std::shared_ptr<const IBoundaryDescriptor>,
const FunctionSymbolId&)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary,
const FunctionSymbolId& pressure_id)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<DirichletBoundaryConditionDescriptor>("pressure", boundary,
pressure_id);
}
));
this->_addBuiltinFunction("velocity",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IBoundaryConditionDescriptor>(std::shared_ptr<const IBoundaryDescriptor>,
const FunctionSymbolId&)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary,
const FunctionSymbolId& velocity_id)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<DirichletBoundaryConditionDescriptor>("velocity", boundary,
velocity_id);
}
));
this->_addBuiltinFunction("external_fsi_velocity",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<
const IBoundaryConditionDescriptor>(std::shared_ptr<const IBoundaryDescriptor>,
const std::shared_ptr<const Socket>&)>>(
[](std::shared_ptr<const IBoundaryDescriptor> boundary,
const std::shared_ptr<const Socket>& socket)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<ExternalBoundaryConditionDescriptor>("external_fsi_velocity",
boundary, socket);
}
));
this->_addBuiltinFunction("glace_solver",
std::make_shared<BuiltinFunctionEmbedder<std::tuple<
std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>>(const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::vector<std::shared_ptr<
const IBoundaryConditionDescriptor>>&,
const double&)>>(
[](const std::shared_ptr<const IDiscreteFunction>& rho,
const std::shared_ptr<const IDiscreteFunction>& u,
const std::shared_ptr<const IDiscreteFunction>& E,
const std::shared_ptr<const IDiscreteFunction>& c,
const std::shared_ptr<const IDiscreteFunction>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list,
const double& dt)
-> std::tuple<std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>> {
return AcousticSolverHandler{AcousticSolverHandler::SolverType::Glace,
rho,
c,
u,
p,
bc_descriptor_list}
.solver()
.apply(dt, rho, u, E);
}
));
this->_addBuiltinFunction("eucclhyd_solver",
std::make_shared<BuiltinFunctionEmbedder<std::tuple<
std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>>(const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::shared_ptr<const IDiscreteFunction>&,
const std::vector<std::shared_ptr<
const IBoundaryConditionDescriptor>>&,
const double&)>>(
[](const std::shared_ptr<const IDiscreteFunction>& rho,
const std::shared_ptr<const IDiscreteFunction>& u,
const std::shared_ptr<const IDiscreteFunction>& E,
const std::shared_ptr<const IDiscreteFunction>& c,
const std::shared_ptr<const IDiscreteFunction>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list,
const double& dt)
-> std::tuple<std::shared_ptr<const IMesh>, std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>,
std::shared_ptr<const IDiscreteFunction>> {
return AcousticSolverHandler{AcousticSolverHandler::SolverType::Eucclhyd,
rho,
c,
u,
p,
bc_descriptor_list}
.solver()
.apply(dt, rho, u, E);
}
));
this
->_addBuiltinFunction("lagrangian",
std::make_shared<BuiltinFunctionEmbedder<
std::shared_ptr<const IDiscreteFunction>(const std::shared_ptr<const IMesh>&,
const std::shared_ptr<const IDiscreteFunction>&)>>(
[](const std::shared_ptr<const IMesh>& mesh,
const std::shared_ptr<const IDiscreteFunction>& v)
-> std::shared_ptr<const IDiscreteFunction> { return shallowCopy(mesh, v); }
));
this->_addBuiltinFunction("acoustic_dt",
std::make_shared<
BuiltinFunctionEmbedder<double(const std::shared_ptr<const IDiscreteFunction>&)>>(
[](const std::shared_ptr<const IDiscreteFunction>& c) -> double { return acoustic_dt(c); }
));
this
->_addBuiltinFunction("cell_volume",
std::make_shared<BuiltinFunctionEmbedder<std::shared_ptr<const IDiscreteFunction>(
const std::shared_ptr<const IMesh>&)>>(
[](const std::shared_ptr<const IMesh>& i_mesh) -> std::shared_ptr<const IDiscreteFunction> {
switch (i_mesh->dimension()) {
case 1: {
constexpr size_t Dimension = 1;
using MeshType = Mesh<Connectivity<Dimension>>;
std::shared_ptr<const MeshType> mesh =
std::dynamic_pointer_cast<const Mesh<Connectivity<Dimension>>>(i_mesh);
return std::make_shared<const DiscreteFunctionP0<
Dimension, double>>(mesh, copy(MeshDataManager::instance().getMeshData(*mesh).Vj()));
}
case 2: {
constexpr size_t Dimension = 2;
using MeshType = Mesh<Connectivity<Dimension>>;
std::shared_ptr<const MeshType> mesh =
std::dynamic_pointer_cast<const Mesh<Connectivity<Dimension>>>(i_mesh);
return std::make_shared<const DiscreteFunctionP0<
Dimension, double>>(mesh, copy(MeshDataManager::instance().getMeshData(*mesh).Vj()));
}
case 3: {
constexpr size_t Dimension = 3;
using MeshType = Mesh<Connectivity<Dimension>>;
std::shared_ptr<const MeshType> mesh =
std::dynamic_pointer_cast<const Mesh<Connectivity<Dimension>>>(i_mesh);
return std::make_shared<const DiscreteFunctionP0<
Dimension, double>>(mesh, copy(MeshDataManager::instance().getMeshData(*mesh).Vj()));
}
default: {
throw UnexpectedError("invalid mesh dimension");
}
}
}
));
MathFunctionRegisterForVh{*this};
}
void
SchemeModule::registerOperators() const
{
BinaryOperatorRegisterForVh{};
UnaryOperatorRegisterForVh{};
}