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::HLLSolver Class Reference

#include <HLLSolver.h>

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

Static Public Member Functions

static RiemannSolverSharedPtr create ()
static RiemannSolverSharedPtr create ()

Static Public Attributes

static std::string solverName

Protected Member Functions

 HLLSolver ()
virtual void v_PointSolve (double rhoL, double rhouL, double rhovL, double rhowL, double EL, double rhoR, double rhouR, double rhovR, double rhowR, double ER, double &rhof, double &rhouf, double &rhovf, double &rhowf, double &Ef)
 HLL Riemann solver.
 HLLSolver ()
virtual void v_PointSolve (double hL, double huL, double hvL, double hR, double huR, double hvR, double &hf, double &huf, double &hvf)
 HLL 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)
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)

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/HLLSolver.h.

Constructor & Destructor Documentation

Nektar::HLLSolver::HLLSolver ( )
protected

Definition at line 46 of file CompressibleFlowSolver/RiemannSolvers/HLLSolver.cpp.

Referenced by create().

Nektar::HLLSolver::HLLSolver ( )
protected

Member Function Documentation

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

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

References HLLSolver().

{
new HLLSolver());
}
static RiemannSolverSharedPtr Nektar::HLLSolver::create ( )
inlinestatic

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

References HLLSolver().

{
new HLLSolver());
}
void Nektar::HLLSolver::v_PointSolve ( double  hL,
double  huL,
double  hvL,
double  hR,
double  huR,
double  hvR,
double &  hf,
double &  huf,
double &  hvf 
)
protectedvirtual

HLL 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/HLLSolver.cpp.

References 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 hstar,fL,fR;
hstar = 0.5 * (cL + cR) + 0.25 * (uL - uR);
hstar *= hstar;
hstar *= (1.0/g);
// 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 (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
{
hf = (SR * hL * uL - SL * hR * uR +
SL * SR * (hR - hL)) / (SR - SL);
fL = uL * uL * hL + 0.5 * g * hL * hL;
fR = uR * uR * hR + 0.5 * g * hR * hR;
huf =(SR * fL - SL * fR +
SL * SR * (hR * uR - hL * uL)) / (SR - SL);
fL = uL * vL * hL;
fR = uR * vR * hR;
hvf =(SR * fL - SL * fR +
SL * SR * (hR * vR - hL * vL)) / (SR - SL);
}
}
void Nektar::HLLSolver::v_PointSolve ( double  rhoL,
double  rhouL,
double  rhovL,
double  rhowL,
double  EL,
double  rhoR,
double  rhouR,
double  rhovR,
double  rhowR,
double  ER,
double &  rhof,
double &  rhouf,
double &  rhovf,
double &  rhowf,
double &  Ef 
)
protectedvirtual

HLL 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 70 of file CompressibleFlowSolver/RiemannSolvers/HLLSolver.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 pressures
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);
// HLL Riemann fluxes (positive case)
if (SL >= 0)
{
rhof = rhouL;
rhouf = rhouL * uL + pL;
rhovf = rhouL * vL;
rhowf = rhouL * wL;
Ef = uL * (EL + pL);
}
// HLL Riemann fluxes (negative case)
else if (SR <= 0)
{
rhof = rhouR;
rhouf = rhouR * uR + pR;
rhovf = rhouR * vR;
rhowf = rhouR * wR;
Ef = uR * (ER + pR);
}
// HLL Riemann fluxes (general case (SL < 0 | SR > 0)
else
{
NekDouble tmp1 = 1.0 / (SR - SL);
NekDouble tmp2 = SR * SL;
rhof = (SR * rhouL - SL * rhouR + tmp2 * (rhoR - rhoL)) * tmp1;
rhouf = (SR * (rhouL * uL + pL) - SL * (rhouR * uR + pR) +
tmp2 * (rhouR - rhouL)) * tmp1;
rhovf = (SR * rhouL * vL - SL * rhouR * vR +
tmp2 * (rhovR - rhovL)) * tmp1;
rhowf = (SR * rhouL * wL - SL * rhouR * wR +
tmp2 * (rhowR - rhowL)) * tmp1;
Ef = (SR * uL * (EL + pL) - SL * uR * (ER + pR) +
tmp2 * (ER - EL)) * tmp1;
}
}

Member Data Documentation

static std::string Nektar::HLLSolver::solverName
static