Skip to content
Snippets Groups Projects
Commit 3514d132 authored by Robert K's avatar Robert K
Browse files

[feature][hllc-bgfix] Added HLLC with bgfix version.

parent b0fb1e41
No related branches found
No related tags found
No related merge requests found
Pipeline #55654 passed
Stage: test
Hide whitespace changes
Inline Side-by-side
dune/fem-dg/operator/fluxes/advection/parameters.hh
+ 8
− 2
View file @ 3514d132
......@@ -42,9 +42,13 @@ namespace Fem
euler_hll_pt,
//! the HLLC flux
euler_hllc,
//! the HLLC flux with potential temperature
euler_hllc_pt,
//! the Harten, Lax and van Leer (HLL) flux with potential temperature formulation
//! and background fix
euler_hll_bgfix,
//! the HLLC with potential temperature and background fix
euler_hllc_bgfix,
//! general flux: Parameter selection is done via parameter file!
euler_general,
......@@ -60,12 +64,14 @@ namespace Fem
//! Contains all known enums for advection fluxes which can be chosen via parameter file.
const Enum _enums[] = { Enum::none, Enum::upwind, Enum::llf,
Enum::euler_llf, Enum::euler_hll, Enum::euler_hll_pt, Enum::euler_hllc,
Enum::euler_llf, Enum::euler_hll, Enum::euler_hll_pt, Enum::euler_hllc, euler_hllc_pt,
euler_hll_bgfix, euler_hllc_bgfix,
Enum::mhd_dw, Enum::mhd_hllem
};
//! Contains all known names of advection fluxes which can be chosen via parameter file.
const std::string _strings[] = { "NONE", "UPWIND" , "LLF",
"EULER-LLF", "EULER-HLL" , "EULER-HLL-PT", "EULER-HLLC",
"EULER-LLF", "EULER-HLL" , "EULER-HLL-PT", "EULER-HLLC", "EULER-HLLC-PT",
"EULER-HLL-BGFIX", "EULER-HLLC-BGFIX",
"MHD-DW", "MHD-HLLEM"
};
//! Number of known advection fluxes which can be chosen via parameter file.
......
dune/fem-dg/operator/fluxes/euler/CMakeLists.txt
+ 1
− 1
View file @ 3514d132
dune_install( eulerflux_impl.hh fluxes.hh llfadv.hh parameters.hh )
dune_install( bgfixflux.hh eulerflux_impl.hh fluxes.hh llfadv.hh parameters.hh )
dune/fem-dg/operator/fluxes/euler/bgfixflux.hh 0 → 100644
+ 102
− 0
View file @ 3514d132
#ifndef DUNE_FEM_DG_EULER_BGFIX_FLUX_HH
#define DUNE_FEM_DG_EULER_BGFIX_FLUX_HH
// system includes
#include <string>
#include <cmath>
#include "../advection/fluxbase.hh"
#include "../advection/fluxes.hh"
namespace Dune
{
namespace Fem
{
/**
* \brief class specialization for a general flux chosen by a parameter file.
*
* The purpose of this class is to allow the selection of an Euler flux
* via an enum given in AdvectionFlux::Enum.
*/
template< class ModelImp, AdvectionFlux::Enum id >
class DGAdvectionFluxBGFix : public DGAdvectionFluxBase< ModelImp, AdvectionFluxParameters >
{
typedef DGAdvectionFluxBase< ModelImp, AdvectionFluxParameters > BaseType;
typedef DGAdvectionFlux< ModelImp, id > NumFluxImplType;
static const int dimRange = ModelImp::dimRange;
typedef typename ModelImp::DomainType DomainType;
typedef typename ModelImp::RangeType RangeType;
typedef typename ModelImp::JacobianRangeType JacobianRangeType;
typedef typename ModelImp::FluxRangeType FluxRangeType;
typedef typename ModelImp::FaceDomainType FaceDomainType;
public:
typedef AdvectionFlux::Enum IdEnum;
typedef typename BaseType::ModelType ModelType;
typedef typename BaseType::ParameterType ParameterType;
using BaseType::model;
/**
* \copydoc DGAdvectionFluxBase::DGAdvectionFluxBase()
*/
template< class ... Args>
DGAdvectionFluxBGFix( Args&&... args )
: BaseType( std::forward<Args>(args)... ),
numFlux_( std::forward<Args>(args)... )
{}
/**
* \copydoc DGAdvectionFluxBase::name()
*/
static std::string name () {
return std::string("(BGFix) + FluxImpl");
}
/**
* \copydoc DGAdvectionFluxBase::numericalFlux()
*/
template< class LocalEvaluation >
inline double
numericalFlux( const LocalEvaluation& left,
const LocalEvaluation& right,
const RangeType& uLeft,
const RangeType& uRight,
const JacobianRangeType& jacLeft,
const JacobianRangeType& jacRight,
RangeType& gLeft,
RangeType& gRight) const
{
// Assumption: background atmosphere is continuous
auto bgLeft = model().background( left.entity(), left.quadraturePoint() );
auto bgRight = model().background( right.entity(), right.quadraturePoint() );
RangeType uLeftTot( uLeft );
RangeType uRightTot( uRight );
uLeftTot += bgLeft;
uRightTot += bgRight;
RangeType gBgLeft, gBgRight;
// numerical flux
const double ws = numFlux_.numericalFlux(left, right, uLeftTot, uRightTot,
jacLeft, jacRight, gLeft, gRight );
// background flux
numFlux_.numericalFlux(left, right, bgLeft, bgRight, jacLeft, jacRight, gBgLeft, gBgRight );
gLeft -= gBgLeft;
gRight -= gBgRight;
return ws;
}
protected:
NumFluxImplType numFlux_;
};
}
}
#endif // file declaration
dune/fem-dg/operator/fluxes/euler/eulerflux_impl.hh
+ 54
− 11
View file @ 3514d132
......@@ -24,7 +24,7 @@ namespace EulerNumFlux
//
////////////////////////////////////////////////////////
typedef enum {LLF, HLL, HLL_PT, HLLC} EulerFluxType;
typedef enum {LLF, HLL, HLL_PT, HLLC, HLLC_PT} EulerFluxType;
const EulerFluxType std_flux_type = HLL;
/**
......@@ -81,7 +81,10 @@ namespace EulerNumFlux
// potential temperature
if constexpr (flux_type == HLL_PT) return num_flux_HLL(Uj, Un, normal, gj, std::true_type());
if constexpr (flux_type == HLLC) return num_flux_HLLC(Uj, Un, normal, gj);
if constexpr (flux_type == HLLC) return num_flux_HLLC(Uj, Un, normal, gj, std::false_type());
// potential temperature version
if constexpr (flux_type == HLLC_PT) return num_flux_HLLC(Uj, Un, normal, gj, std::true_type());
DUNE_THROW( Dune::NotImplemented, "Numerical flux not implemented" );
}
......@@ -284,13 +287,37 @@ namespace EulerNumFlux
}
template <bool pottemp>
FieldType num_flux_HLLC(const FieldType Um[dim+2], const FieldType Up[dim+2],
const FieldType normal[dim], FieldType g[dim+2]) const
const FieldType normal[dim], FieldType g[dim+2],
const std::integral_constant<bool, pottemp> ) const
{
static constexpr bool potentialTemperature = pottemp;
const FieldType rhom = Um[0];
const FieldType rhop = Up[0];
const FieldType Em = Um[dim+1];
const FieldType Ep = Up[dim+1];
FieldType Em = Um[dim+1];
FieldType Ep = Up[dim+1];
// compute Em and Ep differently
if constexpr (potentialTemperature)
{
// j==m and n==p
FieldType pressj, tempj;
FieldType pressn, tempn;
model_.pressureTemperature( Um, pressj, tempj );
model_.pressureTemperature( Up, pressn, tempn );
FieldType Ekinj = 0;
FieldType Ekinn = 0;
for(int i=1; i<dim+1; i++)
{
Ekinj += (0.5/rhom) * Um[i] * Um[i];
Ekinn += (0.5/rhop) * Up[i] * Up[i];
}
Em = pressj/(_gamma-1.) + Ekinj;
Ep = pressn/(_gamma-1.) + Ekinn;
}
FieldType rho_um[dim], rho_up[dim];
rotate( normal, Um+1, rho_um );
......@@ -336,7 +363,11 @@ namespace EulerNumFlux
guj[i] = rho_um[i]*um[0];
guj[0] += pm;
g[dim+1] = (Em+pm)*um[0];
//g[dim+1] = (Em+pm)*um[0];
if constexpr (potentialTemperature)
g[dim+1] = Um[dim+1]*um[0];
else
g[dim+1] = (Em+pm)*um[0];
}
else if (sp <= 0.0)
{
......@@ -346,7 +377,11 @@ namespace EulerNumFlux
guj[i] = rho_up[i]*up[0];
guj[0] += pp;
g[dim+1] = (Ep+pp)*up[0];
//g[dim+1] = (Ep+pp)*up[0];
if constexpr (potentialTemperature)
g[dim+1] = Up[dim+1]*up[0];
else
g[dim+1] = (Ep+pp)*up[0];
}
else
{
......@@ -361,8 +396,12 @@ namespace EulerNumFlux
guj[i] = rho_um[i]*um[0] + rhom*um[i]*tmpm - sm*rho_um[i];
guj[0] += pm + rhom*(u_star-um[0])*tmpm;
g[dim+1] = (Em+pm)*um[0] + Em*(tmpm-sm)
+ tmpm*(u_star-um[0])*( rhom*u_star + pm/(sm-um[0]) );
if constexpr(potentialTemperature)
g[dim+1] = (Um[dim+1]*um[0]) + Um[dim+1]*(tmpm-sm);
else
g[dim+1] = (Em+pm)*um[0] + Em*(tmpm-sm)
+ tmpm*(u_star-um[0])*( rhom*u_star + pm/(sm-um[0]) );
}
else
{
......@@ -372,8 +411,12 @@ namespace EulerNumFlux
guj[i] = rho_up[i]*up[0] + rhop*up[i]*tmpp - sp*rho_up[i];
guj[0] += pp + rhop*(u_star-up[0])*tmpp;
g[dim+1] = (Ep+pp)*up[0] + Ep*(tmpp-sp)
+ tmpp*(u_star-up[0])*( rhop*u_star + pp/(sp-up[0]) );
if constexpr(potentialTemperature)
g[dim+1] = (Up[dim+1]*up[0]) + Up[dim+1]*(tmpp-sp);
else
g[dim+1] = (Ep+pp)*up[0] + Ep*(tmpp-sp)
+ tmpp*(u_star-up[0])*( rhop*u_star + pp/(sp-up[0]) );
}
}
......
dune/fem-dg/operator/fluxes/euler/fluxes.hh
+ 54
− 77
View file @ 3514d132
......@@ -7,6 +7,8 @@
#include "../advection/fluxbase.hh"
#include "../advection/fluxes.hh"
#include "bgfixflux.hh"
#include "eulerflux_impl.hh"
#include "llfadv.hh"
......@@ -31,7 +33,7 @@ namespace Fem
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
static std::string name () { return "LLF (Dennis)"; }
static std::string name () { return "LLF (Euler)"; }
};
/**
......@@ -50,7 +52,7 @@ namespace Fem
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
static std::string name () { return "HLL (Dennis)"; }
static std::string name () { return "HLL (Euler)"; }
};
/**
......@@ -69,7 +71,7 @@ namespace Fem
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
static std::string name () { return "HLL_PT (Dennis)"; }
static std::string name () { return "HLL_PT (Euler)"; }
};
/**
......@@ -88,7 +90,27 @@ namespace Fem
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
static std::string name () { return "HLLC (Dennis)"; }
static std::string name () { return "HLLC (Euler)"; }
};
/**
* \brief class specialization for the HLLC flux in potential temperature
* formulation.
*
* The purpose of this class is to allow the selection of an Euler flux
* via an enum given in AdvectionFlux::Enum.
*/
template< class Model >
class DGAdvectionFlux< Model, AdvectionFlux::Enum::euler_hllc_pt >
: public EulerFluxImpl< Model, EulerNumFlux::EulerFlux<Model,EulerNumFlux::EulerFluxType::HLLC_PT > >
{
typedef EulerFluxImpl< Model, EulerNumFlux::EulerFlux<Model,EulerNumFlux::EulerFluxType::HLLC_PT > > BaseType ;
public:
template< class ... Args>
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
static std::string name () { return "HLLC_PT"; }
};
/**
......@@ -177,89 +199,44 @@ namespace Fem
/**
* \brief class specialization for a general flux chosen by a parameter file.
* \brief class specialization for Euler hll_pt and bgfix.
*
* The purpose of this class is to allow the selection of an Euler flux
* via an enum given in AdvectionFlux::Enum.
*/
template< class ModelImp >
class DGAdvectionFlux< ModelImp, AdvectionFlux::Enum::euler_hll_bgfix >
: public DGAdvectionFluxBase< ModelImp, AdvectionFluxParameters >
template< class Model >
class DGAdvectionFlux< Model, AdvectionFlux::Enum::euler_hll_bgfix >
: public DGAdvectionFluxBGFix< Model, AdvectionFlux::Enum::euler_hll_pt > // pass id of implementation
{
typedef DGAdvectionFluxBase< ModelImp, AdvectionFluxParameters > BaseType;
typedef DGAdvectionFlux< ModelImp, AdvectionFlux::Enum::euler_hll_pt > NumFluxImplType;
static const int dimRange = ModelImp::dimRange;
typedef typename ModelImp::DomainType DomainType;
typedef typename ModelImp::RangeType RangeType;
typedef typename ModelImp::JacobianRangeType JacobianRangeType;
typedef typename ModelImp::FluxRangeType FluxRangeType;
typedef typename ModelImp::FaceDomainType FaceDomainType;
typedef DGAdvectionFluxBGFix< Model, AdvectionFlux::Enum::euler_hll_pt > BaseType;
public:
typedef AdvectionFlux::Enum IdEnum;
typedef typename BaseType::ModelType ModelType;
typedef typename BaseType::ParameterType ParameterType;
using BaseType::model;
/**
* \copydoc DGAdvectionFluxBase::DGAdvectionFluxBase()
*/
template< class ... Args>
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... ),
numFlux_( std::forward<Args>(args)... )
: BaseType( std::forward<Args>(args)... )
{}
using BaseType::name;
};
/**
* \copydoc DGAdvectionFluxBase::name()
*/
static std::string name () { return "AdvectionFlux - Euler (via parameter file)"; }
/**
* \copydoc DGAdvectionFluxBase::numericalFlux()
*/
template< class LocalEvaluation >
inline double
numericalFlux( const LocalEvaluation& left,
const LocalEvaluation& right,
const RangeType& uLeft,
const RangeType& uRight,
const JacobianRangeType& jacLeft,
const JacobianRangeType& jacRight,
RangeType& gLeft,
RangeType& gRight) const
{
// Assumption: background atmosphere is continuous
auto bgLeft = model().background( left.entity(), left.quadraturePoint() );
auto bgRight = model().background( right.entity(), right.quadraturePoint() );
RangeType uLeftTot( uLeft );
RangeType uRightTot( uRight );
uLeftTot += bgLeft;
uRightTot += bgRight;
RangeType gBgLeft, gBgRight;
// numerical flux
const double ws = numFlux_.numericalFlux(left, right, uLeftTot, uRightTot,
jacLeft, jacRight, gLeft, gRight );
// background flux
numFlux_.numericalFlux(left, right, bgLeft, bgRight, jacLeft, jacRight, gBgLeft, gBgRight );
gLeft -= gBgLeft;
gRight -= gBgRight;
return ws;
}
protected:
NumFluxImplType numFlux_;
/**
* \brief class specialization for Euler hllc_pt and bgfix.
*
* The purpose of this class is to allow the selection of an Euler flux
* via an enum given in AdvectionFlux::Enum.
*/
template< class Model >
class DGAdvectionFlux< Model, AdvectionFlux::Enum::euler_hllc_bgfix >
: public DGAdvectionFluxBGFix< Model, AdvectionFlux::Enum::euler_hllc_pt > // pass id of implementation
{
typedef DGAdvectionFluxBGFix< Model, AdvectionFlux::Enum::euler_hllc_pt > BaseType;
public:
template< class ... Args>
DGAdvectionFlux( Args&&... args )
: BaseType( std::forward<Args>(args)... )
{}
using BaseType::name;
};
}
}
} // end namespace Fem
} // end namespace Dune
#endif // file declaration
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment