ryujin 2.1.1 revision 863a4d36dcc743d4e1a9b41cfabd03d0aea57863
initial_state_smooth_wave.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#include <simd.h>
10
11namespace ryujin
12{
13 namespace EulerInitialStates
14 {
26 template <typename Description, int dim, typename Number>
27 class SmoothWave : public InitialState<Description, dim, Number>
28 {
29 public:
31 using View =
32 typename Description::template HyperbolicSystemView<dim, Number>;
33 using state_type = typename View::state_type;
34
35 using ScalarNumber = typename View::ScalarNumber;
36
37 SmoothWave(const HyperbolicSystem &hyperbolic_system,
38 const std::string subsection)
39 : InitialState<Description, dim, Number>("smooth wave", subsection)
40 , hyperbolic_system_(hyperbolic_system)
41 {
42 density_ref_ = 1.;
43 this->add_parameter("reference density",
44 density_ref_,
45 "The material reference density");
46
47 pressure_ref_ = 1.;
48 this->add_parameter("reference pressure",
49 pressure_ref_,
50 "The material reference pressure");
51
52 mach_number_ = 1.0;
53 this->add_parameter("mach number",
54 mach_number_,
55 "Mach number of traveling smooth wave");
56
57 /* These are the x_0 and x_1 parameters from references above. */
58 left_ = 0.1;
59 right_ = 0.3;
60 }
61
62 state_type compute(const dealii::Point<dim> &point, Number t) final
63 {
64 const auto view = hyperbolic_system_.template view<dim, Number>();
65
66 auto point_bar = point;
67 point_bar[0] = point_bar[0] - mach_number_ * t;
68 const auto x = Number(point_bar[0]);
69
70 const Number polynomial = Number(64) *
71 ryujin::fixed_power<3>(x - left_) *
72 ryujin::fixed_power<3>(right_ - x) /
73 ryujin::fixed_power<6>(right_ - left_);
74
75 /* Define density profile */
76 Number rho = density_ref_;
77 if (left_ <= point_bar[0] && point_bar[0] <= right_)
78 rho = density_ref_ + polynomial;
79
80 state_type initial_state;
81 {
82 initial_state[0] = rho;
83 initial_state[1] = mach_number_;
84 initial_state[dim + 1] = pressure_ref_;
85 }
86 return view.from_initial_state(initial_state);
87 }
88
89 private:
90 const HyperbolicSystem &hyperbolic_system_;
91
92 Number density_ref_;
93 Number pressure_ref_;
94 Number mach_number_;
95 Number left_;
96 Number right_;
97 };
98 } // namespace EulerInitialStates
99} // namespace ryujin
SmoothWave(const HyperbolicSystem &hyperbolic_system, const std::string subsection)
typename Description::HyperbolicSystem HyperbolicSystem
state_type compute(const dealii::Point< dim > &point, Number t) final
typename Description::template HyperbolicSystemView< dim, Number > View
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:32