Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | List of all members
Nektar::HLLCSolver Class Reference

#include <HLLCSolver.h>

Inheritance diagram for Nektar::HLLCSolver:
Inheritance graph
[legend]
Collaboration diagram for Nektar::HLLCSolver:
Collaboration graph
[legend]

Static Public Member Functions

static RiemannSolverSharedPtr create ()
static RiemannSolverSharedPtr create ()

Static Public Attributes

static std::string solverName

Protected Member Functions

 HLLCSolver ()
virtual void v_PointSolve (NekDouble rhoL, NekDouble rhouL, NekDouble rhovL, NekDouble rhowL, NekDouble EL, NekDouble rhoR, NekDouble rhouR, NekDouble rhovR, NekDouble rhowR, NekDouble ER, NekDouble &rhof, NekDouble &rhouf, NekDouble &rhovf, NekDouble &rhowf, NekDouble &Ef)
 HLLC Riemann solver.
virtual void v_PointSolveVisc (NekDouble rhoL, NekDouble rhouL, NekDouble rhovL, NekDouble rhowL, NekDouble EL, NekDouble EpsL, NekDouble rhoR, NekDouble rhouR, NekDouble rhovR, NekDouble rhowR, NekDouble ER, NekDouble EpsR, NekDouble &rhof, NekDouble &rhouf, NekDouble &rhovf, NekDouble &rhowf, NekDouble &Ef, NekDouble &Epsf)
 HLLCSolver ()
virtual void v_PointSolve (double hL, double huL, double hvL, double hR, double huR, double hvR, double &hf, double &huf, double &hvf)
 HLLC Riemann solver for the Nonlinear Shallow Water Equations.
- Protected Member Functions inherited from Nektar::NonlinearSWESolver
 NonlinearSWESolver ()
