ryujin 2.1.1 revision 863a4d36dcc743d4e1a9b41cfabd03d0aea57863
initial_state_function.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2023 - 2024 by the ryujin authors
4//
5
6#pragma once
7
9
10#include <deal.II/base/function_parser.h>
11
12namespace ryujin
13{
14 namespace EulerInitialStates
15 {
21 template <typename Description, int dim, typename Number>
22 class Function : public InitialState<Description, dim, Number>
23 {
24 public:
26 using View =
27 typename Description::template HyperbolicSystemView<dim, Number>;
28 using state_type = typename View::state_type;
29
30 Function(const HyperbolicSystem &hyperbolic_system,
31 const std::string subsection)
32 : InitialState<Description, dim, Number>("function", subsection)
33 , hyperbolic_system_(hyperbolic_system)
34 {
35
36 density_expression_ = "1.4";
37 this->add_parameter("density expression",
38 density_expression_,
39 "A function expression describing the density");
40
41 velocity_x_expression_ = "3.0";
42 this->add_parameter(
43 "velocity x expression",
44 velocity_x_expression_,
45 "A function expression describing the x-component of the velocity");
46
47 if constexpr (dim > 1) {
48 velocity_y_expression_ = "0.0";
49 this->add_parameter("velocity y expression",
50 velocity_y_expression_,
51 "A function expression describing the "
52 "y-component of the velocity");
53 }
54
55 if constexpr (dim > 2) {
56 velocity_z_expression_ = "0.0";
57 this->add_parameter("velocity z expression",
58 velocity_z_expression_,
59 "A function expression describing the "
60 "z-component of the velocity");
61 }
62
63 pressure_expression_ = "1.0";
64 this->add_parameter("pressure expression",
65 pressure_expression_,
66 "A function expression describing the pressure");
67
68 /*
69 * Set up the muparser object with the final flux description from
70 * the parameter file:
71 */
72 const auto set_up_muparser = [this] {
73 using FP = dealii::FunctionParser<dim>;
74 /*
75 * This variant of the constructor initializes the function
76 * parser with support for a time-dependent description involving
77 * a variable »t«:
78 */
79 density_function_ = std::make_unique<FP>(density_expression_);
80 velocity_x_function_ = std::make_unique<FP>(velocity_x_expression_);
81 if constexpr (dim > 1)
82 velocity_y_function_ = std::make_unique<FP>(velocity_y_expression_);
83 if constexpr (dim > 2)
84 velocity_z_function_ = std::make_unique<FP>(velocity_z_expression_);
85 pressure_function_ = std::make_unique<FP>(pressure_expression_);
86 };
87
88 set_up_muparser();
89 this->parse_parameters_call_back.connect(set_up_muparser);
90 }
91
92 state_type compute(const dealii::Point<dim> &point, Number t) final
93 {
94 const auto view = hyperbolic_system_.template view<dim, Number>();
95 state_type full_primitive_state;
96
97 density_function_->set_time(t);
98 full_primitive_state[0] = density_function_->value(point);
99
100 velocity_x_function_->set_time(t);
101 full_primitive_state[1] = velocity_x_function_->value(point);
102
103 if constexpr (dim > 1) {
104 velocity_y_function_->set_time(t);
105 full_primitive_state[2] = velocity_y_function_->value(point);
106 }
107 if constexpr (dim > 2) {
108 velocity_z_function_->set_time(t);
109 full_primitive_state[3] = velocity_z_function_->value(point);
110 }
111
112 pressure_function_->set_time(t);
113 full_primitive_state[1 + dim] = pressure_function_->value(point);
114
115 return view.from_primitive_state(full_primitive_state);
116 }
117
118 private:
119 const HyperbolicSystem &hyperbolic_system_;
120
121 std::string density_expression_;
122 std::string velocity_x_expression_;
123 std::string velocity_y_expression_;
124 std::string velocity_z_expression_;
125 std::string pressure_expression_;
126
127 std::unique_ptr<dealii::FunctionParser<dim>> density_function_;
128 std::unique_ptr<dealii::FunctionParser<dim>> velocity_x_function_;
129 std::unique_ptr<dealii::FunctionParser<dim>> velocity_y_function_;
130 std::unique_ptr<dealii::FunctionParser<dim>> velocity_z_function_;
131 std::unique_ptr<dealii::FunctionParser<dim>> pressure_function_;
132 };
133 } // namespace EulerInitialStates
134} // namespace ryujin
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
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:32