12#include <deal.II/base/function_parser.h>
16 namespace ShallowWaterInitialStates
24 template <
typename Description,
int dim,
typename Number>
34 const std::string subsection)
36 , hyperbolic_system_(hyperbolic_system)
39 depth_expression_ =
"1.4";
41 "water elevation expression",
43 "A function expression describing the water elevation. When "
44 "bathymetry is 0, this reduces to the water depth.");
46 bathymetry_expression_ =
"0.";
47 this->add_parameter(
"bathymetry expression",
48 bathymetry_expression_,
49 "A function expression describing the bathymetry");
52 velocity_x_expression_ =
"3.0";
54 "velocity x expression",
55 velocity_x_expression_,
56 "A function expression describing the x-component of the velocity");
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");
70 const auto set_up_muparser = [
this] {
71 using FP = dealii::FunctionParser<dim>;
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_);
85 this->parse_parameters_call_back.connect(set_up_muparser);
90 const auto view = hyperbolic_system_.template view<dim, Number>();
94 const Number z = compute_bathymetry(point);
96 depth_function_->set_time(t);
98 std::max(Number(depth_function_->value(point)) - z, Number(0.));
100 velocity_x_function_->set_time(t);
101 full_primitive[1] = velocity_x_function_->value(point);
103 if constexpr (dim > 1) {
104 velocity_y_function_->set_time(t);
105 full_primitive[2] = velocity_y_function_->value(point);
108 return view.from_primitive_state(full_primitive);
113 initial_precomputed_type
final
116 return {compute_bathymetry(point)};
122 std::string depth_expression_;
123 std::string velocity_x_expression_;
124 std::string velocity_y_expression_;
125 std::string bathymetry_expression_;
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_;
132 DEAL_II_ALWAYS_INLINE
inline Number
133 compute_bathymetry(
const dealii::Point<dim> &point)
const
135 bathymetry_function_->set_time(0.);
136 return bathymetry_function_->value(point);
typename View::state_type state_type
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