68 size_t nDim = pFields[0]->GetCoordim(0);
69 size_t nVariable = pFields.size();
70 size_t nTracePts = pFields[0]->GetTrace()->GetTotPoints();
73 for (
int i = 0; i < nDim; ++i)
79 for (
int i = 0; i < nVariable; ++i)
94 pFields[0]->GetTrace()->GetElmtNormalLength(lengthFwd, lengthBwd);
98 pFields[0]->GetTraceMap()->GetAssemblyCommDG()->PerformExchange(lengthFwd,
102 pFields[0]->PeriodicBwdCopy(lengthFwd, lengthBwd);
110 for (
int i = 0; i < nTracePts; ++i)
115 lengthBwd[i] = lengthFwd[i];
122 Vmath::Vadd(nTracePts, lengthBwd, 1, lengthFwd, 1, lengthsum, 1);
123 Vmath::Vmul(nTracePts, lengthBwd, 1, lengthFwd, 1, lengthmul, 1);
124 Vmath::Vdiv(nTracePts, lengthsum, 1, lengthmul, 1, lengthFwd, 1);
135 for (
int i = 1; i < nVariable; ++i)
137 pFields[i]->GetBwdWeight(tmpBwdWeight, tmpBwdWeightJump);
140 Vmath::Vabs(nTracePts, tmpBwdWeight, 1, tmpBwdWeight, 1);
142 for (
int j = 0; j < nTracePts; ++j)
144 norm += tmpBwdWeight[j];
147 "different BWD for different variable not coded yet");
151 size_t nCoeffs = pFields[0]->GetNcoeffs();
153 for (
int i = 0; i < nVariable; ++i)
161 for (
int nd = 0; nd < nDim; ++nd)
165 for (
int i = 0; i < nVariable; ++i)
174 const size_t nConvectiveFields,
186 for (
int i = 0; i < nConvectiveFields; ++i)
188 fields[i]->BwdTrans(
m_wspDiff[i], outarray[i]);
195 const size_t nConvectiveFields,
205 size_t nDim = fields[0]->GetCoordim(0);
206 size_t nPts = fields[0]->GetTotPoints();
207 size_t nCoeffs = fields[0]->GetNcoeffs();
208 size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
216 for (
int j = 0; j < nDim; ++j)
220 for (
int i = 0; i < nConvectiveFields; ++i)
224 for (
int i = 0; i < nConvectiveFields; ++i)
236 for (
int i = 0; i < nConvectiveFields; ++i)
241 for (
int i = 0; i < nConvectiveFields; ++i)
243 fields[i]->GetFwdBwdTracePhys(inarray[i], vFwd[i], vBwd[i]);
248 for (
int i = 0; i < nConvectiveFields; ++i)
265 for (
int i = 0; i < nonZeroIndex.size(); ++i)
267 int j = nonZeroIndex[i];
268 for (
int k = 0; k < nDim; ++k)
270 tmpFluxIprdct[k] = elmtFlux[k][j];
272 fields[j]->IProductWRTDerivBase(tmpFluxIprdct, outarray[j]);
280 for (
int j = 0; j < nDim; ++j)
287 for (
int j = 0; j < nConvectiveFields; ++j)
297 for (
int i = 0; i < nonZeroIndex.size(); ++i)
299 int j = nonZeroIndex[i];
301 fields[j]->AddTraceIntegral(Traceflux[j], outarray[j]);
302 fields[j]->SetPhysState(
false);
306 elmtFlux, outarray, vFwd, vBwd);
309 for (
int i = 0; i < nonZeroIndex.size(); ++i)
311 int j = nonZeroIndex[i];
313 fields[j]->MultiplyByElmtInvMass(outarray[j], outarray[j]);
327 size_t nDim = fields[0]->GetCoordim(0);
328 for (
int nd = 0; nd < 3; ++nd)
333 size_t nConvectiveFields = inarray.size();
334 for (
int i = 0; i < nConvectiveFields; ++i)
336 for (
int nd = 0; nd < nDim; ++nd)
338 qtmp[nd] = qfield[nd][i];
340 fields[i]->PhysDeriv(inarray[i], qtmp[0], qtmp[1], qtmp[2]);
350 size_t nDim = fields[0]->GetCoordim(0);
373 traceflux3D[0] = TraceFlux;
375 size_t nConvectiveFields = fields.size();
387 const std::size_t nConvectiveFields,
397 size_t nDim = fields[0]->GetCoordim(0);
398 size_t nPts = fields[0]->GetTotPoints();
399 size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
401 for (
int nd = 0; nd < nDim; ++nd)
405 for (
int j = 0; j < nConvectiveFields; ++j)
412 VolumeFlux, traceSymflux, pFwd, pBwd,
416 fields, nonZeroIndex, traceSymflux,
422 const std::size_t nConvectiveFields,
432 size_t nDim = fields[0]->GetCoordim(0);
433 size_t nPts = fields[0]->GetTotPoints();
434 size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
436 for (
int nd = 0; nd < nDim; ++nd)
440 for (
int j = 0; j < nConvectiveFields; ++j)
447 VolumeFlux, traceSymflux, pFwd, pBwd,
451 fields, nonZeroIndex, traceSymflux, outarray);
456 const std::size_t nConvectiveFields,
466 size_t nDim = fields[0]->GetCoordim(0);
469 nonZeroIndex, SymmFlux);
473 const std::size_t nConvectiveFields,
const size_t nDim,
479 size_t nTracePts = solution_jump[nConvectiveFields - 1].size();
484 for (
int nd = 0; nd < nDim; ++nd)
486 for (
int j = 0; j < nonZeroIndexsymm.size(); ++j)
488 int i = nonZeroIndexsymm[j];
490 traceSymflux[nd][i], 1, traceSymflux[nd][i], 1);
496 const std::size_t nConvectiveFields,
const size_t nDim,
const size_t nPts,
497 [[maybe_unused]]
const size_t nTracePts,
503 size_t nCoeffs = outarray[nConvectiveFields - 1].size();
506 for (
int i = 0; i < nDim; ++i)
511 for (
int j = 0; j < nonZeroIndex.size(); ++j)
513 nv = nonZeroIndex[j];
515 for (
int nd = 0; nd < nDim; ++nd)
519 tracelist->MultiplyByQuadratureMetric(tracflux[nd][nv],
522 fields[nv]->AddTraceQuadPhysToField(tracflux[nd][nv],
523 tracflux[nd][nv], tmpfield[nd]);
524 fields[nv]->DivideByQuadratureMetric(tmpfield[nd], tmpfield[nd]);
526 fields[nv]->IProductWRTDerivBase(tmpfield, tmpCoeff);
527 Vmath::Vadd(nCoeffs, tmpCoeff, 1, outarray[nv], 1, outarray[nv], 1);
532 const std::size_t nConvectiveFields,
const size_t nDim,
const size_t nPts,
533 [[maybe_unused]]
const size_t nTracePts,
539 size_t nCoeffs = outarray[nConvectiveFields - 1].size();
543 for (
int i = 0; i < nDim; ++i)
547 for (
int j = 0; j < nonZeroIndex.size(); ++j)
549 int nv = nonZeroIndex[j];
550 for (
int nd = 0; nd < nDim; ++nd)
554 fields[nv]->AddTraceQuadPhysToField(tracflux[nd][nv],
555 tracflux[nd][nv], tmpfield[nd]);
556 fields[nv]->DivideByQuadratureMetric(tmpfield[nd], tmpfield[nd]);
558 fields[nv]->IProductWRTDerivBase(tmpfield, tmpCoeff);
559 fields[nv]->BwdTrans(tmpCoeff, tmpPhysi);
560 Vmath::Vadd(nPts, tmpPhysi, 1, outarray[nv], 1, outarray[nv], 1);
569 std::shared_ptr<LocalRegions::ExpansionVector> traceExp =
571 size_t ntotTrac = (*traceExp).size();
572 int nTracPnt, noffset;
575 fields[0]->GetLocTraceToTraceMap();
578 locTraceToTraceMap->GetLeftRightAdjacentExpId();
580 locTraceToTraceMap->GetLeftRightAdjacentExpFlag();
582 std::shared_ptr<LocalRegions::ExpansionVector> fieldExp =
591 for (
int ntrace = 0; ntrace < ntotTrac; ++ntrace)
593 noffset = tracelist->GetPhys_Offset(ntrace);
594 nTracPnt = tracelist->GetTotPoints(ntrace);
596 factorFwdBwd[0] = 0.0;
597 factorFwdBwd[1] = 0.0;
599 for (
int nlr = 0; nlr < 2; ++nlr)
601 if (LRAdjflag[nlr][ntrace])
604 for (
int nd = 0; nd < spaceDim; nd++)
607 ->GetExp(LRAdjExpid[nlr][ntrace])
608 ->GetBasisNumModes(nd);
609 numModes = std::max(ntmp, numModes);
611 factorFwdBwd[nlr] = (numModes) * (numModes);
615 for (
int np = 0; np < nTracPnt; ++np)
617 factor[noffset + np] = std::max(factorFwdBwd[0], factorFwdBwd[1]);
623 const std::size_t nConvectiveFields,
const size_t nPts,
629 std::vector<NekDouble> vFwdTmp(nConvectiveFields),
630 vBwdTmp(nConvectiveFields), averTmp(nConvectiveFields);
631 for (
size_t p = 0;
p < nPts; ++
p)
634 for (
size_t f = 0; f < nConvectiveFields; ++f)
636 vFwdTmp[f] = vFwd[f][
p];
637 vBwdTmp[f] = vBwd[f][
p];
644 for (
size_t f = 0; f < nConvectiveFields; ++f)
646 aver[f][
p] = averTmp[f];
657 for (
size_t f = 0; f < nConvectiveFields; ++f)
659 for (
size_t p = 0;
p < nPts; ++
p)
666 jump[f][
p] = tmpF * Fweight + tmpB * Bweight;
682 size_t nDim = fields[0]->GetCoordim(0);
683 size_t nPts = fields[0]->GetTotPoints();
684 size_t nTracePts = fields[0]->GetTrace()->GetTotPoints();
685 size_t nConvectiveFields = fields.size();
688 fields[0]->GetTraceMap();
693 for (
int nd = 0; nd < nDim; ++nd)
695 for (
int i = 0; i < nConvectiveFields; ++i)
710 if (fabs(PenaltyFactor2) > 1.0E-12)
720 for (
int nd = 0; nd < nDim; ++nd)
722 for (
int i = 0; i < nConvectiveFields; ++i)
728 fields[i]->GetFwdBwdTracePhys(qfield[nd][i], Fwd, Bwd,
true,
true,
733 for (
size_t p = 0;
p < nTracePts; ++
p)
740 TraceMap->GetAssemblyCommDG()->PerformExchange(
751 ConsVarAveJump(nConvectiveFields, nTracePts, vFwd, vBwd, solution_Aver,
758 for (
size_t p = 0;
p < nTracePts; ++
p)
763 for (
size_t f = 0; f < nConvectiveFields; ++f)
765 NekDouble jumpTmp = solution_jump[f][
p] * PenaltyFactor;
766 for (
size_t d = 0;
d < nDim; ++
d)
778 traceflux, nonZeroIndexflux,
785 const std::size_t nConvectiveFields,
const size_t nDim,
const size_t nPts,
786 const size_t nTracePts,
const NekDouble PenaltyFactor2,
794 std::vector<NekDouble> tmp(nTracePts);
797 for (
int nd1 = 0; nd1 < nDim; ++nd1)
803 for (
int nd = 0; nd < 3; ++nd)
807 for (
int nd2 = 0; nd2 < nDim; ++nd2)
809 qtmp[nd2] = elmt2ndDerv[nd2];
812 for (
int i = 0; i < nTracePts; ++i)
817 for (
int nd1 = 0; nd1 < nDim; ++nd1)
819 for (
int i = 0; i < nConvectiveFields; ++i)
821 fields[i]->PhysDeriv(qfield[nd1][i], qtmp[0], qtmp[1], qtmp[2]);
823 for (
int nd2 = nd1; nd2 < nDim; ++nd2)
826 fields[i]->GetFwdBwdTracePhys(elmt2ndDerv[nd2], Fwd, Bwd,
true,
828 for (
int p = 0;
p < nTracePts; ++
p)
834 numDerivFwd[nd1][i][
p]);
838 for (
int p = 0;
p < nTracePts; ++
p)
840 numDerivBwd[nd2][i][
p] +=
842 numDerivFwd[nd2][i][
p] =
844 numDerivFwd[nd2][i][
p]);
858 const int nConvectiveFields,
862 int nengy = nConvectiveFields - 1;
866 size_t nBndRegions = fields[nengy]->GetBndCondExpansions().size();
867 for (
size_t j = 0; j < nBndRegions; ++j)
869 if (fields[nengy]->GetBndConditions()[j]->GetBoundaryConditionType() ==
876 fields[nengy]->GetBndCondExpansions()[j]->GetExpSize();
877 for (
size_t e = 0; e < nBndEdges; ++e)
879 size_t nBndEdgePts = fields[nengy]
880 ->GetBndCondExpansions()[j]
884 std::size_t id2 = fields[0]->GetTrace()->GetPhys_Offset(
885 fields[0]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt++));
887 if (fields[0]->GetBndConditions()[j]->GetUserDefined() ==
#define ASSERTL0(condition, msg)
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
void AccumulateRegion(std::string, int iolevel=0)
Accumulate elapsed time for a region.
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)
DiffusionFluxCons m_FunctorDiffusionfluxConsTrace
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.
DiffusionFluxCons m_FunctorDiffusionfluxCons
SpecialBndTreat m_SpecialBndTreat
DiffusionSymmFluxCons m_FunctorSymmetricfluxCons
Array< OneD, NekDouble > m_traceNormDirctnElmtLengthRecip
void v_DiffuseCoeffs(const std::size_t nConvectiveFields, 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
static DiffusionSharedPtr create(std::string diffType)
Array< OneD, NekDouble > m_tracBwdWeightAver
Array< OneD, NekDouble > m_traceNormDirctnElmtLength
void ApplyFluxBndConds(const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, Array< OneD, Array< OneD, NekDouble > > &flux)
aplly Neuman boundary conditions on flux Currently only consider WallAdiabatic
void GetPenaltyFactor(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, Array< OneD, NekDouble > &factor)
Get IP penalty factor based on order.
void ConsVarAve(const size_t nConvectiveFields, const T &Bweight, const std::vector< T > &vFwd, const std::vector< T > &vBwd, std::vector< T > &aver)
Calculate the average of conservative variables on traces.
void AddSecondDerivToTrace(const std::size_t nConvectiveFields, const size_t nDim, const size_t nPts, const size_t nTracePts, const NekDouble PenaltyFactor2, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const TensorOfArray3D< NekDouble > &qfield, TensorOfArray3D< NekDouble > &numDerivFwd, TensorOfArray3D< NekDouble > &numDerivBwd)
Add second derivative term to trace jump (for DDG scheme)
void CalcTraceNumFlux(const NekDouble PenaltyFactor2, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, const TensorOfArray3D< NekDouble > &qfield, const Array< OneD, Array< OneD, NekDouble > > &vFwd, const Array< OneD, Array< OneD, NekDouble > > &vBwd, Array< OneD, int > &nonZeroIndexflux, TensorOfArray3D< NekDouble > &traceflux, Array< OneD, Array< OneD, NekDouble > > &solution_Aver, Array< OneD, Array< OneD, NekDouble > > &solution_jump)
Calculate numerical flux on traces.
NekDouble m_IPSymmFluxCoeff
NekDouble m_IP2ndDervCoeff
NekDouble m_IPPenaltyCoeff
void v_DiffuseCalcDerivative(const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd) override
Diffusion Flux, calculate the physical derivatives.
void AddDiffusionSymmFluxToPhys(const std::size_t nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
Array< OneD, Array< OneD, NekDouble > > m_traceJump
void AddSymmFluxIntegralToCoeff(const std::size_t nConvectiveFields, const size_t nDim, const size_t nPts, const size_t nTracePts, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, const int > &nonZeroIndex, TensorOfArray3D< NekDouble > &tracflux, Array< OneD, Array< OneD, NekDouble > > &outarray)
Add symmetric flux integration to field (in coefficient space)
TensorOfArray3D< NekDouble > m_wspNumDerivBwd
Workspace for CallTraceNumFlux.
TensorOfArray3D< NekDouble > m_wspNumDerivFwd
LibUtilities::SessionReaderSharedPtr m_session
Array< OneD, Array< OneD, NekDouble > > m_traceAver
Array< OneD, NekDouble > m_tracBwdWeightJump
void v_Diffuse(const std::size_t nConvectiveFields, 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 ConsVarAveJump(const std::size_t nConvectiveFields, const size_t npnts, const Array< OneD, const Array< OneD, NekDouble > > &vFwd, const Array< OneD, const Array< OneD, NekDouble > > &vBwd, Array< OneD, Array< OneD, NekDouble > > &aver, Array< OneD, Array< OneD, NekDouble > > &jump)
void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields) override
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.
void AddDiffusionSymmFluxToCoeff(const std::size_t nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, TensorOfArray3D< NekDouble > &VolumeFlux, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
void DiffuseTraceSymmFlux(const std::size_t nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &inarray, TensorOfArray3D< NekDouble > &qfield, TensorOfArray3D< NekDouble > &VolumeFlux, TensorOfArray3D< NekDouble > &SymmFlux, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd, Array< OneD, int > &nonZeroIndex)
Calculate symmetric flux on traces interface.
Array< OneD, Array< OneD, NekDouble > > m_wspDiff
Workspace for v_Diffusion.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
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.
void CalcTraceSymFlux(const std::size_t nConvectiveFields, const size_t nDim, const Array< OneD, Array< OneD, NekDouble > > &solution_Aver, Array< OneD, Array< OneD, NekDouble > > &solution_jump, Array< OneD, int > &nonZeroIndexsymm, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &traceSymflux)
Calculate symmetric flux on traces.
void AddSymmFluxIntegralToPhys(const std::size_t nConvectiveFields, const size_t nDim, const size_t nPts, const size_t nTracePts, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, const int > &nonZeroIndex, TensorOfArray3D< NekDouble > &tracflux, Array< OneD, Array< OneD, NekDouble > > &outarray)
Add symmetric flux integration to field (in physical space)
std::shared_ptr< SessionReader > SessionReaderSharedPtr
std::shared_ptr< AssemblyMapDG > AssemblyMapDGSharedPtr
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
std::shared_ptr< LocTraceToTraceMap > LocTraceToTraceMapSharedPtr
static const NekDouble kNekMachineEpsilon
DiffusionFactory & GetDiffusionFactory()
std::vector< double > d(NPUPPER *NPUPPER)
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayOfArray
static Array< OneD, NekDouble > NullNekDouble1DArray
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 Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
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 Vdiv(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 Zero(int n, T *x, const int incx)
Zero vector.
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.
scalarT< T > abs(scalarT< T > in)