85 double rhowL,
double EL,
double rhoR,
double rhouR,
86 double rhovR,
double rhowR,
double ER,
87 double &rhof,
double &rhouf,
double &rhovf,
88 double &rhowf,
double &Ef)
92 RoeKernel(rhoL, rhouL, rhovL, rhowL, EL, rhoR, rhouR, rhovR, rhowR, ER,
93 rhof, rhouf, rhovf, rhowf, Ef, gamma);
104 static auto gamma =
m_params[
"gamma"]();
105 static size_t nVars = fwd.size();
106 static size_t spaceDim = nVars - 2;
112 size_t sizeScalar = fwd[0].size();
113 size_t sizeVec = (sizeScalar / vec_t::width) * vec_t::width;
117 for (; i < sizeVec; i += vec_t::width)
119 vec_t rhoL{}, rhouL{}, rhovL{}, rhowL{}, EL{};
120 vec_t rhoR{}, rhouR{}, rhovR{}, rhowR{}, ER{};
135 else if (spaceDim == 3)
143 vec_t rhof{}, rhouf{}, rhovf{}, rhowf{}, Ef{};
145 RoeKernel(rhoL, rhouL, rhovL, rhowL, EL, rhoR, rhouR, rhovR, rhowR, ER,
146 rhof, rhouf, rhovf, rhowf, Ef, gamma);
156 else if (spaceDim == 3)
165 for (; i < sizeScalar; ++i)
167 NekDouble rhoL{}, rhouL{}, rhovL{}, rhowL{}, EL{};
168 NekDouble rhoR{}, rhouR{}, rhovR{}, rhowR{}, ER{};
173 EL = fwd[spaceDim + 1][i];
176 ER = bwd[spaceDim + 1][i];
183 else if (spaceDim == 3)
191 NekDouble rhof{}, rhouf{}, rhovf{}, rhowf{}, Ef{};
193 RoeKernel(rhoL, rhouL, rhovL, rhowL, EL, rhoR, rhouR, rhovR, rhowR, ER,
194 rhof, rhouf, rhovf, rhowf, Ef, gamma);
199 flux[nVars - 1][i] = Ef;
204 else if (spaceDim == 3)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
void v_PointSolve(ND rhoL, ND rhouL, ND rhovL, ND rhowL, ND EL, ND rhoR, ND rhouR, ND rhovR, ND rhowR, ND ER, ND &rhof, ND &rhouf, ND &rhovf, ND &rhowf, ND &Ef) final
Roe Riemann solver.
void v_ArraySolve(const Array< OneD, const Array< OneD, ND > > &Fwd, const Array< OneD, const Array< OneD, ND > > &Bwd, Array< OneD, Array< OneD, ND > > &flux) final
RoeSolver()
programmatic ctor
static RiemannSolverSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession)
static std::string solverName
std::map< std::string, RSParamFuncType > m_params
Map of parameter function types.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
RiemannSolverFactory & GetRiemannSolverFactory()
void RoeKernel(T &rhoL, T &rhouL, T &rhovL, T &rhowL, T &EL, T &rhoR, T &rhouR, T &rhovR, T &rhowR, T &ER, T &rhof, T &rhouf, T &rhovf, T &rhowf, T &Ef, NekDouble gamma)
tinysimd::simd< NekDouble > vec_t
static constexpr struct tinysimd::is_not_aligned_t is_not_aligned
typename abi< ScalarType, width >::type simd