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

Prepare compilation time reduction

parent 07055986
No related branches found
No related tags found
1 merge request!43Feature/glace boundary conditions
Show whitespace changes
Inline Side-by-side
src/language/algorithms/AcousticSolverAlgorithm.hpp 0 → 100644
+ 222
− 0
View file @ a3c94db1
#ifndef ACOUSTIC_SOLVER_ALGORITHM_HPP
#define ACOUSTIC_SOLVER_ALGORITHM_HPP
#include <language/utils/FunctionSymbolId.hpp>
#include <language/utils/InterpolateItemValue.hpp>
#include <output/VTKWriter.hpp>
#include <scheme/AcousticSolver.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
#include <scheme/IBoundaryDescriptor.hpp>
#include <scheme/PressureBoundaryConditionDescriptor.hpp>
#include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
#include <scheme/VelocityBoundaryConditionDescriptor.hpp>
template <size_t Dimension, AcousticSolverType solver_type>
struct AcousticSolverAlgorithm
{
using ConnectivityType = Connectivity<Dimension>;
using MeshType = Mesh<ConnectivityType>;
using MeshDataType = MeshData<Dimension>;
using UnknownsType = FiniteVolumesEulerUnknowns<MeshType>;
AcousticSolverAlgorithm(std::shared_ptr<const IMesh> i_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list,
const FunctionSymbolId& rho_id,
const FunctionSymbolId& u_id,
const FunctionSymbolId& p_id)
{
std::shared_ptr<const MeshType> mesh = std::dynamic_pointer_cast<const MeshType>(i_mesh);
std::cout << "number of bc descr = " << bc_descriptor_list.size() << '\n';
typename AcousticSolver<MeshType>::BoundaryConditionList bc_list;
{
constexpr ItemType FaceType = [] {
if constexpr (Dimension > 1) {
return ItemType::face;
} else {
return ItemType::node;
}
}();
for (const auto& bc_descriptor : bc_descriptor_list) {
switch (bc_descriptor->type()) {
case IBoundaryConditionDescriptor::Type::symmetry: {
using SymmetryBoundaryCondition = typename AcousticSolver<MeshType>::SymmetryBoundaryCondition;
const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor =
dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == sym_bc_descriptor.boundaryDescriptor()) {
bc_list.push_back(
SymmetryBoundaryCondition{MeshFlatNodeBoundary<MeshType::Dimension>(mesh, ref_face_list)});
}
}
break;
}
case IBoundaryConditionDescriptor::Type::velocity: {
using VelocityBoundaryCondition = typename AcousticSolver<MeshType>::VelocityBoundaryCondition;
const VelocityBoundaryConditionDescriptor& velocity_bc_descriptor =
dynamic_cast<const VelocityBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == velocity_bc_descriptor.boundaryDescriptor()) {
MeshNodeBoundary<Dimension> mesh_node_boundary{mesh, ref_face_list};
const FunctionSymbolId velocity_id = velocity_bc_descriptor.functionSymbolId();
const auto& node_list = mesh_node_boundary.nodeList();
Array<const TinyVector<Dimension>> value_list = InterpolateItemValue<TinyVector<Dimension>(
TinyVector<Dimension>)>::template interpolate<ItemType::node>(velocity_id, mesh->xr(), node_list);
bc_list.push_back(VelocityBoundaryCondition{node_list, value_list});
}
}
break;
}
case IBoundaryConditionDescriptor::Type::pressure: {
using PressureBoundaryCondition = typename AcousticSolver<MeshType>::PressureBoundaryCondition;
const PressureBoundaryConditionDescriptor& pressure_bc_descriptor =
dynamic_cast<const PressureBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == pressure_bc_descriptor.boundaryDescriptor()) {
const auto& face_list = ref_face_list.list();
const FunctionSymbolId pressure_id = pressure_bc_descriptor.functionSymbolId();
Array<const double> face_values = [&] {
if constexpr (Dimension == 1) {
return InterpolateItemValue<double(
TinyVector<Dimension>)>::template interpolate<FaceType>(pressure_id, mesh->xr(), face_list);
} else {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
return InterpolateItemValue<double(
TinyVector<Dimension>)>::template interpolate<FaceType>(pressure_id, mesh_data.xl(), face_list);
}
}();
bc_list.push_back(PressureBoundaryCondition{face_list, face_values});
}
}
break;
}
default: {
std::ostringstream error_msg;
error_msg << *bc_descriptor << " is an invalid boundary condition for acoustic solver";
throw NormalError(error_msg.str());
}
}
}
}
UnknownsType unknowns(*mesh);
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
unknowns.rhoj() =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(rho_id,
mesh_data.xj());
unknowns.pj() =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(p_id, mesh_data.xj());
unknowns.uj() = InterpolateItemValue<TinyVector<Dimension>(
TinyVector<Dimension>)>::template interpolate<ItemType::cell>(u_id, mesh_data.xj());
}
unknowns.gammaj().fill(1.4);
AcousticSolver acoustic_solver(mesh, bc_list, solver_type);
const double tmax = 0.2;
double t = 0;
int itermax = std::numeric_limits<int>::max();
int iteration = 0;
CellValue<double>& rhoj = unknowns.rhoj();
CellValue<double>& ej = unknowns.ej();
CellValue<double>& pj = unknowns.pj();
CellValue<double>& gammaj = unknowns.gammaj();
CellValue<double>& cj = unknowns.cj();
CellValue<TinyVector<Dimension>>& uj = unknowns.uj();
CellValue<double>& Ej = unknowns.Ej();
CellValue<double>& mj = unknowns.mj();
CellValue<double>& inv_mj = unknowns.invMj();
BlockPerfectGas block_eos(rhoj, ej, pj, gammaj, cj);
block_eos.updateEandCFromRhoP();
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
const CellValue<const double> Vj = mesh_data.Vj();
parallel_for(
mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { Ej[j] = ej[j] + 0.5 * (uj[j], uj[j]); });
parallel_for(
mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { mj[j] = rhoj[j] * Vj[j]; });
parallel_for(
mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { inv_mj[j] = 1. / mj[j]; });
}
VTKWriter vtk_writer("mesh_" + std::to_string(Dimension), 0.01);
while ((t < tmax) and (iteration < itermax)) {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
vtk_writer.write(mesh,
{NamedItemValue{"density", rhoj}, NamedItemValue{"velocity", unknowns.uj()},
NamedItemValue{"coords", mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", mesh->connectivity().nodeOwner()}},
t);
double dt = 0.4 * acoustic_solver.acoustic_dt(mesh_data.Vj(), cj);
if (t + dt > tmax) {
dt = tmax - t;
}
std::cout.setf(std::cout.scientific);
std::cout << "iteration " << rang::fg::cyan << std::setw(4) << iteration << rang::style::reset
<< " time=" << rang::fg::green << t << rang::style::reset << " dt=" << rang::fgB::blue << dt
<< rang::style::reset << '\n';
mesh = acoustic_solver.computeNextStep(dt, unknowns);
block_eos.updatePandCFromRhoE();
t += dt;
++iteration;
}
std::cout << rang::style::bold << "Final time=" << rang::fgB::green << t << rang::style::reset << " reached after "
<< rang::fgB::cyan << iteration << rang::style::reset << rang::style::bold << " iterations"
<< rang::style::reset << '\n';
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*mesh);
vtk_writer.write(mesh,
{NamedItemValue{"density", rhoj}, NamedItemValue{"velocity", unknowns.uj()},
NamedItemValue{"coords", mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", mesh->connectivity().nodeOwner()}},
t, true); // forces last output
}
}
};
#endif // ACOUSTIC_SOLVER_ALGORITHM_HPP
src/language/modules/SchemeModule.cpp
+ 7
− 297
View file @ a3c94db1
#include <language/modules/SchemeModule.hpp> #include <language/modules/SchemeModule.hpp>
#include <language/algorithms/AcousticSolverAlgorithm.hpp>
#include <language/utils/BuiltinFunctionEmbedder.hpp> #include <language/utils/BuiltinFunctionEmbedder.hpp>
#include <language/utils/TypeDescriptor.hpp> #include <language/utils/TypeDescriptor.hpp>
#include <mesh/Mesh.hpp> #include <mesh/Mesh.hpp>
...@@ -8,300 +9,9 @@ ...@@ -8,300 +9,9 @@
#include <scheme/IBoundaryDescriptor.hpp> #include <scheme/IBoundaryDescriptor.hpp>
#include <scheme/NamedBoundaryDescriptor.hpp> #include <scheme/NamedBoundaryDescriptor.hpp>
#include <scheme/NumberedBoundaryDescriptor.hpp> #include <scheme/NumberedBoundaryDescriptor.hpp>
#include <scheme/PressureBoundaryConditionDescriptor.hpp>
#include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
#include <scheme/VelocityBoundaryConditionDescriptor.hpp>
#include <memory> #include <memory>
/////////// TEMPORARY
#include <language/utils/PugsFunctionAdapter.hpp>
#include <output/VTKWriter.hpp>
template <typename T>
class InterpolateItemValue;
template <typename OutputType, typename InputType>
class InterpolateItemValue<OutputType(InputType)> : public PugsFunctionAdapter<OutputType(InputType)>
{
static constexpr size_t Dimension = OutputType::Dimension;
using Adapter = PugsFunctionAdapter<OutputType(InputType)>;
public:
template <ItemType item_type>
static inline ItemValue<OutputType, item_type>
interpolate(const FunctionSymbolId& function_symbol_id, const ItemValue<const InputType, item_type>& position)
{
auto& expression = Adapter::getFunctionExpression(function_symbol_id);
auto convert_result = Adapter::getResultConverter(expression.m_data_type);
Array<ExecutionPolicy> context_list = Adapter::getContextList(expression);
using execution_space = typename Kokkos::DefaultExecutionSpace::execution_space;
Kokkos::Experimental::UniqueToken<execution_space, Kokkos::Experimental::UniqueTokenScope::Global> tokens;
const IConnectivity& connectivity = *position.connectivity_ptr();
ItemValue<OutputType, item_type> value(connectivity);
using ItemId = ItemIdT<item_type>;
parallel_for(connectivity.template numberOf<item_type>(), [=, &expression, &tokens](ItemId i) {
const int32_t t = tokens.acquire();
auto& execution_policy = context_list[t];
Adapter::convertArgs(execution_policy.currentContext(), position[i]);
auto result = expression.execute(execution_policy);
value[i] = convert_result(std::move(result));
tokens.release(t);
});
return value;
}
template <ItemType item_type>
static inline Array<OutputType>
interpolate(const FunctionSymbolId& function_symbol_id,
const ItemValue<const InputType, item_type>& position,
const Array<const ItemIdT<item_type>>& list_of_items)
{
auto& expression = Adapter::getFunctionExpression(function_symbol_id);
auto convert_result = Adapter::getResultConverter(expression.m_data_type);
Array<ExecutionPolicy> context_list = Adapter::getContextList(expression);
using execution_space = typename Kokkos::DefaultExecutionSpace::execution_space;
Kokkos::Experimental::UniqueToken<execution_space, Kokkos::Experimental::UniqueTokenScope::Global> tokens;
Array<OutputType> value{list_of_items.size()};
using ItemId = ItemIdT<item_type>;
parallel_for(list_of_items.size(), [=, &expression, &tokens](size_t i_item) {
ItemId item_id = list_of_items[i_item];
const int32_t t = tokens.acquire();
auto& execution_policy = context_list[t];
Adapter::convertArgs(execution_policy.currentContext(), position[item_id]);
auto result = expression.execute(execution_policy);
value[i_item] = convert_result(std::move(result));
tokens.release(t);
});
return value;
}
};
template <size_t Dimension, AcousticSolverType solver_type>
struct AcousticSolverScheme
{
using ConnectivityType = Connectivity<Dimension>;
using MeshType = Mesh<ConnectivityType>;
using MeshDataType = MeshData<Dimension>;
using UnknownsType = FiniteVolumesEulerUnknowns<MeshType>;
std::shared_ptr<const MeshType> m_mesh;
AcousticSolverScheme(std::shared_ptr<const IMesh> i_mesh,
const std::vector<std::shared_ptr<const IBoundaryConditionDescriptor>>& bc_descriptor_list,
const FunctionSymbolId& rho_id,
const FunctionSymbolId& u_id,
const FunctionSymbolId& p_id)
: m_mesh{std::dynamic_pointer_cast<const MeshType>(i_mesh)}
{
std::cout << "number of bc descr = " << bc_descriptor_list.size() << '\n';
typename AcousticSolver<MeshType>::BoundaryConditionList bc_list;
{
constexpr ItemType FaceType = [] {
if constexpr (Dimension > 1) {
return ItemType::face;
} else {
return ItemType::node;
}
}();
for (const auto& bc_descriptor : bc_descriptor_list) {
switch (bc_descriptor->type()) {
case IBoundaryConditionDescriptor::Type::symmetry: {
using SymmetryBoundaryCondition = typename AcousticSolver<MeshType>::SymmetryBoundaryCondition;
const SymmetryBoundaryConditionDescriptor& sym_bc_descriptor =
dynamic_cast<const SymmetryBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < m_mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = m_mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == sym_bc_descriptor.boundaryDescriptor()) {
bc_list.push_back(
SymmetryBoundaryCondition{MeshFlatNodeBoundary<MeshType::Dimension>(m_mesh, ref_face_list)});
}
}
break;
}
case IBoundaryConditionDescriptor::Type::velocity: {
using VelocityBoundaryCondition = typename AcousticSolver<MeshType>::VelocityBoundaryCondition;
const VelocityBoundaryConditionDescriptor& velocity_bc_descriptor =
dynamic_cast<const VelocityBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < m_mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = m_mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == velocity_bc_descriptor.boundaryDescriptor()) {
MeshNodeBoundary<Dimension> mesh_node_boundary{m_mesh, ref_face_list};
const FunctionSymbolId velocity_id = velocity_bc_descriptor.functionSymbolId();
const auto& node_list = mesh_node_boundary.nodeList();
Array<const TinyVector<Dimension>> value_list = InterpolateItemValue<TinyVector<Dimension>(
TinyVector<Dimension>)>::template interpolate<ItemType::node>(velocity_id, m_mesh->xr(), node_list);
bc_list.push_back(VelocityBoundaryCondition{node_list, value_list});
}
}
break;
}
case IBoundaryConditionDescriptor::Type::pressure: {
using PressureBoundaryCondition = typename AcousticSolver<MeshType>::PressureBoundaryCondition;
const PressureBoundaryConditionDescriptor& pressure_bc_descriptor =
dynamic_cast<const PressureBoundaryConditionDescriptor&>(*bc_descriptor);
for (size_t i_ref_face_list = 0;
i_ref_face_list < m_mesh->connectivity().template numberOfRefItemList<FaceType>(); ++i_ref_face_list) {
const auto& ref_face_list = m_mesh->connectivity().template refItemList<FaceType>(i_ref_face_list);
const RefId& ref = ref_face_list.refId();
if (ref == pressure_bc_descriptor.boundaryDescriptor()) {
const auto& face_list = ref_face_list.list();
const FunctionSymbolId pressure_id = pressure_bc_descriptor.functionSymbolId();
Array<const double> face_values = [&] {
if constexpr (Dimension == 1) {
return InterpolateItemValue<double(
TinyVector<Dimension>)>::template interpolate<FaceType>(pressure_id, m_mesh->xr(), face_list);
} else {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
return InterpolateItemValue<double(
TinyVector<Dimension>)>::template interpolate<FaceType>(pressure_id, mesh_data.xl(), face_list);
}
}();
bc_list.push_back(PressureBoundaryCondition{face_list, face_values});
}
}
break;
}
default: {
std::ostringstream error_msg;
error_msg << *bc_descriptor << " is an invalid boundary condition for acoustic solver";
throw NormalError(error_msg.str());
}
}
}
}
UnknownsType unknowns(*m_mesh);
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
unknowns.rhoj() =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(rho_id,
mesh_data.xj());
unknowns.pj() =
InterpolateItemValue<double(TinyVector<Dimension>)>::template interpolate<ItemType::cell>(p_id, mesh_data.xj());
unknowns.uj() = InterpolateItemValue<TinyVector<Dimension>(
TinyVector<Dimension>)>::template interpolate<ItemType::cell>(u_id, mesh_data.xj());
}
unknowns.gammaj().fill(1.4);
AcousticSolver acoustic_solver(m_mesh, bc_list, solver_type);
const double tmax = 0.2;
double t = 0;
int itermax = std::numeric_limits<int>::max();
int iteration = 0;
CellValue<double>& rhoj = unknowns.rhoj();
CellValue<double>& ej = unknowns.ej();
CellValue<double>& pj = unknowns.pj();
CellValue<double>& gammaj = unknowns.gammaj();
CellValue<double>& cj = unknowns.cj();
CellValue<TinyVector<Dimension>>& uj = unknowns.uj();
CellValue<double>& Ej = unknowns.Ej();
CellValue<double>& mj = unknowns.mj();
CellValue<double>& inv_mj = unknowns.invMj();
BlockPerfectGas block_eos(rhoj, ej, pj, gammaj, cj);
block_eos.updateEandCFromRhoP();
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
const CellValue<const double> Vj = mesh_data.Vj();
parallel_for(
m_mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { Ej[j] = ej[j] + 0.5 * (uj[j], uj[j]); });
parallel_for(
m_mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { mj[j] = rhoj[j] * Vj[j]; });
parallel_for(
m_mesh->numberOfCells(), PUGS_LAMBDA(CellId j) { inv_mj[j] = 1. / mj[j]; });
}
VTKWriter vtk_writer("mesh_" + std::to_string(Dimension), 0.01);
while ((t < tmax) and (iteration < itermax)) {
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
vtk_writer.write(m_mesh,
{NamedItemValue{"density", rhoj}, NamedItemValue{"velocity", unknowns.uj()},
NamedItemValue{"coords", m_mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", m_mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", m_mesh->connectivity().nodeOwner()}},
t);
double dt = 0.4 * acoustic_solver.acoustic_dt(mesh_data.Vj(), cj);
if (t + dt > tmax) {
dt = tmax - t;
}
std::cout.setf(std::cout.scientific);
std::cout << "iteration " << rang::fg::cyan << std::setw(4) << iteration << rang::style::reset
<< " time=" << rang::fg::green << t << rang::style::reset << " dt=" << rang::fgB::blue << dt
<< rang::style::reset << '\n';
m_mesh = acoustic_solver.computeNextStep(dt, unknowns);
block_eos.updatePandCFromRhoE();
t += dt;
++iteration;
}
std::cout << rang::style::bold << "Final time=" << rang::fgB::green << t << rang::style::reset << " reached after "
<< rang::fgB::cyan << iteration << rang::style::reset << rang::style::bold << " iterations"
<< rang::style::reset << '\n';
{
MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(*m_mesh);
vtk_writer.write(m_mesh,
{NamedItemValue{"density", rhoj}, NamedItemValue{"velocity", unknowns.uj()},
NamedItemValue{"coords", m_mesh->xr()}, NamedItemValue{"xj", mesh_data.xj()},
NamedItemValue{"cell_owner", m_mesh->connectivity().cellOwner()},
NamedItemValue{"node_owner", m_mesh->connectivity().nodeOwner()}},
t, true); // forces last output
}
}
};
SchemeModule::SchemeModule() SchemeModule::SchemeModule()
{ {
this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IBoundaryDescriptor>>); this->_addTypeDescriptor(ast_node_data_type_from<std::shared_ptr<const IBoundaryDescriptor>>);
...@@ -376,15 +86,15 @@ SchemeModule::SchemeModule() ...@@ -376,15 +86,15 @@ SchemeModule::SchemeModule()
const FunctionSymbolId& rho_id, const FunctionSymbolId& u_id, const FunctionSymbolId& p_id) -> void { const FunctionSymbolId& rho_id, const FunctionSymbolId& u_id, const FunctionSymbolId& p_id) -> void {
switch (p_mesh->dimension()) { switch (p_mesh->dimension()) {
case 1: { case 1: {
AcousticSolverScheme<1, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<1, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
case 2: { case 2: {
AcousticSolverScheme<2, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<2, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
case 3: { case 3: {
AcousticSolverScheme<3, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<3, AcousticSolverType::Glace>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
default: { default: {
...@@ -406,15 +116,15 @@ SchemeModule::SchemeModule() ...@@ -406,15 +116,15 @@ SchemeModule::SchemeModule()
const FunctionSymbolId& rho_id, const FunctionSymbolId& u_id, const FunctionSymbolId& p_id) -> void { const FunctionSymbolId& rho_id, const FunctionSymbolId& u_id, const FunctionSymbolId& p_id) -> void {
switch (p_mesh->dimension()) { switch (p_mesh->dimension()) {
case 1: { case 1: {
AcousticSolverScheme<1, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<1, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
case 2: { case 2: {
AcousticSolverScheme<2, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<2, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
case 3: { case 3: {
AcousticSolverScheme<3, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id}; AcousticSolverAlgorithm<3, AcousticSolverType::Eucclhyd>{p_mesh, bc_descriptor_list, rho_id, u_id, p_id};
break; break;
} }
default: { default: {
......
src/language/utils/InterpolateItemValue.hpp 0 → 100644
+ 82
− 0
View file @ a3c94db1
#ifndef INTERPOLATE_ITEM_VALUE_HPP
#define INTERPOLATE_ITEM_VALUE_HPP
#include <language/utils/PugsFunctionAdapter.hpp>
#include <mesh/ItemType.hpp>
#include <mesh/ItemValue.hpp>
template <typename T>
class InterpolateItemValue;
template <typename OutputType, typename InputType>
class InterpolateItemValue<OutputType(InputType)> : public PugsFunctionAdapter<OutputType(InputType)>
{
static constexpr size_t Dimension = OutputType::Dimension;
using Adapter = PugsFunctionAdapter<OutputType(InputType)>;
public:
template <ItemType item_type>
static inline ItemValue<OutputType, item_type>
interpolate(const FunctionSymbolId& function_symbol_id, const ItemValue<const InputType, item_type>& position)
{
auto& expression = Adapter::getFunctionExpression(function_symbol_id);
auto convert_result = Adapter::getResultConverter(expression.m_data_type);
Array<ExecutionPolicy> context_list = Adapter::getContextList(expression);
using execution_space = typename Kokkos::DefaultExecutionSpace::execution_space;
Kokkos::Experimental::UniqueToken<execution_space, Kokkos::Experimental::UniqueTokenScope::Global> tokens;
const IConnectivity& connectivity = *position.connectivity_ptr();
ItemValue<OutputType, item_type> value(connectivity);
using ItemId = ItemIdT<item_type>;
parallel_for(connectivity.template numberOf<item_type>(), [=, &expression, &tokens](ItemId i) {
const int32_t t = tokens.acquire();
auto& execution_policy = context_list[t];
Adapter::convertArgs(execution_policy.currentContext(), position[i]);
auto result = expression.execute(execution_policy);
value[i] = convert_result(std::move(result));
tokens.release(t);
});
return value;
}
template <ItemType item_type>
static inline Array<OutputType>
interpolate(const FunctionSymbolId& function_symbol_id,
const ItemValue<const InputType, item_type>& position,
const Array<const ItemIdT<item_type>>& list_of_items)
{
auto& expression = Adapter::getFunctionExpression(function_symbol_id);
auto convert_result = Adapter::getResultConverter(expression.m_data_type);
Array<ExecutionPolicy> context_list = Adapter::getContextList(expression);
using execution_space = typename Kokkos::DefaultExecutionSpace::execution_space;
Kokkos::Experimental::UniqueToken<execution_space, Kokkos::Experimental::UniqueTokenScope::Global> tokens;
Array<OutputType> value{list_of_items.size()};
using ItemId = ItemIdT<item_type>;
parallel_for(list_of_items.size(), [=, &expression, &tokens](size_t i_item) {
ItemId item_id = list_of_items[i_item];
const int32_t t = tokens.acquire();
auto& execution_policy = context_list[t];
Adapter::convertArgs(execution_policy.currentContext(), position[item_id]);
auto result = expression.execute(execution_policy);
value[i_item] = convert_result(std::move(result));
tokens.release(t);
});
return value;
}
};
#endif // INTERPOLATE_ITEM_VALUE_HPP
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment