Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::UpwindPulseSolver Class Reference

#include <UpwindPulseSolver.h>

Inheritance diagram for Nektar::UpwindPulseSolver:
[legend]

Static Public Member Functions

static RiemannSolverSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession)
 

Static Public Attributes

static std::string solverName
 

Protected Member Functions

 UpwindPulseSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 
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)
 
void RiemannSolverUpwind (NekDouble AL, NekDouble uL, NekDouble AR, NekDouble uR, NekDouble &Aflux, NekDouble &uflux, NekDouble A0, NekDouble beta, NekDouble n, NekDouble alpha=0.5)
 
- Protected Member Functions inherited from Nektar::SolverUtils::RiemannSolver
SOLVER_UTILS_EXPORT RiemannSolver ()
 
SOLVER_UTILS_EXPORT RiemannSolver (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual SOLVER_UTILS_EXPORT ~RiemannSolver ()
 
SOLVER_UTILS_EXPORT void GenerateRotationMatrices (const Array< OneD, const Array< OneD, NekDouble >> &normals)
 Generate rotation matrices for 3D expansions. More...
 
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. More...
 
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. More...
 
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. More...
 
SOLVER_UTILS_EXPORT bool CheckScalars (std::string name)
 Determine whether a scalar has been defined in m_scalars. More...
 
SOLVER_UTILS_EXPORT bool CheckVectors (std::string name)
 Determine whether a vector has been defined in m_vectors. More...
 
SOLVER_UTILS_EXPORT bool CheckParams (std::string name)
 Determine whether a parameter has been defined in m_params. More...
 
SOLVER_UTILS_EXPORT bool CheckAuxScal (std::string name)
 Determine whether a scalar has been defined in m_auxScal. More...
 
SOLVER_UTILS_EXPORT bool CheckAuxVec (std::string name)
 Determine whether a vector has been defined in m_auxVec. More...
 
virtual SOLVER_UTILS_EXPORT void v_CalcFluxJacobian (const int nDim, const Array< OneD, const Array< OneD, NekDouble >> &Fwd, const Array< OneD, const Array< OneD, NekDouble >> &Bwd, const Array< OneD, const Array< OneD, NekDouble >> &normals, DNekBlkMatSharedPtr &FJac, DNekBlkMatSharedPtr &BJac)
 

Protected Attributes

LibUtilities::SessionReaderSharedPtr m_session
 
int m_nVariables
 
Array< OneD, MultiRegions::ExpListSharedPtrm_vessels
 
PulseWavePressureAreaSharedPtr m_pressureArea
 
- Protected Attributes inherited from Nektar::SolverUtils::RiemannSolver
bool m_requiresRotation
 Indicates whether the Riemann solver requires a rotation to be applied to the velocity fields. More...
 
std::map< std::string, RSScalarFuncTypem_scalars
 Map of scalar function types. More...
 
std::map< std::string, RSVecFuncTypem_vectors
 Map of vector function types. More...
 
std::map< std::string, RSParamFuncTypem_params
 Map of parameter function types. More...
 
std::map< std::string, RSScalarFuncTypem_auxScal
 Map of auxiliary scalar function types. More...
 
std::map< std::string, RSVecFuncTypem_auxVec
 Map of auxiliary vector function types. More...
 
Array< OneD, Array< OneD, NekDouble > > m_rotMat
 Rotation matrices for each trace quadrature point. More...
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_rotStorage
 Rotation storage. More...
 

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. More...
 
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)
 
void SetAuxVec (std::string name, RSVecFuncType fp)
 
std::map< std::string, RSScalarFuncType > & GetScalars ()
 
std::map< std::string, RSVecFuncType > & GetVectors ()
 
std::map< std::string, RSParamFuncType > & GetParams ()
 
SOLVER_UTILS_EXPORT void CalcFluxJacobian (const int nDim, const Array< OneD, const Array< OneD, NekDouble >> &Fwd, const Array< OneD, const Array< OneD, NekDouble >> &Bwd, DNekBlkMatSharedPtr &FJac, DNekBlkMatSharedPtr &BJac)
 Calculate the flux jacobian of Fwd and Bwd. More...
 
- Public Attributes inherited from Nektar::SolverUtils::RiemannSolver
int m_spacedim
 

Detailed Description

Definition at line 47 of file UpwindPulseSolver.h.

Constructor & Destructor Documentation

◆ UpwindPulseSolver()

Nektar::UpwindPulseSolver::UpwindPulseSolver ( const LibUtilities::SessionReaderSharedPtr pSession)
protected

Definition at line 44 of file UpwindPulseSolver.cpp.

46  : RiemannSolver(pSession), m_session(pSession)
47 {
48 }
SOLVER_UTILS_EXPORT RiemannSolver()
LibUtilities::SessionReaderSharedPtr m_session

Member Function Documentation

◆ create()

static RiemannSolverSharedPtr Nektar::UpwindPulseSolver::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Definition at line 50 of file UpwindPulseSolver.h.

52  {
53  return RiemannSolverSharedPtr(new UpwindPulseSolver(pSession));
54  }
UpwindPulseSolver(const LibUtilities::SessionReaderSharedPtr &pSession)
std::shared_ptr< RiemannSolver > RiemannSolverSharedPtr
A shared pointer to an EquationSystem object.

◆ RiemannSolverUpwind()

void Nektar::UpwindPulseSolver::RiemannSolverUpwind ( NekDouble  AL,
NekDouble  uL,
NekDouble  AR,
NekDouble  uR,
NekDouble Aflux,
NekDouble uflux,
NekDouble  A0,
NekDouble  beta,
NekDouble  n,
NekDouble  alpha = 0.5 
)
protected

Riemann solver for upwinding at an interface between two elements. Uses the characteristic variables for calculating the upwinded state \((A_u, u_u)\) from the left \((A_L, u_L)\) and right state \((A_R, u_R)\). Returns the upwinded flux $\mathbf{F}^u$ needed for the weak formulation (1). Details can be found in "Pulse wave propagation in the human vascular system", section 3.3

Definition at line 93 of file UpwindPulseSolver.cpp.

98 {
99  NekDouble W1 = 0.0;
100  NekDouble W2 = 0.0;
101  NekDouble IL = 0.0;
102  NekDouble IR = 0.0;
103  NekDouble Au = 0.0;
104  NekDouble uu = 0.0;
105  NekDouble cL = 0.0;
106  NekDouble cR = 0.0;
107  NekDouble P = 0.0;
108 
109  ASSERTL1(CheckParams("rho"), "rho not defined.");
110  NekDouble rho = m_params["rho"]();
111  NekDouble nDomains = m_params["domains"]();
112 
113  m_nVariables = m_session->GetVariables().size();
114 
115  m_vessels =
116  Array<OneD, MultiRegions::ExpListSharedPtr>(m_nVariables * nDomains);
117 
118  if (m_session->DefinesSolverInfo("PressureArea"))
119  {
121  m_session->GetSolverInfo("PressureArea"), m_vessels, m_session);
122  }
123  else
124  {
126  "Beta", m_vessels, m_session);
127  }
128 
129  // Compute the wave speeds to check dynamics are sub-sonic
130  m_pressureArea->GetC(cL, beta, AL, A0, alpha);
131  m_pressureArea->GetC(cR, beta, AR, A0, alpha);
132  ASSERTL1(fabs(cL + cR) > fabs(uL + uR), "Conditions are not sub-sonic");
133 
134  /*
135  Calculate the characteristics. The use of the normal here allows
136  for the definition of the characteristics (and hence the left
137  and right state) to be inverted if n is in the -ve
138  x-direction. This means we end up with the positive
139  defintion of the flux which has to therefore be multiplied
140  by the normal at the end of the method. This is a bit of a
141  mind twister but is efficient from a coding perspective.
142  */
143  m_pressureArea->GetCharIntegral(IL, beta, AL, A0, alpha);
144  m_pressureArea->GetCharIntegral(IR, beta, AR, A0, alpha);
145  W1 = uL + n * IL;
146  W2 = uR - n * IR;
147 
148  // Calculate conservative variables from characteristics
149  m_pressureArea->GetAFromChars(Au, n * W1, n * W2, beta, A0, alpha);
150  m_pressureArea->GetUFromChars(uu, W1, W2);
151 
152  // Pressure for the energy flux
153  m_pressureArea->GetPressure(P, beta, Au, A0, 0, 0, alpha);
154 
155  // Compute the fluxes multiplied by the normal
156  Aflux = Au * uu * n;
157  uflux = (uu * uu / 2 + P / rho) * n;
158 }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
SOLVER_UTILS_EXPORT bool CheckParams(std::string name)
Determine whether a parameter has been defined in m_params.
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
PulseWavePressureAreaSharedPtr m_pressureArea
Array< OneD, MultiRegions::ExpListSharedPtr > m_vessels
@ beta
Gauss Radau pinned at x=-1,.
Definition: PointsType.h:61
PressureAreaFactory & GetPressureAreaFactory()
double NekDouble

References ASSERTL1, Nektar::LibUtilities::beta, Nektar::SolverUtils::RiemannSolver::CheckParams(), Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::GetPressureAreaFactory(), m_nVariables, Nektar::SolverUtils::RiemannSolver::m_params, m_pressureArea, m_session, m_vessels, and Nektar::LibUtilities::P.

Referenced by v_Solve().

◆ v_Solve()

void Nektar::UpwindPulseSolver::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 
)
protectedvirtual

Calculates the third term of the weak form (1): numerical flux at boundary \( \left[ \mathbf{\psi}^{\delta} \cdot \{ \mathbf{F}^u - \mathbf{F}(\mathbf{U}^{\delta}) \} \right]_{x_e^l}^{x_eû} \)

Implements Nektar::SolverUtils::RiemannSolver.

Definition at line 56 of file UpwindPulseSolver.cpp.

60 {
61  int i;
62  int nTracePts = Fwd[0].size();
63 
64  ASSERTL1(CheckScalars("A0"), "A0 not defined.");
65  const Array<OneD, NekDouble> &A0 = m_scalars["A0"]();
66 
67  ASSERTL1(CheckScalars("beta"), "beta not defined.");
68  const Array<OneD, NekDouble> &beta = m_scalars["beta"]();
69 
70  const Array<OneD, NekDouble> &alpha = m_scalars["alpha"]();
71 
72  ASSERTL1(CheckScalars("N"), "N not defined.");
73  const Array<OneD, NekDouble> &N = m_scalars["N"]();
74 
75  for (i = 0; i < nTracePts; ++i)
76  {
77  RiemannSolverUpwind(Fwd[0][i], Fwd[1][i], Bwd[0][i], Bwd[1][i],
78  flux[0][i], flux[1][i], A0[i], beta[i], N[i],
79  alpha[i]);
80  }
81 }
SOLVER_UTILS_EXPORT bool CheckScalars(std::string name)
Determine whether a scalar has been defined in m_scalars.
std::map< std::string, RSScalarFuncType > m_scalars
Map of scalar function types.
void RiemannSolverUpwind(NekDouble AL, NekDouble uL, NekDouble AR, NekDouble uR, NekDouble &Aflux, NekDouble &uflux, NekDouble A0, NekDouble beta, NekDouble n, NekDouble alpha=0.5)

References ASSERTL1, Nektar::LibUtilities::beta, Nektar::SolverUtils::RiemannSolver::CheckScalars(), Nektar::SolverUtils::RiemannSolver::m_scalars, and RiemannSolverUpwind().

Member Data Documentation

◆ m_nVariables

int Nektar::UpwindPulseSolver::m_nVariables
protected

Definition at line 62 of file UpwindPulseSolver.h.

Referenced by RiemannSolverUpwind().

◆ m_pressureArea

PulseWavePressureAreaSharedPtr Nektar::UpwindPulseSolver::m_pressureArea
protected

Definition at line 64 of file UpwindPulseSolver.h.

Referenced by RiemannSolverUpwind().

◆ m_session

LibUtilities::SessionReaderSharedPtr Nektar::UpwindPulseSolver::m_session
protected

Definition at line 61 of file UpwindPulseSolver.h.

Referenced by RiemannSolverUpwind().

◆ m_vessels

Array<OneD, MultiRegions::ExpListSharedPtr> Nektar::UpwindPulseSolver::m_vessels
protected

Definition at line 63 of file UpwindPulseSolver.h.

Referenced by RiemannSolverUpwind().

◆ solverName

std::string Nektar::UpwindPulseSolver::solverName
static
Initial value:
=
"UpwindPulse", UpwindPulseSolver::create, "UpwindPulseSolver")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
static RiemannSolverSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession)
RiemannSolverFactory & GetRiemannSolverFactory()

Definition at line 56 of file UpwindPulseSolver.h.