38#include <boost/algorithm/string/predicate.hpp>
63 std::size_t nDim = pFields[0]->GetCoordim(0);
64 std::size_t nTracePts = pFields[0]->GetTrace()->GetTotPoints();
67 for (std::size_t i = 0; i < nDim; ++i)
75 const std::size_t nConvectiveFields,
82 std::size_t nCoeffs = fields[0]->GetNcoeffs();
85 for (std::size_t i = 0; i < nConvectiveFields; ++i)
93 for (std::size_t i = 0; i < nConvectiveFields; ++i)
95 fields[i]->BwdTrans(tmp[i], outarray[i]);
100 const std::size_t nConvectiveFields,
107 std::size_t nDim = fields[0]->GetCoordim(0);
108 std::size_t nPts = fields[0]->GetTotPoints();
109 std::size_t nCoeffs = fields[0]->GetNcoeffs();
110 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
115 for (std::size_t j = 0; j < nDim; ++j)
118 for (std::size_t i = 0; i < nConvectiveFields; ++i)
125 for (std::size_t i = 0; i < nConvectiveFields; ++i)
134 for (std::size_t j = 0; j < nDim; ++j)
137 for (std::size_t i = 0; i < nConvectiveFields; ++i)
148 for (std::size_t i = 0; i < nConvectiveFields; ++i)
150 for (std::size_t j = 0; j < nDim; ++j)
152 qdbase[j] = viscTensor[j][i];
154 fields[i]->IProductWRTDerivBase(qdbase, tmp);
157 fields[i]->AddTraceIntegral(traceflux[i], tmp);
158 fields[i]->SetPhysState(
false);
159 fields[i]->MultiplyByElmtInvMass(tmp, outarray[i]);
170 std::size_t nConvectiveFields = fields.size();
171 std::size_t nDim = fields[0]->GetCoordim(0);
172 std::size_t nCoeffs = fields[0]->GetNcoeffs();
173 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
177 for (std::size_t j = 0; j < nDim; ++j)
180 for (std::size_t i = 0; i < nConvectiveFields; ++i)
188 for (std::size_t j = 0; j < nDim; ++j)
190 for (std::size_t i = 0; i < nConvectiveFields; ++i)
192 fields[i]->IProductWRTDerivBase(j, inarray[i], tmp);
194 fields[i]->AddTraceIntegral(flux[j][i], tmp);
195 fields[i]->SetPhysState(
false);
196 fields[i]->MultiplyByElmtInvMass(tmp, tmp);
197 fields[i]->BwdTrans(tmp, qfield[j][i]);
231 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
232 std::size_t nvariables = fields.size();
233 std::size_t nDim = fields[0]->GetCoordim(0);
242 for (std::size_t i = 0; i < nvariables; ++i)
248 fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
260 if (fields[0]->GetBndCondExpansions().size())
265 for (std::size_t j = 0; j < nDim; ++j)
275 const std::size_t var,
282 std::size_t nBndRegions = fields[var]->GetBndCondExpansions().size();
284 for (std::size_t i = 0; i < nBndRegions; ++i)
286 if (fields[var]->GetBndConditions()[i]->GetBoundaryConditionType() ==
293 std::size_t nBndEdges =
294 fields[var]->GetBndCondExpansions()[i]->GetExpSize();
297 for (std::size_t e = 0; e < nBndEdges; ++e)
299 std::size_t nBndEdgePts = fields[var]
300 ->GetBndCondExpansions()[i]
305 fields[var]->GetBndCondExpansions()[i]->GetPhys_Offset(e);
307 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
308 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
312 fields[var]->GetBndConditions()[i]->GetUserDefined(),
315 fields[var]->GetBndConditions()[i]->GetUserDefined(),
318 fields[var]->GetBndConditions()[i]->GetUserDefined(),
321 fields[var]->GetBndConditions()[i]->GetUserDefined(),
324 Vmath::Vcopy(nBndEdgePts, &Fwd[id2], 1, &penaltyflux[id2], 1);
328 ->GetBndConditions()[i]
329 ->GetBoundaryConditionType() ==
334 &(fields[var]->GetBndCondExpansions()[i]->GetPhys())[id1],
335 1, &penaltyflux[id2], 1);
338 else if ((fields[var]->GetBndConditions()[i])
339 ->GetBoundaryConditionType() ==
342 Vmath::Vcopy(nBndEdgePts, &Fwd[id2], 1, &penaltyflux[id2], 1);
358 std::size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
359 std::size_t nvariables = fields.size();
360 std::size_t nDim = qfield.size();
371 for (std::size_t i = 0; i < nvariables; ++i)
374 fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
379 for (std::size_t j = 0; j < nDim; ++j)
382 fields[i]->GetFwdBwdTracePhys(qfield[j][i], qFwd, qBwd);
391 Vmath::Vadd(nTracePts, uterm, 1, qfluxtemp, 1, qfluxtemp, 1);
394 if (fields[0]->GetBndCondExpansions().size())
403 Vmath::Vadd(nTracePts, qfluxtemp, 1, qflux[i], 1, qflux[i], 1);
415 const std::size_t var,
const std::size_t dir,
421 std::size_t nBndRegions = fields[var]->GetBndCondExpansions().size();
424 for (std::size_t i = 0; i < nBndRegions; ++i)
426 if (fields[var]->GetBndConditions()[i]->GetBoundaryConditionType() ==
431 std::size_t nBndEdges =
432 fields[var]->GetBndCondExpansions()[i]->GetExpSize();
435 for (std::size_t e = 0; e < nBndEdges; ++e)
437 std::size_t nBndEdgePts = fields[var]
438 ->GetBndCondExpansions()[i]
443 fields[var]->GetBndCondExpansions()[i]->GetPhys_Offset(e);
445 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
446 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
450 fields[var]->GetBndConditions()[i]->GetUserDefined(),
453 fields[var]->GetBndConditions()[i]->GetUserDefined(),
456 fields[var]->GetBndConditions()[i]->GetUserDefined(),
459 fields[var]->GetBndConditions()[i]->GetUserDefined(),
467 ->GetBndConditions()[i]
468 ->GetBoundaryConditionType() ==
472 &qFwd[id2], 1, &penaltyflux[id2], 1);
475 else if ((fields[var]->GetBndConditions()[i])
476 ->GetBoundaryConditionType() ==
481 &(fields[var]->GetBndCondExpansions()[i]->GetPhys())[id1],
482 1, &penaltyflux[id2], 1);
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
SOLVER_UTILS_EXPORT void DiffuseCalcDerivative(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray)
SOLVER_UTILS_EXPORT void DiffuseVolumeFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, int > &nonZeroIndex=NullInt1DArray)
Diffusion Volume FLux.
SOLVER_UTILS_EXPORT void DiffuseTraceFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, Array< OneD, NekDouble > > &TraceFlux, const Array< OneD, Array< OneD, NekDouble > > &pFwd=NullNekDoubleArrayOfArray, const Array< OneD, Array< OneD, NekDouble > > &pBwd=NullNekDoubleArrayOfArray, Array< OneD, int > &nonZeroIndex=NullInt1DArray)
Diffusion term Trace Flux.
DiffusionFluxVecCB m_fluxVector
void v_Diffuse(const std::size_t nConvective, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
void ApplyScalarBCs(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const std::size_t var, const Array< OneD, const NekDouble > &ufield, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &penaltyflux)
void v_DiffuseVolumeFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, int > &nonZeroIndex) override
Diffusion Volume Flux.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
static DiffusionSharedPtr create(std::string diffType)
void ApplyVectorBCs(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const std::size_t var, const std::size_t dir, const Array< OneD, const NekDouble > &qfield, const Array< OneD, const NekDouble > &qFwd, const Array< OneD, const NekDouble > &qBwd, Array< OneD, NekDouble > &penaltyflux)
void v_DiffuseCalcDerivative(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
Diffusion Flux, calculate the physical derivatives.
void NumFluxforScalar(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &ufield, TensorOfArray3D< NekDouble > &uflux, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
void v_DiffuseCoeffs(const std::size_t nConvective, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
LibUtilities::SessionReaderSharedPtr m_session
std::string m_shockCaptureType
void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields) override
void NumFluxforVector(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &ufield, TensorOfArray3D< NekDouble > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
Build the numerical flux for the 2nd order derivatives todo: add variable coeff and h dependence to p...
NekDouble m_C11
Coefficient of penalty term.
void v_DiffuseTraceFlux(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfields, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, Array< OneD, NekDouble > > &TraceFlux, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, int > &nonZeroIndex) override
Diffusion term Trace Flux.
std::shared_ptr< SessionReader > SessionReaderSharedPtr
DiffusionFactory & GetDiffusionFactory()
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayOfArray
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
void Neg(int n, T *x, const int incx)
Negate x = -x.
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
void Zero(int n, T *x, const int incx)
Zero vector.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.