8#include <compile_time_options.h>
12#include <deal.II/base/function_parser.h>
16 namespace EulerInitialStates
23 template <
typename Description,
int dim,
typename Number>
33 const std::string subsection)
35 , hyperbolic_system_(hyperbolic_system)
38 density_expression_ =
"1.4";
39 this->add_parameter(
"density expression",
41 "A function expression describing the density");
43 velocity_x_expression_ =
"3.0";
45 "velocity x expression",
46 velocity_x_expression_,
47 "A function expression describing the x-component of the velocity");
49 if constexpr (dim > 1) {
50 velocity_y_expression_ =
"0.0";
51 this->add_parameter(
"velocity y expression",
52 velocity_y_expression_,
53 "A function expression describing the "
54 "y-component of the velocity");
57 if constexpr (dim > 2) {
58 velocity_z_expression_ =
"0.0";
59 this->add_parameter(
"velocity z expression",
60 velocity_z_expression_,
61 "A function expression describing the "
62 "z-component of the velocity");
65 pressure_expression_ =
"1.0";
66 this->add_parameter(
"pressure expression",
68 "A function expression describing the pressure");
74 const auto set_up_muparser = [
this] {
75 using FP = dealii::FunctionParser<dim>;
81 density_function_ = std::make_unique<FP>(density_expression_);
82 velocity_x_function_ = std::make_unique<FP>(velocity_x_expression_);
83 if constexpr (dim > 1)
84 velocity_y_function_ = std::make_unique<FP>(velocity_y_expression_);
85 if constexpr (dim > 2)
86 velocity_z_function_ = std::make_unique<FP>(velocity_z_expression_);
87 pressure_function_ = std::make_unique<FP>(pressure_expression_);
91 this->parse_parameters_call_back.connect(set_up_muparser);
96 const auto view = hyperbolic_system_.template view<dim, Number>();
99 density_function_->set_time(t);
100 full_primitive_state[0] = density_function_->value(point);
102 velocity_x_function_->set_time(t);
103 full_primitive_state[1] = velocity_x_function_->value(point);
105 if constexpr (dim > 1) {
106 velocity_y_function_->set_time(t);
107 full_primitive_state[2] = velocity_y_function_->value(point);
109 if constexpr (dim > 2) {
110 velocity_z_function_->set_time(t);
111 full_primitive_state[3] = velocity_z_function_->value(point);
114 pressure_function_->set_time(t);
115 full_primitive_state[1 + dim] = pressure_function_->value(point);
117 return view.from_initial_state(full_primitive_state);
123 std::string density_expression_;
124 std::string velocity_x_expression_;
125 std::string velocity_y_expression_;
126 std::string velocity_z_expression_;
127 std::string pressure_expression_;
129 std::unique_ptr<dealii::FunctionParser<dim>> density_function_;
130 std::unique_ptr<dealii::FunctionParser<dim>> velocity_x_function_;
131 std::unique_ptr<dealii::FunctionParser<dim>> velocity_y_function_;
132 std::unique_ptr<dealii::FunctionParser<dim>> velocity_z_function_;
133 std::unique_ptr<dealii::FunctionParser<dim>> pressure_function_;
typename Description::HyperbolicSystem HyperbolicSystem
state_type compute(const dealii::Point< dim > &point, Number t) final
Function(const HyperbolicSystem &hyperbolic_system, const std::string subsection)
typename Description::template HyperbolicSystemView< dim, Number > View
typename View::state_type state_type
Euler::HyperbolicSystem HyperbolicSystem