8#include <compile_time_options.h>
14#include <deal.II/base/function_parser.h>
18 namespace ShallowWaterInitialStates
27 template <
typename Description,
int dim,
typename Number>
38 const std::string subsection)
40 , hyperbolic_system_(hyperbolic_system)
41 , geotiff_reader_(subsection +
"/geotiff")
43 height_expression_ =
"1.4";
45 "water height expression",
47 "A function expression describing the initial total water height");
49 velocity_expression_ =
"0.0";
51 "velocity expression",
53 "A function expression describing the initial velocity");
55 const auto set_up = [
this] {
56 using FP = dealii::FunctionParser<dim>;
62 height_function_ = std::make_unique<FP>(height_expression_);
63 velocity_function_ = std::make_unique<FP>(velocity_expression_);
67 this->parse_parameters_call_back.connect(set_up);
72 const auto z = geotiff_reader_.compute_height(point);
74 dealii::Tensor<1, 2, Number> primitive;
76 height_function_->set_time(t);
77 primitive[0] = std::max(0., height_function_->value(point) - z);
79 velocity_function_->set_time(t);
80 primitive[1] = velocity_function_->value(point);
82 const auto view = hyperbolic_system_.template view<dim, Number>();
83 return view.from_initial_state(primitive);
88 initial_precomputed_type
final
91 return {
static_cast<Number
>(geotiff_reader_.compute_height(point))};
100 std::string height_expression_;
101 std::string velocity_expression_;
105 std::unique_ptr<dealii::FunctionParser<dim>> height_function_;
106 std::unique_ptr<dealii::FunctionParser<dim>> velocity_function_;
typename Description::HyperbolicSystem HyperbolicSystem
auto initial_precomputations(const dealii::Point< dim > &point) -> typename InitialState< Description, dim, Number >::initial_precomputed_type final
typename View::state_type state_type
state_type compute(const dealii::Point< dim > &point, Number t) final
GeoTIFF(const HyperbolicSystem &hyperbolic_system, const std::string subsection)
typename Description::template HyperbolicSystemView< dim, Number > View
Euler::HyperbolicSystem HyperbolicSystem