From 28070fe6ad63b9b5cfc5f22424a0f906e2c342e7 Mon Sep 17 00:00:00 2001
From: Felix Schindler <felix.schindler@wwu.de>
Date: Wed, 8 Aug 2018 00:24:04 +0200
Subject: [PATCH] [test.inviscid-compressible-flow] add another 1d euler test
case
---
dune/gdt/test/instationary-eocstudies/base.hh | 2 +-
...mpressible_flow__euler_1d__explicit__fv.cc | 165 +++++++++++++++---
2 files changed, 146 insertions(+), 21 deletions(-)
diff --git a/dune/gdt/test/instationary-eocstudies/base.hh b/dune/gdt/test/instationary-eocstudies/base.hh
index 3f2911fef..28136a53f 100644
--- a/dune/gdt/test/instationary-eocstudies/base.hh
+++ b/dune/gdt/test/instationary-eocstudies/base.hh
@@ -299,7 +299,7 @@ protected:
return time_points;
}
- const double T_end_;
+ double T_end_;
const size_t num_refinements_;
const size_t num_additional_refinements_for_reference_;
const std::function<void(const DiscreteBochnerFunction<V, GV, m>&, const std::string&)> visualize_;
diff --git a/dune/gdt/test/inviscid-compressible-flow/inviscid_compressible_flow__euler_1d__explicit__fv.cc b/dune/gdt/test/inviscid-compressible-flow/inviscid_compressible_flow__euler_1d__explicit__fv.cc
index 9b013da8e..b597786e1 100644
--- a/dune/gdt/test/inviscid-compressible-flow/inviscid_compressible_flow__euler_1d__explicit__fv.cc
+++ b/dune/gdt/test/inviscid-compressible-flow/inviscid_compressible_flow__euler_1d__explicit__fv.cc
@@ -16,6 +16,8 @@
#include <dune/xt/grid/grids.hh>
#include <dune/xt/grid/walker.hh>
#include <dune/xt/functions/base/sliced.hh>
+#include <dune/xt/functions/interfaces/element-functions.hh>
+#include <dune/xt/functions/lambda/function.hh>
#include <dune/gdt/local/integrands/identity.hh>
#include <dune/gdt/test/instationary-eocstudies/hyperbolic-fv.hh>
@@ -40,6 +42,7 @@ protected:
using BaseType::m;
using typename BaseType::E;
using typename BaseType::F;
+ using typename BaseType::R;
using typename BaseType::DF;
using typename BaseType::GV;
using typename BaseType::GP;
@@ -78,7 +81,7 @@ public:
std::vector<std::string> quantities() const override final
{
- return {"rel mass conserv error"};
+ return {"rel mass conserv error"}; // <- This is on purpose, see column header formatting in ConvergenceStudy.
}
virtual std::map<std::string, std::map<std::string, double>>
@@ -124,8 +127,35 @@ protected:
DF make_initial_values(const S& space) override final
{
- return Problem::access().template make_initial_values<V>(space);
- }
+ if (boundary_treatment == "inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right") {
+ const auto& euler_tools = Problem::access().euler_tools;
+ return interpolate<V>(0,
+ [&](const auto& /*xx*/, const auto& /*param*/) {
+ FieldVector<R, m> primitive_variables(0.);
+ // density
+ primitive_variables[0] = 0.5;
+ // velocity
+ for (size_t ii = 0; ii < d; ++ii)
+ primitive_variables[1 + ii] = 0.;
+ // pressure
+ primitive_variables[m - 1] = 0.4;
+ return euler_tools.to_conservative(primitive_variables);
+ },
+ space);
+ } else
+ return Problem::access().template make_initial_values<V>(space);
+ } // ... make_initial_values(...)
+
+ double estimate_dt(const S& space) override final
+ {
+ if (boundary_treatment == "inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right")
+ return estimate_dt_for_hyperbolic_system(space.grid_view(),
+ make_initial_values(space),
+ flux(),
+ /*boundary_data_range=*/{{0.5, 0., 0.4}, {1.5, 0.5, 0.4}});
+ else
+ return estimate_dt_for_hyperbolic_system(space.grid_view(), make_initial_values(space), flux());
+ } // ... estimate_dt(...)
GP make_initial_grid() override final
{
@@ -146,17 +176,18 @@ protected:
boundary_info = std::make_unique<XT::Grid::NormalBasedBoundaryInfo<XT::Grid::extract_intersection_t<GV>>>();
boundary_info->register_new_normal({-1}, new XT::Grid::ImpermeableBoundary());
boundary_info->register_new_normal({1}, new XT::Grid::ImpermeableBoundary());
- XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> impermeable_wall_filter(*boundary_info,
- new XT::Grid::ImpermeableBoundary());
+ const XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> impermeable_wall_filter(
+ *boundary_info, new XT::Grid::ImpermeableBoundary());
auto op = std::make_unique<AdvectionFvOperator<GV, V, m>>(space.grid_view(),
*self.numerical_flux_,
space,
space,
/*periodicity_restriction=*/impermeable_wall_filter);
// the actual handling of impermeable walls
- op->append([&](const auto& u,
+ op->append(/*numerical_boundary_flux=*/[&](
+ const auto& u,
const auto& n,
- const auto& /*param*/ = {}) { return euler_tools.flux_at_impermeable_walls(u, n); },
+ const auto& /*param*/) { return euler_tools.flux_at_impermeable_walls(u, n); },
{},
impermeable_wall_filter);
return op;
@@ -164,8 +195,8 @@ protected:
boundary_info = std::make_unique<XT::Grid::NormalBasedBoundaryInfo<XT::Grid::extract_intersection_t<GV>>>();
boundary_info->register_new_normal({-1}, new XT::Grid::ImpermeableBoundary());
boundary_info->register_new_normal({1}, new XT::Grid::ImpermeableBoundary());
- XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> impermeable_wall_filter(*boundary_info,
- new XT::Grid::ImpermeableBoundary());
+ const XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> impermeable_wall_filter(
+ *boundary_info, new XT::Grid::ImpermeableBoundary());
auto op = std::make_unique<AdvectionFvOperator<GV, V, m>>(space.grid_view(),
*self.numerical_flux_,
space,
@@ -173,11 +204,104 @@ protected:
/*periodicity_restriction=*/impermeable_wall_filter);
// the actual handling of impermeable walls, see [DF2015, p. 415, (8.66 - 8.67)]
op->append(
+ /*boundary_extrapolation=*/
[&](const auto& intersection,
const auto& xx_in_reference_intersection_coordinates,
- const auto& /*f*/,
+ const auto& /*flux*/,
const auto& u,
- const auto& /*mu*/ = {}) {
+ const auto& /*param*/) {
+ const auto normal = intersection.unitOuterNormal(xx_in_reference_intersection_coordinates);
+ const auto rho = euler_tools.density_from_conservative(u);
+ auto velocity = euler_tools.velocity_from_conservative(u);
+ velocity -= normal * 2. * (velocity * normal);
+ const auto pressure = euler_tools.pressure_from_conservative(u);
+ return euler_tools.to_conservative(XT::Common::hstack(rho, velocity, pressure));
+ },
+ {},
+ impermeable_wall_filter);
+ return op;
+ } else if (boundary_treatment == "inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right") {
+ periodic_density_variation_ =
+ std::unique_ptr<XT::Functions::LambdaFunction<d, m>>(new XT::Functions::LambdaFunction<d, m>(
+ /*order=*/0,
+ [&](const auto& /*xx*/, const auto& param) {
+ FieldVector<R, m> primitive_variables(0.);
+ const auto t = param.get("_t").at(0);
+ // density
+ primitive_variables[0] = 1. + 0.5 * std::sin(2. * M_PI * (t - 0.25));
+ // velocity
+ for (size_t ii = 0; ii < d; ++ii)
+ primitive_variables[1 + ii] = 0.25 + 0.25 * std::sin(2. * M_PI * (t - 0.25));
+ // pressure
+ primitive_variables[m - 1] = 0.4;
+ return euler_tools.to_conservative(primitive_variables);
+ },
+ /*name=*/"periodic_density_variation",
+ /*parameter_type=*/{"_t", 1}));
+ local_periodic_density_variation_ = periodic_density_variation_->template as_grid_function<E>().local_function();
+ // inflow/outflow left, impermeable wall right
+ boundary_info = std::make_unique<XT::Grid::NormalBasedBoundaryInfo<XT::Grid::extract_intersection_t<GV>>>();
+ boundary_info->register_new_normal({-1}, new XT::Grid::InflowOutflowBoundary());
+ boundary_info->register_new_normal({1}, new XT::Grid::ImpermeableBoundary());
+ const XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> impermeable_wall_filter(
+ *boundary_info, new XT::Grid::ImpermeableBoundary());
+ const XT::Grid::ApplyOn::CustomBoundaryIntersections<GV> inflow_outflow_filter(
+ *boundary_info, new XT::Grid::InflowOutflowBoundary());
+ auto op = std::make_unique<AdvectionFvOperator<GV, V, m>>(
+ space.grid_view(),
+ *self.numerical_flux_,
+ space,
+ space,
+ /*periodicity_restriction=*/*(inflow_outflow_filter || impermeable_wall_filter));
+ // the actual handling of inflow/outflow, see [DF2015, p. 421, (8.88)]
+ // (supposedly this does not work well for slow flows!)
+ const auto heuristic_euler_inflow_outflow_treatment = [&](const auto& intersection,
+ const auto& xx_in_reference_intersection_coordinates,
+ const auto& /*flux*/,
+ const auto& u,
+ const auto& param) {
+ using RangeType = XT::Common::FieldVector<R, m>;
+ // evaluate boundary values
+ const auto element = intersection.inside();
+ const auto xx_in_reference_element_coordinates =
+ intersection.geometryInInside().global(xx_in_reference_intersection_coordinates);
+ local_periodic_density_variation_->bind(element);
+ const RangeType bv = local_periodic_density_variation_->evaluate(xx_in_reference_element_coordinates, param);
+ // determine flow regime
+ const auto a = euler_tools.speed_of_sound_from_conservative(u);
+ const auto velocity = euler_tools.velocity_from_conservative(u);
+ const auto normal = intersection.unitOuterNormal(xx_in_reference_intersection_coordinates);
+ const auto flow_speed = velocity * normal;
+ // compute v
+ if (flow_speed < -a) {
+ // supersonic inlet
+ return bv;
+ } else if (!(flow_speed > 0)) {
+ // subsonic inlet
+ const auto rho_outer = euler_tools.density_from_conservative(bv);
+ const auto v_outer = euler_tools.velocity_from_conservative(bv);
+ const auto p_inner = euler_tools.pressure_from_conservative(u);
+ return euler_tools.to_conservative(XT::Common::hstack(rho_outer, v_outer, p_inner));
+ } else if (flow_speed < a) {
+ // subsonic outlet
+ const auto rho_inner = euler_tools.density_from_conservative(u);
+ const auto v_inner = euler_tools.velocity_from_conservative(u);
+ const auto p_outer = euler_tools.pressure_from_conservative(bv);
+ return euler_tools.to_conservative(XT::Common::hstack(rho_inner, v_inner, p_outer));
+ } else {
+ // supersonic outlet
+ return RangeType(u);
+ }
+ }; // ... heuristic_euler_inflow_outflow_treatment(...)
+ op->append(heuristic_euler_inflow_outflow_treatment, {}, inflow_outflow_filter);
+ // the actual handling of impermeable walls, see above
+ op->append(
+ /*boundary_extrapolation=*/
+ [&](const auto& intersection,
+ const auto& xx_in_reference_intersection_coordinates,
+ const auto& /*flux*/,
+ const auto& u,
+ const auto& /*param*/) {
const auto normal = intersection.unitOuterNormal(xx_in_reference_intersection_coordinates);
const auto rho = euler_tools.density_from_conservative(u);
auto velocity = euler_tools.velocity_from_conservative(u);
@@ -203,6 +327,8 @@ protected:
size_t visualization_steps_;
std::string boundary_treatment;
std::unique_ptr<XT::Grid::NormalBasedBoundaryInfo<XT::Grid::extract_intersection_t<GV>>> boundary_info;
+ std::unique_ptr<XT::Functions::LambdaFunction<d, m>> periodic_density_variation_;
+ std::unique_ptr<typename XT::Functions::ElementFunctionInterface<E, m>> local_periodic_density_variation_;
}; // class InviscidCompressibleFlowEulerExplicitFvTest
@@ -226,12 +352,11 @@ TEST_F(InviscidCompressibleFlow1dEulerExplicitFvTest, impermeable_walls_by_invis
this->boundary_treatment = "impermeable_walls_by_inviscid_mirror_treatment";
this->run();
}
-
-
-#if 0
-// TEST_F(InviscidCompressibleFlowEuler1dExplicitFV,
-// inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right)
-//{
-// this->inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right(10*100);
-//}
-#endif // 0
+TEST_F(InviscidCompressibleFlow1dEulerExplicitFvTest,
+ inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right)
+{
+ this->T_end_ = 2; // We need more time to hit the right wall
+ this->set_numerical_flux("vijayasundaram");
+ this->boundary_treatment = "inflow_from_the_left_by_heuristic_euler_treatment_impermeable_wall_right";
+ this->run();
+}
--
GitLab