ryujin 2.1.1 revision 9dcb748690310d6a540ebb8b066d1a0834fc7604
initial_state_unsteady_vortex.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 {
20 template <int dim, typename Number, typename state_type>
21 class UnsteadyVortex : public InitialState<dim, Number, state_type, 1>
22 {
23 public:
24 UnsteadyVortex(const HyperbolicSystem &hyperbolic_system,
25 const std::string subsection)
26 : InitialState<dim, Number, state_type, 1>("unsteady vortex",
27 subsection)
28 , hyperbolic_system(hyperbolic_system)
29 {
30 depth_ = 1.0;
31 this->add_parameter("reference depth", depth_, "Reference water depth");
32
33 mach_number_ = 2.0;
34 this->add_parameter(
35 "mach number", mach_number_, "Mach number of unsteady vortex");
36
37 beta_ = 0.1;
38 this->add_parameter("beta", beta_, "vortex strength beta");
39 }
40
41 state_type compute(const dealii::Point<dim> &point, Number t) final
42 {
43 const auto gravity = hyperbolic_system.gravity();
44
45 dealii::Point<2> point_bar;
46 point_bar[0] = point[0] - mach_number_ * t;
47 point_bar[1] = point[1];
48
49 const Number r_square = Number(point_bar.norm_square());
50
51 /* We assume r0 = 1 */
52 const Number factor =
53 beta_ / Number(2. * M_PI) *
54 exp(Number(0.5) - Number(0.5) * r_square /* /r0^ 2*/);
55
56 const Number h =
57 depth_ - Number(1. / (2. * gravity /* * r0^2 */)) * factor * factor;
58 const Number u = mach_number_ - factor * Number(point_bar[1]);
59 const Number v = factor * Number(point_bar[0]);
60
61
62 if constexpr (dim == 2)
63 return hyperbolic_system.template expand_state<dim>(
65 else {
66 AssertThrow(false, dealii::ExcNotImplemented());
67 __builtin_trap();
68 }
69 }
70
71 private:
72 const HyperbolicSystem &hyperbolic_system;
73
74 Number depth_;
75 Number mach_number_;
76 Number beta_;
77 };
78
79 } // namespace ShallowWater
80} // 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
UnsteadyVortex(const HyperbolicSystem &hyperbolic_system, const std::string subsection)
Euler::HyperbolicSystem HyperbolicSystem