57 : StdExpansion(Ba.GetNumModes(), 1, Ba),
58 StdExpansion1D(Ba.GetNumModes(), Ba), StdRegions::StdSegExp(Ba),
61 std::bind(&
SegExp::CreateMatrix, this,
std::placeholders::_1),
62 std::string(
"SegExpMatrix")),
63 m_staticCondMatrixManager(
std::bind(&
Expansion::CreateStaticCondMatrix,
64 this,
std::placeholders::_1),
65 std::string(
"SegExpStaticCondMatrix"))
74 : StdExpansion(S), StdExpansion1D(S), StdRegions::StdSegExp(S),
76 m_staticCondMatrixManager(S.m_staticCondMatrixManager)
100 int nquad0 =
m_base[0]->GetNumPoints();
116 ival = StdSegExp::v_Integral(tmp);
150 int nquad0 =
m_base[0]->GetNumPoints();
161 Vmath::Vmul(nquad0, &gmat[0][0], 1, &diff[0], 1, &out_d0[0], 1);
166 Vmath::Vmul(nquad0, &gmat[1][0], 1, &diff[0], 1, &out_d1[0], 1);
171 Vmath::Vmul(nquad0, &gmat[2][0], 1, &diff[0], 1, &out_d2[0], 1);
178 Vmath::Smul(nquad0, gmat[0][0], diff, 1, out_d0, 1);
183 Vmath::Smul(nquad0, gmat[1][0], diff, 1, out_d1, 1);
188 Vmath::Smul(nquad0, gmat[2][0], diff, 1, out_d2, 1);
204 int nquad0 =
m_base[0]->GetNumPoints();
205 int coordim =
m_geom->GetCoordim();
243 int nquad0 =
m_base[0]->GetNumPoints();
246 int coordim =
m_geom->GetCoordim();
258 cout <<
"der_n" << endl;
259 for (
int k = 0; k < coordim; ++k)
267 for (
int i = 0; i < nquad0; i++)
269 cout <<
"nx= " << normals[0][i] <<
" ny=" << normals[1][i]
273 "normal vectors do not exist: check if a"
274 "boundary region is defined as I ");
276 Vmath::Vmul(nquad0, normals[0], 1, inarray_d0, 1, out_dn_tmp, 1);
277 Vmath::Vadd(nquad0, out_dn_tmp, 1, out_dn, 1, out_dn, 1);
279 Vmath::Vmul(nquad0, normals[1], 1, inarray_d1, 1, out_dn_tmp, 1);
280 Vmath::Vadd(nquad0, out_dn_tmp, 1, out_dn, 1, out_dn, 1);
282 for (
int i = 0; i < nquad0; i++)
284 cout <<
"deps/dx =" << inarray_d0[i]
285 <<
" deps/dy=" << inarray_d1[i] << endl;
315 ASSERTL1(
false,
"input dir is out of range");
349 if (
m_base[0]->Collocation())
365 out = (*matsys) * in;
373 if (
m_base[0]->Collocation())
403 "Cannot use FwdTrans_BndConstrained with these points.");
408 ASSERTL0(
false,
"This type of FwdTrans is not defined"
409 "for this expansion type");
412 fill(outarray.data(), outarray.data() +
m_ncoeffs, 0.0);
440 Blas::Dgemv(
'N', nInteriorDofs, nInteriorDofs, matsys->Scale(),
441 &((matsys->GetOwnedMatrix())->GetPtr())[0],
442 nInteriorDofs, tmp1.data() + offset, 1, 0.0,
443 outarray.data() + offset, 1);
509 int nquad0 =
m_base[0]->GetNumPoints();
522 StdSegExp::v_IProductWRTBase(base, tmp, outarray, coll_check);
529 ASSERTL1(dir < 3,
"input dir is out of range");
531 "input dir is out of range");
533 int nquad =
m_base[0]->GetNumPoints();
541 Vmath::Vmul(nquad, gmat[dir], 1, inarray, 1, tmp1, 1);
545 Vmath::Smul(nquad, gmat[dir][0], inarray, 1, tmp1, 1);
555 int nq =
m_base[0]->GetNumPoints();
563 Vmath::Vmul(nq, &Fx[0], 1, &normals[0][0], 1, &Fn[0], 1);
564 Vmath::Vvtvp(nq, &Fy[0], 1, &normals[1][0], 1, &Fn[0], 1, &Fn[0], 1);
590 return StdExpansion1D::v_PhysEvaluate(Lcoord, physvals);
599 m_geom->GetLocCoords(coord, Lcoord);
601 return StdExpansion1D::v_PhysEvaluate(Lcoord, physvals);
607 std::array<NekDouble, 3> &firstOrderDerivs)
611 m_geom->GetLocCoords(coord, Lcoord);
612 return StdSegExp::v_PhysEvalFirstDeriv(Lcoord, inarray, firstOrderDerivs);
618 std::array<NekDouble, 3> &firstOrderDerivs,
619 std::array<NekDouble, 6> &secondOrderDerivs)
623 m_geom->GetLocCoords(coord, Lcoord);
624 return StdSegExp::v_PhysEvalFirstSecondDeriv(
625 Lcoord, inarray, firstOrderDerivs, secondOrderDerivs);
633 ASSERTL1(Lcoords[0] >= -1.0 && Lcoords[0] <= 1.0,
634 "Local coordinates are not in region [-1,1]");
637 for (i = 0; i <
m_geom->GetCoordim(); ++i)
639 coords[i] =
m_geom->GetCoord(i, Lcoords);
655 int nquad =
m_base[0]->GetNumPoints();
662 outarray = inarray[0];
665 outarray = inarray[nquad - 1];
680 Blas::Ddot(nquad, mat_gauss->GetOwnedMatrix()->GetPtr().data(), 1,
695 outarray[0] = result;
701 int nquad =
m_base[0]->GetNumPoints();
703 ASSERTL1(vertex == 0 || vertex == 1,
"Vertex value should be 0 or 1");
707 map[0] = vertex == 0 ? 0 : nquad - 1;
721 if (&inarray[0] != &outarray[0])
736 m_base[0]->GetBasisKey());
742 m_base[0]->GetPointsKey());
766 const NekDouble *data,
const std::vector<unsigned int> &nummodes,
767 const int mode_offset,
NekDouble *coeffs,
768 [[maybe_unused]] std::vector<LibUtilities::BasisType> &fromType)
774 int fillorder = min((
int)nummodes[mode_offset],
m_ncoeffs);
803 ASSERTL0(
false,
"basis is either not set up or not hierarchicial");
820 for (i = 0; i < vCoordDim; ++i)
826 size_t nbnd = vertex;
839 for (i = 0; i < vCoordDim; ++i)
845 for (i = 0; i < vCoordDim; ++i)
851 ASSERTL0(
false,
"point is out of range (point < 2)");
856 for (i = 0; i < vCoordDim; ++i)
858 vert += normal[i][0] * normal[i][0];
860 vert = 1.0 /
sqrt(vert);
864 for (i = 0; i < vCoordDim; ++i)
866 Vmath::Smul(nqe, vert, normal[i], 1, normal[i], 1);
880 int nquad =
m_base[0]->GetNumPoints();
890 switch (
m_geom->GetCoordim())
899 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
900 dPhysValuesdx.data(), 1);
904 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
912 PhysDeriv(physValues, dPhysValuesdx, dPhysValuesdy);
917 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
918 dPhysValuesdx.data(), 1);
919 Vmath::Vvtvp(nquad, &gmat[1][0], 1, dPhysValuesdy.data(), 1,
920 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
924 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
925 Blas::Daxpy(nquad, gmat[1][0], dPhysValuesdy.data(), 1,
926 dPhysValuesdx.data(), 1);
935 PhysDeriv(physValues, dPhysValuesdx, dPhysValuesdy, dPhysValuesdz);
940 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
941 dPhysValuesdx.data(), 1);
942 Vmath::Vvtvp(nquad, &gmat[1][0], 1, dPhysValuesdy.data(), 1,
943 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
944 Vmath::Vvtvp(nquad, &gmat[2][0], 1, dPhysValuesdz.data(), 1,
945 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
949 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
950 Blas::Daxpy(nquad, gmat[1][0], dPhysValuesdy.data(), 1,
951 dPhysValuesdx.data(), 1);
952 Blas::Daxpy(nquad, gmat[2][0], dPhysValuesdz.data(), 1,
953 dPhysValuesdx.data(), 1);
958 ASSERTL0(
false,
"Wrong number of dimensions");
970 StdExpansion::LaplacianMatrixOp_MatFree(k1, k2, inarray, outarray, mkey);
977 int nquad =
m_base[0]->GetNumPoints();
992 switch (
m_geom->GetCoordim())
1001 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
1002 dPhysValuesdx.data(), 1);
1006 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
1014 PhysDeriv(physValues, dPhysValuesdx, dPhysValuesdy);
1019 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
1020 dPhysValuesdx.data(), 1);
1021 Vmath::Vvtvp(nquad, &gmat[1][0], 1, dPhysValuesdy.data(), 1,
1022 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
1026 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
1027 Blas::Daxpy(nquad, gmat[1][0], dPhysValuesdy.data(), 1,
1028 dPhysValuesdx.data(), 1);
1037 PhysDeriv(physValues, dPhysValuesdx, dPhysValuesdy, dPhysValuesdz);
1042 Vmath::Vmul(nquad, &gmat[0][0], 1, dPhysValuesdx.data(), 1,
1043 dPhysValuesdx.data(), 1);
1044 Vmath::Vvtvp(nquad, &gmat[1][0], 1, dPhysValuesdy.data(), 1,
1045 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
1046 Vmath::Vvtvp(nquad, &gmat[2][0], 1, dPhysValuesdz.data(), 1,
1047 dPhysValuesdx.data(), 1, dPhysValuesdx.data(), 1);
1051 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.data(), 1);
1052 Blas::Daxpy(nquad, gmat[1][0], dPhysValuesdy.data(), 1,
1053 dPhysValuesdx.data(), 1);
1054 Blas::Daxpy(nquad, gmat[2][0], dPhysValuesdz.data(), 1,
1055 dPhysValuesdx.data(), 1);
1060 ASSERTL0(
false,
"Wrong number of dimensions");
1098 return tmp->GetStdMatrix(mkey);
1108 "Geometric information is not set up");
1118 goto UseLocRegionsMatrix;
1123 goto UseStdRegionsMatrix;
1144 goto UseStdRegionsMatrix;
1156 goto UseLocRegionsMatrix;
1168 "Cannot call eWeakDeriv2 in a "
1169 "coordinate system which is not at "
1170 "least two-dimensional");
1175 "Cannot call eWeakDeriv2 in a "
1176 "coordinate system which is not "
1177 "three-dimensional");
1201 goto UseLocRegionsMatrix;
1205 int coordim =
m_geom->GetCoordim();
1207 for (
int i = 0; i < coordim; ++i)
1213 goto UseStdRegionsMatrix;
1226 int rows = LapMat.GetRows();
1227 int cols = LapMat.GetColumns();
1233 (*helm) = LapMat + factor * MassMat;
1274 coords[0] = (vertex == 0) ? -1.0 : 1.0;
1276 m_Ix =
m_base[0]->GetI(coords);
1282 UseLocRegionsMatrix:
1288 UseStdRegionsMatrix:
1317 returnval = StdSegExp::v_GenMatrix(mkey);
1339 ASSERTL1(&inarray[0] != &outarray[0],
1340 "inarray and outarray can not be the same");
1345 outarray[0] = inarray[1];
1346 outarray[1] = inarray[0];
1350 outarray[m] = sgn * inarray[m];
1358 outarray[
m_ncoeffs - 1 - m] = inarray[m];
1362 ASSERTL0(
false,
"This basis is not allowed in this method");
1383 out = (*matsys) * in;
#define ASSERTL0(condition, msg)
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed to...
Describes the specification for a Basis.
Defines a specification for a set of points.
DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey) override
std::map< int, NormalVector > m_traceNormals
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
SpatialDomains::GeometrySharedPtr GetGeom() const
SpatialDomains::GeometrySharedPtr m_geom
ExpansionSharedPtr GetLeftAdjacentElementExp() const
void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
int GetLeftAdjacentElementTrace() const
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
void v_LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Forward transform from physical quadrature space stored in inarray and evaluate the expansion coeffic...
void v_DropLocMatrix(const MatrixKey &mkey) override
StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const override
void v_GetTracePhysVals(const int edge, const StdRegions::StdExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient) override
void v_PhysDeriv_s(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_ds) override
Evaluate the derivative along a line: . The derivative is calculated performing the product .
void v_GetVertexPhysVals(const int vertex, const Array< OneD, const NekDouble > &inarray, NekDouble &outarray) override
NekDouble v_PhysEvalFirstSecondDeriv(const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs, std::array< NekDouble, 6 > &secondOrderDerivs) override
NekDouble v_StdPhysEvaluate(const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals) override
void ReverseCoeffsAndSign(const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Reverse the coefficients in a boundary interior expansion this routine is of use when we need the seg...
SegExp(const LibUtilities::BasisKey &Ba, const SpatialDomains::Geometry1DSharedPtr &geom)
Constructor using BasisKey class for quadrature points and order definition.
void v_NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray) override
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const override
void v_FwdTransBndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
Inner product of inarray over region with respect to the expansion basis (this)->_Base[0] and return ...
DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey) override
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray) override
Integrate the physical point list inarray over region and return the value.
void v_HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey) override
void v_PhysDeriv_n(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dn) override
Evaluate the derivative normal to a line: . The derivative is calculated performing the product .
void v_ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType) override
Unpack data from input file assuming it comes from.
void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray) override
Evaluate the derivative at the physical quadrature points given by inarray and return in outarray.
void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3) override
DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix(const MatrixKey &mkey) override
void MultiplyByElmtInvMass(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
int v_NumBndryCoeffs() const override
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coord, const Array< OneD, const NekDouble > &physvals) override
void v_GetTracePhysMap(const int vertex, Array< OneD, int > &map) override
int v_NumDGBndryCoeffs() const override
DNekMatSharedPtr v_CreateStdMatrix(const StdRegions::StdMatrixKey &mkey) override
void v_SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
void v_DropLocStaticCondMatrix(const MatrixKey &mkey) override
const Array< OneD, const NekDouble > & v_GetPhysNormals() override
DNekScalMatSharedPtr v_GetLocMatrix(const MatrixKey &mkey) override
void v_GetCoord(const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords) override
NekDouble v_PhysEvalFirstDeriv(const Array< OneD, NekDouble > &coord, const Array< OneD, const NekDouble > &inarray, std::array< NekDouble, 3 > &firstOrderDerivs) override
void v_ComputeTraceNormal(const int vertex) override
void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray) override
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
void PhysTensorDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Evaluate the derivative at the physical quadrature points given by inarray and return in outarray.
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
const LibUtilities::PointsKeyVector GetPointsKeys() const
void MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
int GetVertexMap(const int localVertexId, bool useCoeffPacking=false)
void NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Backward transformation from coefficient space to physical space.
DNekMatSharedPtr GenMatrix(const StdMatrixKey &mkey)
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1=NullNekDouble1DArray, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Array< OneD, LibUtilities::BasisSharedPtr > m_base
LibUtilities::ShapeType GetShapeType() const
const VarCoeffMap & GetVarCoeffs() const
MatrixType GetMatrixType() const
const ConstFactorMap & GetConstFactors() const
NekDouble GetConstFactor(const ConstFactorType &factor) const
static void Dgemv(const char &trans, const int &m, const int &n, const double &alpha, const double *a, const int &lda, const double *x, const int &incx, const double &beta, double *y, const int &incy)
BLAS level 2: Matrix vector multiply y = alpha A x plus beta y where A[m x n].
static void Dscal(const int &n, const double &alpha, double *x, const int &incx)
BLAS level 1: x = alpha x.
static double Ddot(const int &n, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: output = .
static void Daxpy(const int &n, const double &alpha, const double *x, const int &incx, const double *y, const int &incy)
BLAS level 1: y = alpha x plus y.
void Interp1D(const BasisKey &fbasis0, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, Array< OneD, NekDouble > &to)
this function interpolates a 1D function evaluated at the quadrature points of the basis fbasis0 to ...
std::vector< PointsKey > PointsKeyVector
@ eGaussLobattoLegendre
1D Gauss-Lobatto-Legendre quadrature points
@ ePolyEvenlySpaced
1D Evenly-spaced points using Lagrange polynomial
@ eGaussGaussLegendre
1D Gauss-Gauss-Legendre quadrature points
@ eModified_B
Principle Modified Functions .
@ eGauss_Lagrange
Lagrange Polynomials using the Gauss points.
@ eGLL_Lagrange
Lagrange for SEM basis .
@ eModified_A
Principle Modified Functions .
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
GeomType
Indicates the type of element geometry.
@ eRegular
Geometry is straight-sided with constant geometric factors.
@ eNoGeomType
No type defined.
@ eMovingRegular
Currently unused.
@ eDeformed
Geometry is curved or has non-constant factors.
std::shared_ptr< Geometry1D > Geometry1DSharedPtr
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::map< ConstFactorType, NekDouble > ConstFactorMap
std::shared_ptr< StdSegExp > StdSegExpSharedPtr
StdRegions::ConstFactorMap factors
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayOfArray
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
static Array< OneD, NekDouble > NullNekDouble1DArray
std::shared_ptr< DNekMat > DNekMatSharedPtr
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 Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
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 Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
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.
scalarT< T > sqrt(scalarT< T > in)