Select Git revision
SchemeModule.cpp
-
clovis schoeck authored
Erreur oscillation fin de l'onde de détente Déplacement de l'onde de détente, choc et discontinuité de contact lrosque dt est modifié
clovis schoeck authoredErreur oscillation fin de l'onde de détente Déplacement de l'onde de détente, choc et discontinuité de contact lrosque dt est modifié
SchemeModule.cpp 46.80 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/IZoneDescriptor.hpp>
#include <mesh/Mesh.hpp>
#include <mesh/MeshData.hpp>
#include <mesh/MeshDataManager.hpp>
#include <mesh/MeshRandomizer.hpp>
#include <mesh/MeshSmoother.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/DiscreteFunctionVariant.hpp>
#include <scheme/DiscreteFunctionVectorIntegrator.hpp>
#include <scheme/DiscreteFunctionVectorInterpoler.hpp>
#include <scheme/ExternalBoundaryConditionDescriptor.hpp>
#include <scheme/FixedBoundaryConditionDescriptor.hpp>
#include <scheme/FluxingAdvectionSolver.hpp>
#include <scheme/FourierBoundaryConditionDescriptor.hpp>
#include <scheme/FreeBoundaryConditionDescriptor.hpp>
#include <scheme/HyperelasticSolver.hpp>
#include <scheme/HypoelasticSolver.hpp>
#include <scheme/IBoundaryConditionDescriptor.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 DiscreteFunctionVariant>>);
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::function(
[]() -> std::shared_ptr<const IDiscreteFunctionDescriptor> {
return std::make_shared<DiscreteFunctionDescriptorP0>();
}
));
this->_addBuiltinFunction("P0Vector", std::function(
[]() -> std::shared_ptr<const IDiscreteFunctionDescriptor> {
return std::make_shared<DiscreteFunctionDescriptorP0Vector>();
}
));
this->_addBuiltinFunction("Gauss", std::function(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussQuadratureDescriptor>(degree);
}
));
this->_addBuiltinFunction("GaussLobatto", std::function(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussLobattoQuadratureDescriptor>(degree);
}
));
this->_addBuiltinFunction("GaussLegendre", std::function(
[](uint64_t degree) -> std::shared_ptr<const IQuadratureDescriptor> {
return std::make_shared<GaussLegendreQuadratureDescriptor>(degree);
}
));
this->_addBuiltinFunction("integrate",
std::function(
[](std::shared_ptr<const IMesh> mesh,
const std::vector<std::shared_ptr<const IZoneDescriptor>>& integration_zone_list,
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 DiscreteFunctionVariant> {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionVectorIntegrator{mesh, integration_zone_list, quadrature_descriptor,
discrete_function_descriptor, function_id_list}
.integrate());
}
));
this->_addBuiltinFunction("integrate",
std::function(
[](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 DiscreteFunctionVariant> {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionVectorIntegrator{mesh, quadrature_descriptor,
discrete_function_descriptor, function_id_list}
.integrate());
}
));
this->_addBuiltinFunction("integrate",
std::function(
[](std::shared_ptr<const IMesh> mesh,
const std::vector<std::shared_ptr<const IZoneDescriptor>>& integration_zone_list,
std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
const FunctionSymbolId& function_id)
-> std::shared_ptr<const DiscreteFunctionVariant> {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionIntegrator{mesh, integration_zone_list, quadrature_descriptor,
function_id}
.integrate());
}
));
this->_addBuiltinFunction("integrate",
std::function(
[](std::shared_ptr<const IMesh> mesh,
std::shared_ptr<const IQuadratureDescriptor> quadrature_descriptor,
const FunctionSymbolId& function_id)
-> std::shared_ptr<const DiscreteFunctionVariant> {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionIntegrator{mesh, quadrature_descriptor, function_id}.integrate());
}
));
this->_addBuiltinFunction("interpolate",
std::function(
[](std::shared_ptr<const IMesh> mesh,
const std::vector<std::shared_ptr<const IZoneDescriptor>>& interpolation_zone_list,
std::shared_ptr<const IDiscreteFunctionDescriptor> discrete_function_descriptor,
const std::vector<FunctionSymbolId>& function_id_list)
-> std::shared_ptr<const DiscreteFunctionVariant> {
switch (discrete_function_descriptor->type()) {
case DiscreteFunctionType::P0: {
if (function_id_list.size() != 1) {
throw NormalError("invalid function descriptor type");
}
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionInterpoler{mesh, interpolation_zone_list,
discrete_function_descriptor, function_id_list[0]}
.interpolate());
}
case DiscreteFunctionType::P0Vector: {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionVectorInterpoler{mesh, interpolation_zone_list,
discrete_function_descriptor, function_id_list}
.interpolate());
}
default: {
throw NormalError("invalid function descriptor type");
}
}
}
));
this->_addBuiltinFunction("interpolate",
std::function(
[](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 DiscreteFunctionVariant> {
switch (discrete_function_descriptor->type()) {
case DiscreteFunctionType::P0: {
if (function_id_list.size() != 1) {
throw NormalError("invalid function descriptor type");
}
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionInterpoler{mesh, discrete_function_descriptor, function_id_list[0]}
.interpolate());
}
case DiscreteFunctionType::P0Vector: {
return std::make_shared<DiscreteFunctionVariant>(
DiscreteFunctionVectorInterpoler{mesh, discrete_function_descriptor,
function_id_list}
.interpolate());
}
default: {
throw NormalError("invalid function descriptor type");
}
}
}
));
this->_addBuiltinFunction("randomizeMesh",
std::function(
[](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::function(
[](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("smoothMesh", std::function(
[](std::shared_ptr<const IMesh> p_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list) -> std::shared_ptr<const IMesh> {
MeshSmootherHandler handler;
return handler.getSmoothedMesh(p_mesh, bc_descriptor_list);
}
));
this->_addBuiltinFunction("smoothMesh",
std::function(
[](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> {
MeshSmootherHandler handler;
return handler.getSmoothedMesh(p_mesh, bc_descriptor_list, function_symbol_id);
}
));
this->_addBuiltinFunction("smoothMesh",
std::function(
[](std::shared_ptr<const IMesh> p_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list,
const std::vector<std::shared_ptr<const IZoneDescriptor>>& smoothing_zone_list)
-> std::shared_ptr<const IMesh> {
MeshSmootherHandler handler;
return handler.getSmoothedMesh(p_mesh, bc_descriptor_list, smoothing_zone_list);
}
));
this->_addBuiltinFunction("smoothMesh", std::function(
[](std::shared_ptr<const IMesh> p_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>&
discrete_function_variant_list) -> std::shared_ptr<const IMesh> {
MeshSmootherHandler handler;
return handler.getSmoothedMesh(p_mesh, bc_descriptor_list,
discrete_function_variant_list);
}
));
this->_addBuiltinFunction("fixed", std::function(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<FixedBoundaryConditionDescriptor>(boundary);
}
));
this->_addBuiltinFunction("axis", std::function(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<AxisBoundaryConditionDescriptor>(boundary);
}
));
this->_addBuiltinFunction("symmetry", std::function(
[](std::shared_ptr<const IBoundaryDescriptor> boundary)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<SymmetryBoundaryConditionDescriptor>(boundary);
}
));
this->_addBuiltinFunction("pressure", std::function(
[](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("normalstress",
std::function(
[](std::shared_ptr<const IBoundaryDescriptor> boundary,
const FunctionSymbolId& normal_stress_id)
-> std::shared_ptr<const IBoundaryConditionDescriptor> {
return std::make_shared<DirichletBoundaryConditionDescriptor>("normal-stress", boundary,
normal_stress_id);
}
));
this->_addBuiltinFunction("velocity", std::function(
[](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::function(
[](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_fluxes", std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& c,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return AcousticSolverHandler{getCommonMesh({rho, c, u, p})}
.solver()
.compute_fluxes(AcousticSolverHandler::SolverType::Glace, rho, c, u,
p, bc_descriptor_list);
}
));
this->_addBuiltinFunction("glace_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& c,
const std::shared_ptr<const DiscreteFunctionVariant>& 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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return AcousticSolverHandler{getCommonMesh({rho, u, E, c, p})}
.solver()
.apply(AcousticSolverHandler::SolverType::Glace, dt, rho, u, E, c, p,
bc_descriptor_list);
}
));
this->_addBuiltinFunction("eucclhyd_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& c,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return AcousticSolverHandler{getCommonMesh({rho, c, u, p})}
.solver()
.compute_fluxes(AcousticSolverHandler::SolverType::Eucclhyd, rho, c, u, p,
bc_descriptor_list);
}
));
this->_addBuiltinFunction("eucclhyd_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& c,
const std::shared_ptr<const DiscreteFunctionVariant>& 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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return AcousticSolverHandler{getCommonMesh({rho, u, E, c, p})}
.solver()
.apply(AcousticSolverHandler::SolverType::Eucclhyd, dt, rho, u, E, c, p,
bc_descriptor_list);
}
));
this->_addBuiltinFunction("apply_acoustic_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho, //
const std::shared_ptr<const DiscreteFunctionVariant>& u, //
const std::shared_ptr<const DiscreteFunctionVariant>& E, //
const std::shared_ptr<const ItemValueVariant>& ur, //
const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, //
const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return AcousticSolverHandler{getCommonMesh({rho, u, E})} //
.solver()
.apply_fluxes(dt, rho, u, E, ur, Fjr);
}
));
this->_addBuiltinFunction("hyperelastic_eucclhyd_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return HyperelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigma})}
.solver()
.compute_fluxes(HyperelasticSolverHandler::SolverType::Eucclhyd, rho, aL, aT, u,
sigma, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hyperelastic_eucclhyd_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& CG,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG, aL, aT, sigma})}
.solver()
.apply(HyperelasticSolverHandler::SolverType::Eucclhyd, dt, rho, u, E, CG, aL, aT,
sigma, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hyperelastic_glace_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return HyperelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigma})}
.solver()
.compute_fluxes(HyperelasticSolverHandler::SolverType::Glace, //
rho, aL, aT, u, sigma, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hyperelastic_glace_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& CG,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& sigma,
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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG, aL, aT, sigma})}
.solver()
.apply(HyperelasticSolverHandler::SolverType::Glace, dt, rho, u, E, CG, aL, aT, sigma,
bc_descriptor_list);
}
));
this->_addBuiltinFunction("apply_hyperelastic_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho, //
const std::shared_ptr<const DiscreteFunctionVariant>& u, //
const std::shared_ptr<const DiscreteFunctionVariant>& E, //
const std::shared_ptr<const DiscreteFunctionVariant>& CG, //
const std::shared_ptr<const ItemValueVariant>& ur, //
const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, //
const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HyperelasticSolverHandler{getCommonMesh({rho, u, E, CG})} //
.solver()
.apply_fluxes(dt, rho, u, E, CG, ur, Fjr);
}
));
this->_addBuiltinFunction("hypoelastic_eucclhyd_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& sigmad,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return HypoelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigmad, p})}
.solver()
.compute_fluxes(HypoelasticSolverHandler::SolverType::Eucclhyd, rho, aL, aT, u,
sigmad, p, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hypoelastic_eucclhyd_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& sigmad,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::shared_ptr<const DiscreteFunctionVariant>& mu,
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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HypoelasticSolverHandler{getCommonMesh({rho, u, E, aL, aT, sigmad, p, mu})}
.solver()
.apply(HypoelasticSolverHandler::SolverType::Eucclhyd, dt, rho, u, E, sigmad, aL, aT,
p, mu, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hypoelastic_glace_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& sigmad,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>&
bc_descriptor_list)
-> std::tuple<std::shared_ptr<const ItemValueVariant>,
std::shared_ptr<const SubItemValuePerItemVariant>> {
return HypoelasticSolverHandler{getCommonMesh({rho, aL, aT, u, sigmad, p})}
.solver()
.compute_fluxes(HypoelasticSolverHandler::SolverType::Glace, //
rho, aL, aT, u, sigmad, p, bc_descriptor_list);
}
));
this->_addBuiltinFunction("hypoelastic_glace_solver",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho,
const std::shared_ptr<const DiscreteFunctionVariant>& u,
const std::shared_ptr<const DiscreteFunctionVariant>& E,
const std::shared_ptr<const DiscreteFunctionVariant>& aL,
const std::shared_ptr<const DiscreteFunctionVariant>& aT,
const std::shared_ptr<const DiscreteFunctionVariant>& sigmad,
const std::shared_ptr<const DiscreteFunctionVariant>& p,
const std::shared_ptr<const DiscreteFunctionVariant>& mu,
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 DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HypoelasticSolverHandler{getCommonMesh({rho, u, E, aL, aT, sigmad, mu, p})}
.solver()
.apply(HypoelasticSolverHandler::SolverType::Glace, dt, rho, u, E, sigmad, aL, aT, p,
mu, bc_descriptor_list);
}
));
this->_addBuiltinFunction("apply_hypoelastic_fluxes",
std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& rho, //
const std::shared_ptr<const DiscreteFunctionVariant>& u, //
const std::shared_ptr<const DiscreteFunctionVariant>& E, //
const std::shared_ptr<const DiscreteFunctionVariant>& sigmad, //
const std::shared_ptr<const DiscreteFunctionVariant>& mu, //
const std::shared_ptr<const ItemValueVariant>& ur, //
const std::shared_ptr<const SubItemValuePerItemVariant>& Fjr, //
const double& dt) -> std::tuple<std::shared_ptr<const IMesh>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>,
std::shared_ptr<const DiscreteFunctionVariant>> {
return HypoelasticSolverHandler{getCommonMesh({rho, u, E, sigmad, mu})} //
.solver()
.apply_fluxes(dt, rho, u, E, sigmad, mu, ur, Fjr);
}
));
this->_addBuiltinFunction("lagrangian",
std::function(
[](const std::shared_ptr<const IMesh>& mesh,
const std::shared_ptr<const DiscreteFunctionVariant>& v)
-> std::shared_ptr<const DiscreteFunctionVariant> { return shallowCopy(mesh, v); }
));
this->_addBuiltinFunction("acoustic_dt", std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
return acoustic_dt(c);
}
));
this->_addBuiltinFunction("cell_volume",
std::function(
[](const std::shared_ptr<const IMesh>& i_mesh)
-> std::shared_ptr<const DiscreteFunctionVariant> {
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<DiscreteFunctionVariant>(
DiscreteFunctionP0(mesh, 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<DiscreteFunctionVariant>(
DiscreteFunctionP0(mesh, 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<DiscreteFunctionVariant>(
DiscreteFunctionP0(mesh, MeshDataManager::instance().getMeshData(*mesh).Vj()));
}
default: {
throw UnexpectedError("invalid mesh dimension");
}
}
}
));
this->_addBuiltinFunction("hyperelastic_dt", std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
return hyperelastic_dt(c);
}
));
this->_addBuiltinFunction("hypoelastic_dt", std::function(
[](const std::shared_ptr<const DiscreteFunctionVariant>& c) -> double {
return hypoelastic_dt(c);
}
));
this->_addBuiltinFunction("fluxing_advection",
std::function(
[](const std::shared_ptr<const IMesh> new_mesh,
const std::vector<std::shared_ptr<const DiscreteFunctionVariant>>& remapped_variables)
-> std::vector<std::shared_ptr<const DiscreteFunctionVariant>> {
return advectByFluxing(new_mesh, remapped_variables);
}
));
MathFunctionRegisterForVh{*this};
}
void
SchemeModule::registerOperators() const
{
BinaryOperatorRegisterForVh{};
UnaryOperatorRegisterForVh{};
}