ryujin 2.1.1 revision ef0fcd4010d109b860652ece4a7b8963fb7d46b1
geometry_tank.h
Go to the documentation of this file.
1//
2// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
3// Copyright (C) 2024 by the ryujin authors
4//
5
6#pragma once
7
8#include <compile_time_options.h>
9
11
12namespace ryujin
13{
14 namespace GridGenerator
15 {
29 template <int dim, int spacedim, template <int, int> class Triangulation>
30 void wavetank(Triangulation<dim, spacedim> &,
31 const double /*reservoir_length*/,
32 const double /*reservoir_width*/,
33 const double /*flume_length*/,
34 const double /*flume_width*/)
35 {
36 AssertThrow(false, dealii::ExcNotImplemented());
37 __builtin_trap();
38 }
39
40
41#ifndef DOXYGEN
42 template <template <int, int> class Triangulation>
43 void wavetank(Triangulation<2, 2> &triangulation,
44 const double reservoir_length,
45 const double reservoir_width,
46 const double flume_length,
47 const double flume_width)
48 {
49 using namespace dealii;
50
51 dealii::Triangulation<2, 2> res1, res2, res3, flume, final;
52
53 const double tolerance = 1.e-8;
54
55 Assert(reservoir_width - flume_width > tolerance,
56 dealii::ExcInternalError());
57
58 /* We split the reservoir into three triangulations and subdivide to
59 * get close to uniform refinement */
60
61 const double diff = (reservoir_width - flume_width) / 2.;
62 unsigned int sub_x = int(std::round(reservoir_length * 100.));
63 unsigned int sub_y = int(std::round(diff * 100.));
64
65 GridGenerator::subdivided_hyper_rectangle(
66 res1,
67 {sub_x, sub_y},
68 Point<2>(-reservoir_length, -reservoir_width / 2.),
69 Point<2>(0, -flume_width / 2.));
70
71 GridGenerator::subdivided_hyper_rectangle(
72 res3,
73 {sub_x, sub_y},
74 Point<2>(-reservoir_length, flume_width / 2.),
75 Point<2>(0, reservoir_width / 2.));
76
77 sub_y = int(std::round(flume_width * 100.));
78
79 GridGenerator::subdivided_hyper_rectangle(
80 res2,
81 {sub_x, sub_y},
82 Point<2>(-reservoir_length, -flume_width / 2.),
83 Point<2>(0, flume_width / 2.));
84
85 sub_x = int(std::round(flume_length * 100.));
86
87 GridGenerator::subdivided_hyper_rectangle(
88 flume,
89 {sub_x, sub_y},
90 Point<2>(0., -flume_width / 2.),
91 Point<2>(flume_length, flume_width / 2.));
92
93 final.set_mesh_smoothing(triangulation.get_mesh_smoothing());
94 GridGenerator::merge_triangulations(
95 {&res1, &res2, &res3, &flume}, final, tolerance);
96 triangulation.copy_triangulation(final);
97
98 /*
99 * Set boundary ids:
100 */
101
102 for (auto cell : triangulation.active_cell_iterators()) {
103 for (auto f : cell->face_indices()) {
104 const auto face = cell->face(f);
105
106 if (!face->at_boundary())
107 continue;
108
109 /*
110 * We want slip everywhere except right edge of tank.
111 */
112
113 face->set_boundary_id(Boundary::slip);
114
115 const auto center = face->center();
116 if (center[0] > flume_length - tolerance)
117 face->set_boundary_id(Boundary::dynamic);
118
119 } /*f*/
120 } /*cell*/
121 }
122#endif
123 } /* namespace GridGenerator */
124
125
126 namespace Geometries
127 {
133 template <int dim>
134 class WaveTank : public Geometry<dim>
135 {
136 public:
137 WaveTank(const std::string subsection)
138 : Geometry<dim>("wave tank", subsection)
139 {
140 reservoir_length_ = 157. / 100.;
141 this->add_parameter("reservoir length",
142 reservoir_length_,
143 "length of water reservoir [meters]");
144
145 reservoir_width_ = 8.1 / 100.;
146 this->add_parameter("reservoir width",
147 reservoir_width_,
148 "width of water reservoir [meters]");
149
150 flume_length_ = 600.78 / 100.;
151 this->add_parameter(
152 "flume length", flume_length_, "length of flume [meters]");
153
154 flume_width_ = 24 / 100.;
155 this->add_parameter(
156 "flume width", flume_width_, "width of flume [meters]");
157 }
158
160 typename dealii::Triangulation<dim> &triangulation) final
161 {
162 GridGenerator::wavetank(triangulation,
163 reservoir_length_,
164 reservoir_width_,
165 flume_length_,
166 flume_width_);
167 }
168
169 private:
170 double reservoir_length_;
171 double reservoir_width_;
172 double flume_length_;
173 double flume_width_;
174 };
175 } /* namespace Geometries */
176} /* namespace ryujin */
WaveTank(const std::string subsection)
void create_triangulation(typename dealii::Triangulation< dim > &triangulation) final
void wavetank(Triangulation< dim, spacedim > &, const double, const double, const double, const double)
Definition: geometry_tank.h:30