8#include <compile_time_options.h>
14 namespace GridGenerator
29 template <
int dim,
int spacedim,
template <
int,
int>
class Triangulation>
36 AssertThrow(
false, dealii::ExcNotImplemented());
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)
49 using namespace dealii;
51 dealii::Triangulation<2, 2> res1, res2, res3, flume,
final;
53 const double tolerance = 1.e-8;
55 Assert(reservoir_width - flume_width > tolerance,
56 dealii::ExcInternalError());
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.));
65 GridGenerator::subdivided_hyper_rectangle(
68 Point<2>(-reservoir_length, -reservoir_width / 2.),
69 Point<2>(0, -flume_width / 2.));
71 GridGenerator::subdivided_hyper_rectangle(
74 Point<2>(-reservoir_length, flume_width / 2.),
75 Point<2>(0, reservoir_width / 2.));
77 sub_y = int(std::round(flume_width * 100.));
79 GridGenerator::subdivided_hyper_rectangle(
82 Point<2>(-reservoir_length, -flume_width / 2.),
83 Point<2>(0, flume_width / 2.));
85 sub_x = int(std::round(flume_length * 100.));
87 GridGenerator::subdivided_hyper_rectangle(
90 Point<2>(0., -flume_width / 2.),
91 Point<2>(flume_length, flume_width / 2.));
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);
102 for (
auto cell : triangulation.active_cell_iterators()) {
103 for (
auto f : cell->face_indices()) {
104 const auto face = cell->face(f);
106 if (!face->at_boundary())
115 const auto center = face->center();
116 if (center[0] > flume_length - tolerance)
138 :
Geometry<dim>(
"wave tank", subsection)
140 reservoir_length_ = 157. / 100.;
141 this->add_parameter(
"reservoir length",
143 "length of water reservoir [meters]");
145 reservoir_width_ = 8.1 / 100.;
146 this->add_parameter(
"reservoir width",
148 "width of water reservoir [meters]");
150 flume_length_ = 600.78 / 100.;
152 "flume length", flume_length_,
"length of flume [meters]");
154 flume_width_ = 24 / 100.;
156 "flume width", flume_width_,
"width of flume [meters]");
160 typename dealii::Triangulation<dim> &triangulation)
final
170 double reservoir_length_;
171 double reservoir_width_;
172 double flume_length_;
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)