virtual void v_Solve (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
virtual void v_ArraySolve (const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
- Protected Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT RiemannSolver ()
void GenerateRotationMatrices (const Array< OneD, const Array< OneD, NekDouble > > &normals)
 Generate rotation matrices for 3D expansions.
void FromToRotation (Array< OneD, const NekDouble > &from, Array< OneD, const NekDouble > &to, NekDouble *mat)
 A function for creating a rotation matrix that rotates a vector from into another vector to.
SOLVER_UTILS_EXPORT void rotateToNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Rotate a vector field to trace normal.
SOLVER_UTILS_EXPORT void rotateFromNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Rotate a vector field from trace normal.
SOLVER_UTILS_EXPORT bool CheckScalars (std::string name)
 Determine whether a scalar has been defined in m_scalars.
SOLVER_UTILS_EXPORT bool CheckVectors (std::string name)
 Determine whether a vector has been defined in m_vectors.
SOLVER_UTILS_EXPORT bool CheckParams (std::string name)
 Determine whether a parameter has been defined in m_params.
SOLVER_UTILS_EXPORT bool CheckAuxScal (std::string name)
 Determine whether a scalar has been defined in m_auxScal.
SOLVER_UTILS_EXPORT bool CheckAuxVec (std::string name)
 Determine whether a vector has been defined in m_auxVec.
- Protected Member Functions inherited from Nektar::CompressibleSolver
 CompressibleSolver ()
virtual void v_Solve (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
virtual void v_ArraySolve (const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT void Solve (const int nDim, const Array< OneD, const Array< OneD, NekDouble > > &Fwd, const Array< OneD, const Array< OneD, NekDouble > > &Bwd, Array< OneD, Array< OneD, NekDouble > > &flux)
 Perform the Riemann solve given the forwards and backwards spaces.
template<typename FuncPointerT , typename ObjectPointerT >
void SetScalar (std::string name, FuncPointerT func, ObjectPointerT obj)
void SetScalar (std::string name, RSScalarFuncType fp)
template<typename FuncPointerT , typename ObjectPointerT >
void SetVector (std::string name, FuncPointerT func, ObjectPointerT obj)
void SetVector (std::string name, RSVecFuncType fp)
template<typename FuncPointerT , typename ObjectPointerT >
void SetParam (std::string name, FuncPointerT func, ObjectPointerT obj)
void SetParam (std::string name, RSParamFuncType fp)
template<typename FuncPointerT , typename ObjectPointerT >
void SetAuxScal (std::string name, FuncPointerT func, ObjectPointerT obj)
template<typename FuncPointerT , typename ObjectPointerT >
void SetAuxVec (std::string name, FuncPointerT func, ObjectPointerT obj)
std::map< std::string,
RSScalarFuncType > & 
GetScalars ()
std::map< std::string,
RSVecFuncType > & 
GetVectors ()
std::map< std::string,
RSParamFuncType > & 
GetParams ()
- Public Attributes inherited from Nektar::SolverUtils::RiemannSolver
int m_spacedim
- Protected Attributes inherited from Nektar::NonlinearSWESolver
bool m_pointSolve
- Protected Attributes inherited from Nektar::CompressibleSolver
bool m_pointSolve

Detailed Description

Definition at line 43 of file CompressibleFlowSolver/RiemannSolvers/HLLCSolver.h.

Constructor & Destructor Documentation

Nektar::HLLCSolver::HLLCSolver ( )
protected

Definition at line 44 of file CompressibleFlowSolver/RiemannSolvers/HLLCSolver.cpp.

Referenced by create().

Nektar::HLLCSolver::HLLCSolver ( )
protected

Member Function Documentation

static RiemannSolverSharedPtr Nektar::HLLCSolver::create ( )
inlinestatic

Definition at line 46 of file CompressibleFlowSolver/RiemannSolvers/HLLCSolver.h.

References HLLCSolver().

{
}
static RiemannSolverSharedPtr Nektar::HLLCSolver::create ( )
inlinestatic

Definition at line 46 of file ShallowWaterSolver/RiemannSolvers/HLLCSolver.h.

References HLLCSolver().

{
new HLLCSolver());
}
void Nektar::HLLCSolver::v_PointSolve ( NekDouble  rhoL,
NekDouble  rhouL,
NekDouble  rhovL,
NekDouble  rhowL,
NekDouble  EL,
NekDouble  rhoR,
NekDouble  rhouR,
NekDouble  rhovR,
NekDouble  rhowR,
NekDouble  ER,
NekDouble rhof,
NekDouble rhouf,
NekDouble rhovf,
NekDouble rhowf,
NekDouble Ef 
)
protectedvirtual

HLLC Riemann solver.

Parameters
rhoLDensity left state.
rhoRDensity right state.
rhouLx-momentum component left state.
rhouRx-momentum component right state.
rhovLy-momentum component left state.
rhovRy-momentum component right state.
rhowLz-momentum component left state.
rhowRz-momentum component right state.
ELEnergy left state.
EREnergy right state.
rhofComputed Riemann flux for density.
rhoufComputed Riemann flux for x-momentum component
rhovfComputed Riemann flux for y-momentum component
rhowfComputed Riemann flux for z-momentum component
EfComputed Riemann flux for energy.

Reimplemented from Nektar::CompressibleSolver.

Definition at line 68 of file CompressibleFlowSolver/RiemannSolvers/HLLCSolver.cpp.

References Nektar::SolverUtils::RiemannSolver::m_params.

{
static NekDouble gamma = m_params["gamma"]();
// Left and Right velocities
NekDouble uL = rhouL / rhoL;
NekDouble vL = rhovL / rhoL;
NekDouble wL = rhowL / rhoL;
NekDouble uR = rhouR / rhoR;
NekDouble vR = rhovR / rhoR;
NekDouble wR = rhowR / rhoR;
// Left and right pressure, sound speed and enthalpy.
NekDouble pL = (gamma - 1.0) *
(EL - 0.5 * (rhouL * uL + rhovL * vL + rhowL * wL));
NekDouble pR = (gamma - 1.0) *
(ER - 0.5 * (rhouR * uR + rhovR * vR + rhowR * wR));
NekDouble cL = sqrt(gamma * pL / rhoL);
NekDouble cR = sqrt(gamma * pR / rhoR);
NekDouble hL = (EL + pL) / rhoL;
NekDouble hR = (ER + pR) / rhoR;
// Square root of rhoL and rhoR.
NekDouble srL = sqrt(rhoL);
NekDouble srR = sqrt(rhoR);
NekDouble srLR = srL + srR;
// Velocity Roe averages
NekDouble uRoe = (srL * uL + srR * uR) / srLR;
NekDouble vRoe = (srL * vL + srR * vR) / srLR;
NekDouble wRoe = (srL * wL + srR * wR) / srLR;
NekDouble hRoe = (srL * hL + srR * hR) / srLR;
NekDouble cRoe = sqrt((gamma - 1.0)*(hRoe - 0.5 *
(uRoe * uRoe + vRoe * vRoe + wRoe * wRoe)));
// Maximum wave speeds
NekDouble SL = std::min(uL-cL, uRoe-cRoe);
NekDouble SR = std::max(uR+cR, uRoe+cRoe);
// HLLC Riemann fluxes (positive case)
if (SL >= 0)
{
rhof = rhouL;
rhouf = rhouL * uL + pL;
rhovf = rhouL * vL;
rhowf = rhouL * wL;
Ef = uL * (EL + pL);
}
// HLLC Riemann fluxes (negative case)
else if (SR <= 0)
{
rhof = rhouR;
rhouf = rhouR * uR + pR;
rhovf = rhouR * vR;
rhowf = rhouR * wR;
Ef = uR * (ER + pR);
}
// HLLC Riemann fluxes (general case (SL < 0 | SR > 0)
else
{
NekDouble SM = (pR - pL + rhouL * (SL - uL) - rhouR * (SR - uR)) /
(rhoL * (SL - uL) - rhoR * (SR - uR));
NekDouble rhoML = rhoL * (SL - uL) / (SL - SM);
NekDouble rhouML = rhoML * SM;
NekDouble rhovML = rhoML * vL;
NekDouble rhowML = rhoML * wL;
NekDouble EML = rhoML * (EL / rhoL +
(SM - uL) * (SM + pL / (rhoL * (SL - uL))));
NekDouble rhoMR = rhoR * (SR - uR) / (SR - SM);
NekDouble rhouMR = rhoMR * SM;
NekDouble rhovMR = rhoMR * vR;
NekDouble rhowMR = rhoMR * wR;
NekDouble EMR = rhoMR * (ER / rhoR +
(SM - uR) * (SM + pR / (rhoR * (SR - uR))));
if (SL < 0.0 && SM >= 0.0)
{
rhof = rhouL + SL * (rhoML - rhoL);
rhouf = rhouL * uL + pL + SL * (rhouML - rhouL);
rhovf = rhouL * vL + SL * (rhovML - rhovL);
rhowf = rhouL * wL + SL * (rhowML - rhowL);
Ef = uL * (EL + pL) + SL * (EML - EL);
}
else if(SM < 0.0 && SR > 0.0)
{
rhof = rhouR + SR * (rhoMR - rhoR);
rhouf = rhouR * uR + pR + SR * (rhouMR - rhouR);
rhovf = rhouR * vR + SR * (rhovMR - rhovR);
rhowf = rhouR * wR + SR * (rhowMR - rhowR);
Ef = uR * (ER + pR) + SR * (EMR - ER);
}
}
}
void Nektar::HLLCSolver::v_PointSolve ( double  hL,
double  huL,
double  hvL,
double  hR,
double  huR,
double  hvR,
double &  hf,
double &  huf,
double &  hvf 
)
protectedvirtual

HLLC Riemann solver for the Nonlinear Shallow Water Equations.

Parameters
hLWater depth left state.
hRWater depth right state.
huLx-momentum component left state.
huRx-momentum component right state.
hvLy-momentum component left state.
hvRy-momentum component right state.
hfComputed Riemann flux for density.
hufComputed Riemann flux for x-momentum component
hvfComputed Riemann flux for y-momentum component

Reimplemented from Nektar::NonlinearSWESolver.

Definition at line 64 of file ShallowWaterSolver/RiemannSolvers/HLLCSolver.cpp.

References ASSERTL0, and Nektar::SolverUtils::RiemannSolver::m_params.

{
static NekDouble g = m_params["gravity"]();
// Left and Right velocities
NekDouble uL = huL / hL;
NekDouble vL = hvL / hL;
NekDouble uR = huR / hR;
NekDouble vR = hvR / hR;
// Left and right wave speeds
NekDouble cL = sqrt(g * hL);
NekDouble cR = sqrt(g * hR);
// the two-rarefaction wave assumption
NekDouble hC,huC,hvC,SL,SR,hstar,ustar,Sstar;
hstar = 0.5*(cL + cR) + 0.25*(uL - uR);
hstar *= hstar;
hstar *= (1.0/g);
ustar = 0.5*(uL + uR) + cL - cR;
// Compute SL
if (hstar > hL)
SL = uL - cL * sqrt(0.5*((hstar*hstar + hstar*hL)/(hL*hL)));
else
SL = uL - cL;
// Compute SR
if (hstar > hR)
SR = uR + cR * sqrt(0.5*((hstar*hstar + hstar*hR)/(hR*hR)));
else
SR = uR + cR;
if (fabs(hR*(uR-SR)-hL*(uL-SL)) <= 1.0e-10)
Sstar = ustar;
else
Sstar = (SL*hR*(uR-SR)-SR*hL*(uL-SL))/(hR*(uR-SR)-hL*(uL-SL));
if (SL >= 0)
{
hf = hL * uL;
huf = uL * uL * hL + 0.5 * g * hL * hL;
hvf = hL * uL * vL;
}
else if (SR <= 0)
{
hf = hR * uR;
huf = uR * uR * hR + 0.5 * g * hR * hR;
hvf = hR * uR *vR;
}
else if ((SL < 0) && (Sstar >= 0))
{
hC = hL * ((SL - uL) / (SL - Sstar));
huC = hC * Sstar;
hvC = hC * vL;
hf = hL*uL + SL * (hC - hL);
huf = (uL*uL*hL+0.5*g*hL*hL) + SL * (huC - hL*uL);
hvf = (uL*vL*hL) + SL * (hvC - hL*vL);
}
else if ((SR > 0) && (Sstar <= 0))
{
hC = hR * ((SR - uR) / (SR - Sstar));
huC = hC * Sstar;
hvC = hC * vR;
hf = hR*uR + SR * (hC - hR);
huf = (uR*uR*hR+0.5*g*hR*hR) + SR * (huC - hR*uR);
hvf = (uR*vR*hR) + SR * (hvC - hR*vR);
}
else
{
ASSERTL0(false,"Error in HLLC solver -- non physical combination of SR, SL and Sstar");
}
}
void Nektar::HLLCSolver::v_PointSolveVisc ( NekDouble  rhoL,
NekDouble  rhouL,
NekDouble  rhovL,
NekDouble  rhowL,
NekDouble  EL,
NekDouble  EpsL,
NekDouble  rhoR,
NekDouble  rhouR,
NekDouble  rhovR,
NekDouble  rhowR,
NekDouble  ER,
NekDouble  EpsR,
NekDouble rhof,
NekDouble rhouf,
NekDouble rhovf,
NekDouble rhowf,
NekDouble Ef,
NekDouble Epsf 
)
protectedvirtual

Reimplemented from Nektar::CompressibleSolver.

Definition at line 166 of file CompressibleFlowSolver/RiemannSolvers/HLLCSolver.cpp.

References Nektar::SolverUtils::RiemannSolver::m_params.

{
static NekDouble gamma = m_params["gamma"]();
// Left and Right velocities
NekDouble uL = rhouL / rhoL;
NekDouble vL = rhovL / rhoL;
NekDouble wL = rhowL / rhoL;
NekDouble uR = rhouR / rhoR;
NekDouble vR = rhovR / rhoR;
NekDouble wR = rhowR / rhoR;
// Left and right pressure, sound speed and enthalpy.
NekDouble pL = (gamma - 1.0) *
(EL - 0.5 * (rhouL * uL + rhovL * vL + rhowL * wL));
NekDouble pR = (gamma - 1.0) *
(ER - 0.5 * (rhouR * uR + rhovR * vR + rhowR * wR));
NekDouble cL = sqrt(gamma * pL / rhoL);
NekDouble cR = sqrt(gamma * pR / rhoR);
NekDouble hL = (EL + pL) / rhoL;
NekDouble hR = (ER + pR) / rhoR;
// Square root of rhoL and rhoR.
NekDouble srL = sqrt(rhoL);
NekDouble srR = sqrt(rhoR);
NekDouble srLR = srL + srR;
// Velocity Roe averages
NekDouble uRoe = (srL * uL + srR * uR) / srLR;
NekDouble vRoe = (srL * vL + srR * vR) / srLR;
NekDouble wRoe = (srL * wL + srR * wR) / srLR;
NekDouble hRoe = (srL * hL + srR * hR) / srLR;
NekDouble cRoe = sqrt((gamma - 1.0)*(hRoe - 0.5 *
(uRoe * uRoe + vRoe * vRoe + wRoe * wRoe)));
// Maximum wave speeds
NekDouble SL = std::min(uL-cL, uRoe-cRoe);
NekDouble SR = std::max(uR+cR, uRoe+cRoe);
// HLLC Riemann fluxes (positive case)
if (SL >= 0)
{
rhof = rhouL;
rhouf = rhouL * uL + pL;
rhovf = rhouL * vL;
rhowf = rhouL * wL;
Ef = uL * (EL + pL);
Epsf = 0.0;
}
// HLLC Riemann fluxes (negative case)
else if (SR <= 0)
{
rhof = rhouR;
rhouf = rhouR * uR + pR;
rhovf = rhouR * vR;
rhowf = rhouR * wR;
Ef = uR * (ER + pR);
Epsf = 0.0;
}
// HLLC Riemann fluxes (general case (SL < 0 | SR > 0)
else
{
NekDouble SM = (pR - pL + rhouL * (SL - uL) - rhouR * (SR - uR)) /
(rhoL * (SL - uL) - rhoR * (SR - uR));
NekDouble rhoML = rhoL * (SL - uL) / (SL - SM);
NekDouble rhouML = rhoML * SM;
NekDouble rhovML = rhoML * vL;
NekDouble rhowML = rhoML * wL;
NekDouble EML = rhoML * (EL / rhoL +
(SM - uL) * (SM + pL / (rhoL * (SL - uL))));
NekDouble EpsML = EpsL * (SL - uL) / (SL - SM);
NekDouble rhoMR = rhoR * (SR - uR) / (SR - SM);
NekDouble rhouMR = rhoMR * SM;
NekDouble rhovMR = rhoMR * vR;
NekDouble rhowMR = rhoMR * wR;
NekDouble EMR = rhoMR * (ER / rhoR +
(SM - uR) * (SM + pR / (rhoR * (SR - uR))));
NekDouble EpsMR = EpsR * (SL - uR) / (SL - SM);
if (SL < 0.0 && SM >= 0.0)
{
rhof = rhouL + SL * (rhoML - rhoL);
rhouf = rhouL * uL + pL + SL * (rhouML - rhouL);
rhovf = rhouL * vL + SL * (rhovML - rhovL);
rhowf = rhouL * wL + SL * (rhowML - rhowL);
Ef = uL * (EL + pL) + SL * (EML - EL);
Epsf = 0.0 + SL * (EpsML - EpsL);
}
else if(SM < 0.0 && SR > 0.0)
{
rhof = rhouR + SR * (rhoMR - rhoR);
rhouf = rhouR * uR + pR + SR * (rhouMR - rhouR);
rhovf = rhouR * vR + SR * (rhovMR - rhovR);
rhowf = rhouR * wR + SR * (rhowMR - rhowR);
Ef = uR * (ER + pR) + SR * (EMR - ER);
Epsf = 0.0 + SR * (EpsMR - EpsR);
}
}
}

Member Data Documentation

static std::string Nektar::HLLCSolver::solverName
static