ryujin 2.1.1 revision 350e54cc11f3d21282bddcf3e6517944dbc508bf
initial_state_function.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0
3// [LANL Copyright Statement]
4// Copyright (C) 2023 - 2024 by the ryujin authors
5// Copyright (C) 2023 - 2024 by Triad National Security, LLC
6//
7
8#pragma once
9
11
12#include <deal.II/base/function_parser.h>
13
14namespace ryujin
15{
16 namespace ShallowWaterInitialStates
17 {
24 template <typename Description, int dim, typename Number>
25 class Function : public InitialState<Description, dim, Number>
26 {
27 public:
29 using View =
30 typename Description::template HyperbolicSystemView<dim, Number>;
31 using state_type = typename View::state_type;
32
33 Function(const HyperbolicSystem &hyperbolic_system,
34 const std::string subsection)
35 : InitialState<Description, dim, Number>("function", subsection)
36 , hyperbolic_system_(hyperbolic_system)
37 {
38
39 depth_expression_ = "1.4";
40 this->add_parameter(
41 "water elevation expression",
42 depth_expression_,
43 "A function expression describing the water elevation. When "
44 "bathymetry is 0, this reduces to the water depth.");
45
46 bathymetry_expression_ = "0.";
47 this->add_parameter("bathymetry expression",
48 bathymetry_expression_,
49 "A function expression describing the bathymetry");
50
51
52 velocity_x_expression_ = "3.0";
53 this->add_parameter(
54 "velocity x expression",
55 velocity_x_expression_,
56 "A function expression describing the x-component of the velocity");
57
58 if constexpr (dim > 1) {
59 velocity_y_expression_ = "0.0";
60 this->add_parameter("velocity y expression",
61 velocity_y_expression_,
62 "A function expression describing the "
63 "y-component of the velocity");
64 }
65
66 /*
67 * Set up the muparser object with the final flux description from
68 * the parameter file:
69 */
70 const auto set_up_muparser = [this] {
71 using FP = dealii::FunctionParser<dim>;
72 /*
73 * This variant of the constructor initializes the function
74 * parser with support for a time-dependent description involving
75 * a variable »t«:
76 */
77 depth_function_ = std::make_unique<FP>(depth_expression_);
78 velocity_x_function_ = std::make_unique<FP>(velocity_x_expression_);
79 if constexpr (dim > 1)
80 velocity_y_function_ = std::make_unique<FP>(velocity_y_expression_);
81 bathymetry_function_ = std::make_unique<FP>(bathymetry_expression_);
82 };
83
84 set_up_muparser();
85 this->parse_parameters_call_back.connect(set_up_muparser);
86 }
87
88 state_type compute(const dealii::Point<dim> &point, Number t) final
89 {
90 const auto view = hyperbolic_system_.template view<dim, Number>();
91 state_type full_primitive;
92
93 /* Compute bathymetry */
94 const Number z = compute_bathymetry(point);
95
96 depth_function_->set_time(t);
97 full_primitive[0] =
98 std::max(Number(depth_function_->value(point)) - z, Number(0.));
99
100 velocity_x_function_->set_time(t);
101 full_primitive[1] = velocity_x_function_->value(point);
102
103 if constexpr (dim > 1) {
104 velocity_y_function_->set_time(t);
105 full_primitive[2] = velocity_y_function_->value(point);
106 }
107
108 return view.from_primitive_state(full_primitive);
109 }
110
111 auto initial_precomputations(const dealii::Point<dim> &point) ->
113 initial_precomputed_type final
114 {
115 /* Compute bathymetry: */
116 return {compute_bathymetry(point)};
117 }
118
119 private:
120 const HyperbolicSystem &hyperbolic_system_;
121
122 std::string depth_expression_;
123 std::string velocity_x_expression_;
124 std::string velocity_y_expression_;
125 std::string bathymetry_expression_;
126
127 std::unique_ptr<dealii::FunctionParser<dim>> depth_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>> bathymetry_function_;
131
132 DEAL_II_ALWAYS_INLINE inline Number
133 compute_bathymetry(const dealii::Point<dim> &point) const
134 {
135 bathymetry_function_->set_time(0.);
136 return bathymetry_function_->value(point);
137 }
138 };
139 } // namespace ShallowWaterInitialStates
140} // namespace ryujin
typename Description::template HyperbolicSystemView< dim, Number > View
auto initial_precomputations(const dealii::Point< dim > &point) -> typename InitialState< Description, dim, Number >::initial_precomputed_type final
state_type compute(const dealii::Point< dim > &point, Number t) final
typename Description::HyperbolicSystem HyperbolicSystem
Function(const HyperbolicSystem &hyperbolic_system, const std::string subsection)
Euler::HyperbolicSystem HyperbolicSystem
Definition: description.h:32