9#include <initial_state.h>
13 namespace ShallowWater
21 template <
int dim,
typename Number,
typename state_type>
26 const std::string subsection)
29 , hyperbolic_system(hyperbolic_system)
31 dealii::ParameterAcceptor::parse_parameters_call_back.connect(
35 this->add_parameter(
"time initial",
37 "Time at which initial state is prescribed");
40 this->add_parameter(
"left water depth",
42 "Depth of water to the left of pseudo-dam (x<0)");
47 AssertThrow(t_initial_ > 0.,
48 dealii::ExcMessage(
"Expansion must be computed at an "
49 "initial time greater than 0."));
54 const auto g = hyperbolic_system.gravity();
55 const Number x = point[0];
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;
61 const Number tmp = aL - x / (2. * (t + t_initial_));
63 const Number h_expansion = 4. / (9. * g) * tmp * tmp;
64 const Number v_expansion = 2. / 3. * (x / (t + t_initial_) + aL);
67 return hyperbolic_system.template expand_state<dim>(
69 {left_depth, Number(0.)}});
71 return hyperbolic_system.template expand_state<dim>(
73 {h_expansion, h_expansion * v_expansion}});
75 return hyperbolic_system.template expand_state<dim>(
77 {Number(0.), Number(0.)}});
typename HyperbolicSystemView::state_type state_type
dealii::Tensor< 1, problem_dimension< dim >, Number > state_type
void parse_parameters_callback()
state_type compute(const dealii::Point< dim > &point, Number t) final
RitterDamBreak(const HyperbolicSystem &hyperbolic_system, const std::string subsection)