ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
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
8#include <compile_time_options.h>
9
11
12#include <deal.II/base/function_parser.h>
13
14namespace ryujin
15{
16 namespace EulerInitialStates
17 {
23 template <typename Description, int dim, typename Number>
24 class Function : public InitialState<Description, dim, Number>
25 {
26 public:
28 using View =
29 typename Description::template HyperbolicSystemView<dim, Number>;
30 using state_type = typename View::state_type;
31
32 Function(const HyperbolicSystem &hyperbolic_system,
33 const std::string subsection)
34 : InitialState<Description, dim, Number>("function", subsection)
35 , hyperbolic_system_(hyperbolic_system)
36 {
37
38 density_expression_ = "1.4";
39 this->add_parameter("density expression",
40 density_expression_,
41 "A function expression describing the density");
42
43 velocity_x_expression_ = "3.0";
44 this->add_parameter(
45 "velocity x expression",
46 velocity_x_expression_,
47 "A function expression describing the x-component of the velocity");
48
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");
55 }
56
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");
63 }
64
65 pressure_expression_ = "1.0";
66 this->add_parameter("pressure expression",
67 pressure_expression_,
68 "A function expression describing the pressure");
69
70 /*
71 * Set up the muparser object with the final flux description from
72 * the parameter file:
73 */
74 const auto set_up_muparser = [this] {
75 using FP = dealii::FunctionParser<dim>;
76 /*
77 * This variant of the constructor initializes the function
78 * parser with support for a time-dependent description involving
79 * a variable »t«:
80 */
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_);
88 };
89
90 set_up_muparser();
91 this->parse_parameters_call_back.connect(set_up_muparser);
92 }
93
94 state_type compute(const dealii::Point<dim> &point, Number t) final
95 {
96 const auto view = hyperbolic_system_.template view<dim, Number>();
97 state_type full_primitive_state;
98
99 density_function_->set_time(t);
100 full_primitive_state[0] = density_function_->value(point);
101
102 velocity_x_function_->set_time(t);
103 full_primitive_state[1] = velocity_x_function_->value(point);
104
105 if constexpr (dim > 1) {
106 velocity_y_function_->set_time(t);
107 full_primitive_state[2] = velocity_y_function_->value(point);
108 }
109 if constexpr (dim > 2) {
110 velocity_z_function_->set_time(t);
111 full_primitive_state[3] = velocity_z_function_->value(point);
112 }
113
114 pressure_function_->set_time(t);
115 full_primitive_state[1 + dim] = pressure_function_->value(point);
116
117 return view.from_initial_state(full_primitive_state);
118 }
119
120 private:
121 const HyperbolicSystem &hyperbolic_system_;
122
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_;
128
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_;
134 };
135 } // namespace EulerInitialStates
136} // 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:34