ryujin 2.1.1 revision 9dcb748690310d6a540ebb8b066d1a0834fc7604
initial_state_ritter_dam_break.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: MIT
3// Copyright (C) 2020 - 2023 by the ryujin authors
4//
5
6#pragma once
7
8#include "hyperbolic_system.h"
9#include <initial_state.h>
10
11namespace ryujin
12{
13 namespace ShallowWater
14 {
21 template <int dim, typename Number, typename state_type>
22 class RitterDamBreak : public InitialState<dim, Number, state_type, 1>
23 {
24 public:
25 RitterDamBreak(const HyperbolicSystem &hyperbolic_system,
26 const std::string subsection)
27 : InitialState<dim, Number, state_type, 1>("ritter dam break",
28 subsection)
29 , hyperbolic_system(hyperbolic_system)
30 {
31 dealii::ParameterAcceptor::parse_parameters_call_back.connect(
33
34 t_initial_ = 0.1;
35 this->add_parameter("time initial",
36 t_initial_,
37 "Time at which initial state is prescribed");
38
39 left_depth = 0.005;
40 this->add_parameter("left water depth",
41 left_depth,
42 "Depth of water to the left of pseudo-dam (x<0)");
43 }
44
46 {
47 AssertThrow(t_initial_ > 0.,
48 dealii::ExcMessage("Expansion must be computed at an "
49 "initial time greater than 0."));
50 }
51
52 state_type compute(const dealii::Point<dim> &point, Number t) final
53 {
54 const auto g = hyperbolic_system.gravity();
55 const Number x = point[0];
56
57 const Number aL = std::sqrt(g * left_depth);
58 const Number xA = -(t + t_initial_) * aL;
59 const Number xB = Number(2.) * (t + t_initial_) * aL;
60
61 const Number tmp = aL - x / (2. * (t + t_initial_));
62
63 const Number h_expansion = 4. / (9. * g) * tmp * tmp;
64 const Number v_expansion = 2. / 3. * (x / (t + t_initial_) + aL);
65
66 if (x <= xA)
67 return hyperbolic_system.template expand_state<dim>(
69 {left_depth, Number(0.)}});
70 else if (x <= xB)
71 return hyperbolic_system.template expand_state<dim>(
73 {h_expansion, h_expansion * v_expansion}});
74 else
75 return hyperbolic_system.template expand_state<dim>(
77 {Number(0.), Number(0.)}});
78 }
79
80 /* Default bathymetry of 0 */
81
82 private:
83 const HyperbolicSystem &hyperbolic_system;
84
85 Number t_initial_;
86
87 Number left_depth;
88 };
89
90 } // namespace ShallowWater
91} // namespace ryujin
typename HyperbolicSystemView::state_type state_type
dealii::Tensor< 1, problem_dimension< dim >, Number > state_type
state_type compute(const dealii::Point< dim > &point, Number t) final
RitterDamBreak(const HyperbolicSystem &hyperbolic_system, const std::string subsection)