44 namespace LocalRegions
60 StdExpansion(Ba.GetNumModes(), 1, Ba),
61 StdExpansion1D(Ba.GetNumModes(), Ba),
62 StdRegions::StdSegExp(Ba),
66 boost::bind(&
SegExp::CreateMatrix, this, _1),
67 std::string(
"SegExpMatrix")),
68 m_staticCondMatrixManager(
69 boost::bind(&
SegExp::CreateStaticCondMatrix, this, _1),
70 std::string(
"SegExpStaticCondMatrix"))
82 StdRegions::StdSegExp(S),
85 m_matrixManager(S.m_matrixManager),
86 m_staticCondMatrixManager(S.m_staticCondMatrixManager)
120 int nquad0 =
m_base[0]->GetNumPoints();
136 ival = StdSegExp::v_Integral(tmp);
171 int nquad0 =
m_base[0]->GetNumPoints();
180 if(out_d0.num_elements())
186 if(out_d1.num_elements())
192 if(out_d2.num_elements())
200 if(out_d0.num_elements())
206 if(out_d1.num_elements())
212 if(out_d2.num_elements())
232 int nquad0 =
m_base[0]->GetNumPoints();
233 int coordim =
m_geom->GetCoordim();
272 int nquad0 =
m_base[0]->GetNumPoints();
275 int coordim =
m_geom->GetCoordim();
288 for(
int k=0; k<coordim; ++k)
295 for(
int i=0; i<nquad0; i++)
297 cout<<
"nx= "<<normals[0][i]<<
" ny="<<normals[1][i]<<endl;
300 "normal vectors do not exist: check if a"
301 "boundary region is defined as I ");
303 Vmath::Vmul(nquad0,normals[0],1,inarray_d0,1,out_dn_tmp,1);
304 Vmath::Vadd(nquad0,out_dn_tmp,1,out_dn,1,out_dn,1);
306 Vmath::Vmul(nquad0,normals[1],1,inarray_d1,1,out_dn_tmp,1);
307 Vmath::Vadd(nquad0,out_dn_tmp,1,out_dn,1,out_dn,1);
309 for(
int i=0; i<nquad0; i++)
311 cout<<
"deps/dx ="<<inarray_d0[i]<<
" deps/dy="<<inarray_d1[i]<<endl;
344 ASSERTL1(
false,
"input dir is out of range");
380 if (
m_base[0]->Collocation())
404 if(
m_base[0]->Collocation())
434 ASSERTL0(
false,
"This type of FwdTrans is not defined"
435 "for this expansion type");
438 fill(outarray.get(), outarray.get()+
m_ncoeffs, 0.0 );
445 inarray[
m_base[0]->GetNumPoints()-1];
468 Blas::Dgemv(
'N',nInteriorDofs,nInteriorDofs,
470 &((matsys->GetOwnedMatrix())->GetPtr())[0],
471 nInteriorDofs,tmp1.get()+offset,1,0.0,
472 outarray.get()+offset,1);
544 int nquad0 =
m_base[0]->GetNumPoints();
558 StdSegExp::v_IProductWRTBase(base,tmp,outarray,coll_check);
567 int nquad =
m_base[0]->GetNumPoints();
583 Vmath::Smul(nquad, gmat[0][0], inarray, 1, tmp1, 1);
595 Vmath::Smul(nquad, gmat[1][0], inarray, 1, tmp1, 1);
608 Vmath::Smul(nquad, gmat[2][0], inarray, 1, tmp1, 1);
614 ASSERTL1(
false,
"input dir is out of range");
626 int nq =
m_base[0]->GetNumPoints();
635 Vmath::Vmul (nq, &Fx[0], 1, &normals[0][0], 1, &Fn[0], 1);
636 Vmath::Vvtvp(nq, &Fy[0], 1, &normals[1][0], 1, &Fn[0], 1, &Fn[0], 1);
663 return StdSegExp::v_PhysEvaluate(Lcoord,physvals);
673 m_geom->GetLocCoords(coord,Lcoord);
675 return StdSegExp::v_PhysEvaluate(Lcoord, physvals);
685 ASSERTL1(Lcoords[0] >= -1.0&& Lcoords[0] <= 1.0,
686 "Local coordinates are not in region [-1,1]");
689 for(i = 0; i <
m_geom->GetCoordim(); ++i)
691 coords[i] =
m_geom->GetCoord(i,Lcoords);
709 int nquad =
m_base[0]->GetNumPoints();
716 outarray = inarray[0];
719 outarray = inarray[nquad - 1];
734 outarray =
Blas::Ddot(nquad, mat_gauss->GetOwnedMatrix()
735 ->GetPtr().get(), 1, &inarray[0], 1);
749 outarray[0] = result;
771 if (&inarray[0] != &outarray[0])
785 return m_geom->GetPorient(point);
797 return m_geom->GetCoordim();
843 const std::vector<unsigned int > &nummodes,
844 const int mode_offset,
851 int fillorder = min((
int) nummodes[mode_offset],
m_ncoeffs);
862 nummodes[mode_offset],
865 p0,data,
m_base[0]->GetPointsKey(), coeffs);
873 nummodes[mode_offset],
876 p0,data,
m_base[0]->GetPointsKey(), coeffs);
881 "basis is either not set up or not hierarchicial");
900 for (i = 0; i < vCoordDim; ++i)
914 for(i = 0; i < vCoordDim; ++i)
920 for(i = 0; i < vCoordDim; ++i)
927 "point is out of range (point < 2)");
932 for (i =0 ; i < vCoordDim; ++i)
934 vert += normal[i][0]*normal[i][0];
936 vert = 1.0/sqrt(vert);
937 for (i = 0; i < vCoordDim; ++i)
939 Vmath::Smul(nqe, vert, normal[i], 1, normal[i], 1);
954 int nquad =
m_base[0]->GetNumPoints();
964 switch (
m_geom->GetCoordim())
974 &gmat[0][0],1,dPhysValuesdx.get(),1,
975 dPhysValuesdx.get(),1);
980 gmat[0][0], dPhysValuesdx.get(), 1);
988 PhysDeriv(physValues,dPhysValuesdx,dPhysValuesdy);
994 &gmat[0][0],1,dPhysValuesdx.get(),1,
995 dPhysValuesdx.get(),1);
997 &gmat[1][0],1,dPhysValuesdy.get(),1,
998 dPhysValuesdx.get(),1,
999 dPhysValuesdx.get(),1);
1004 gmat[0][0], dPhysValuesdx.get(), 1);
1006 gmat[1][0], dPhysValuesdy.get(), 1,
1007 dPhysValuesdx.get(), 1);
1017 dPhysValuesdy,dPhysValuesdz);
1023 &gmat[0][0], 1, dPhysValuesdx.get(), 1,
1024 dPhysValuesdx.get(),1);
1026 &gmat[1][0], 1, dPhysValuesdy.get(), 1,
1027 dPhysValuesdx.get(),1,
1028 dPhysValuesdx.get(),1);
1030 &gmat[2][0], 1, dPhysValuesdz.get(), 1,
1031 dPhysValuesdx.get(),1,
1032 dPhysValuesdx.get(),1);
1036 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.get(), 1);
1038 gmat[1][0], dPhysValuesdy.get(), 1,
1039 dPhysValuesdx.get(), 1);
1041 gmat[2][0], dPhysValuesdz.get(), 1,
1042 dPhysValuesdx.get(), 1);
1047 ASSERTL0(
false,
"Wrong number of dimensions");
1060 int nquad =
m_base[0]->GetNumPoints();
1076 switch (
m_geom->GetCoordim())
1086 &gmat[0][0],1,dPhysValuesdx.get(),1,
1087 dPhysValuesdx.get(),1);
1091 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.get(), 1);
1099 PhysDeriv(physValues, dPhysValuesdx, dPhysValuesdy);
1105 &gmat[0][0], 1, dPhysValuesdx.get(), 1,
1106 dPhysValuesdx.get(), 1);
1108 &gmat[1][0], 1, dPhysValuesdy.get(), 1,
1109 dPhysValuesdx.get(), 1,
1110 dPhysValuesdx.get(), 1);
1114 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.get(), 1);
1116 gmat[1][0], dPhysValuesdy.get(), 1,
1117 dPhysValuesdx.get(), 1);
1127 dPhysValuesdy, dPhysValuesdz);
1133 &gmat[0][0], 1, dPhysValuesdx.get(), 1,
1134 dPhysValuesdx.get(), 1);
1136 &gmat[1][0], 1, dPhysValuesdy.get(), 1,
1137 dPhysValuesdx.get(), 1,
1138 dPhysValuesdx.get(), 1);
1140 &gmat[2][0], 1, dPhysValuesdz.get(), 1,
1141 dPhysValuesdx.get(), 1,
1142 dPhysValuesdx.get(), 1);
1146 Blas::Dscal(nquad, gmat[0][0], dPhysValuesdx.get(), 1);
1148 gmat[1][0], dPhysValuesdy.get(), 1,
1149 dPhysValuesdx.get(), 1);
1151 gmat[2][0], dPhysValuesdz.get(),
1152 1, dPhysValuesdx.get(), 1);
1157 ASSERTL0(
false,
"Wrong number of dimensions");
1162 Blas::Daxpy(
m_ncoeffs, lambda, wsp.get(), 1, outarray.get(), 1);
1194 return tmp->GetStdMatrix(mkey);
1206 "Geometric information is not set up");
1216 goto UseLocRegionsMatrix;
1221 goto UseStdRegionsMatrix;
1242 goto UseStdRegionsMatrix;
1254 goto UseLocRegionsMatrix;
1266 "Cannot call eWeakDeriv2 in a "
1267 "coordinate system which is not at "
1268 "least two-dimensional");
1273 "Cannot call eWeakDeriv2 in a "
1274 "coordinate system which is not "
1275 "three-dimensional");
1299 goto UseLocRegionsMatrix;
1303 int coordim =
m_geom->GetCoordim();
1305 for (
int i = 0; i < coordim; ++i)
1311 goto UseStdRegionsMatrix;
1327 int rows = LapMat.GetRows();
1328 int cols = LapMat.GetColumns();
1334 (*helm) = LapMat + factor*MassMat;
1378 coords[0] = (vertex == 0) ? -1.0 : 1.0;
1380 m_Ix =
m_base[0]->GetI(coords);
1386 UseLocRegionsMatrix:
1393 UseStdRegionsMatrix:
1424 returnval = StdSegExp::v_GenMatrix(mkey);
1439 "Geometric information is not set up");
1445 exp_size[0] = nbdry;
1461 goto UseLocRegionsMatrix;
1467 goto UseLocRegionsMatrix;
1472 factor = mat->Scale();
1473 goto UseStdRegionsMatrix;
1476 UseStdRegionsMatrix:
1484 returnval->SetBlock(0,0,Atmp =
1486 factor,Asubmat = mat->GetBlock(0,0)));
1487 returnval->SetBlock(0,1,Atmp =
1489 one,Asubmat = mat->GetBlock(0,1)));
1490 returnval->SetBlock(1,0,Atmp =
1492 factor,Asubmat = mat->GetBlock(1,0)));
1493 returnval->SetBlock(1,1,Atmp =
1495 invfactor,Asubmat = mat->GetBlock(1,1)));
1498 UseLocRegionsMatrix:
1518 for (i = 0; i < nbdry; ++i)
1520 for (j = 0; j < nbdry; ++j)
1522 (*A)(i,j) = mat(bmap[i],bmap[j]);
1525 for (j = 0; j < nint; ++j)
1527 (*B)(i,j) = mat(bmap[i],imap[j]);
1531 for (i = 0; i < nint; ++i)
1533 for (j = 0; j < nbdry; ++j)
1535 (*C)(i,j) = mat(imap[i],bmap[j]);
1538 for (j = 0; j < nint; ++j)
1540 (*D)(i,j) = mat(imap[i],imap[j]);
1549 (*A) = (*A) - (*B)*(*C);
1555 AllocateSharedPtr(factor,A));
1557 AllocateSharedPtr(one,B));
1559 AllocateSharedPtr(factor,C));
1561 AllocateSharedPtr(invfactor,D));
1587 ASSERTL1(&inarray[0] != &outarray[0],
1588 "inarray and outarray can not be the same");
1593 outarray[0] = inarray[1];
1594 outarray[1] = inarray[0];
1598 outarray[m] = sgn*inarray[m];
1606 outarray[m_ncoeffs-1-m] = inarray[m];
1610 ASSERTL0(
false,
"This basis is not allowed in this method");
const LibUtilities::PointsKeyVector GetPointsKeys() const
virtual DNekScalBlkMatSharedPtr v_GetLocStaticCondMatrix(const MatrixKey &mkey)
LibUtilities::ShapeType DetShapeType() const
This function returns the shape of the expansion domain.
NekDouble GetConstFactor(const ConstFactorType &factor) const
DNekMatSharedPtr GenMatrix(const StdMatrixKey &mkey)
#define ASSERTL0(condition, msg)
const NormalVector & GetEdgeNormal(const int edge) const
const LibUtilities::BasisSharedPtr & GetBasis(int dir) const
This function gets the shared point to basis in the dir direction.
const ConstFactorMap & GetConstFactors() const
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
virtual const LibUtilities::BasisSharedPtr & v_GetBasis(int dir) const
virtual void v_GetCoord(const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords)
const VarCoeffMap & GetVarCoeffs() const
std::vector< PointsKey > PointsKeyVector
MatrixType GetMatrixType() const
void v_DropLocStaticCondMatrix(const MatrixKey &mkey)
void MassMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdMatrixKey &mkey)
static Array< OneD, NekDouble > NullNekDouble1DArray
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
LibUtilities::NekManager< MatrixKey, DNekScalBlkMat, MatrixKey::opLess > m_staticCondMatrixManager
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
virtual NekDouble v_StdPhysEvaluate(const Array< OneD, const NekDouble > &Lcoord, const Array< OneD, const NekDouble > &physvals)
std::map< int, NormalVector > m_vertexNormals
virtual void v_HelmholtzMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Forward transform from physical quadrature space stored in inarray and evaluate the expansion coeffic...
int GetNumPoints(const int dir) const
This function returns the number of quadrature points in the dir direction.
virtual void v_PhysDeriv_n(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_dn)
Evaluate the derivative normal to a line: . The derivative is calculated performing the product ...
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
Principle Modified Functions .
void NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, Array< OneD, NekDouble > &outarray)
virtual void v_ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int mode_offset, NekDouble *coeffs)
Unpack data from input file assuming it comes from.
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Lagrange Polynomials using the Gauss points .
virtual void v_GetTracePhysVals(const int edge, const StdRegions::StdExpansionSharedPtr &EdgeExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
LibUtilities::ShapeType GetShapeType() const
std::map< ConstFactorType, NekDouble > ConstFactorMap
virtual int v_GetNcoeffs(void) const
virtual StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const
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.
SpatialDomains::GeometrySharedPtr m_geom
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Inner product of inarray over region with respect to the expansion basis (this)->_Base[0] and return ...
virtual void v_PhysDeriv_s(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_ds)
Evaluate the derivative along a line: . The derivative is calculated performing the product ...
boost::shared_ptr< StdSegExp > StdSegExpSharedPtr
boost::shared_ptr< DNekMat > DNekMatSharedPtr
virtual int v_GetNumPoints(const int dir) const
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)
virtual void v_NormVectorIProductWRTBase(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
1D Gauss-Gauss-Legendre quadrature points
Expansion2DSharedPtr GetLeftAdjacentElementExp() const
int NumBndryCoeffs(void) const
DNekBlkMatSharedPtr GetStdStaticCondMatrix(const StdMatrixKey &mkey)
int GetVertexMap(const int localVertexId, bool useCoeffPacking=false)
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
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...
virtual SpatialDomains::GeomType v_MetricInfoType()
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
virtual StdRegions::Orientation v_GetPorient(int point)
Principle Modified Functions .
virtual void v_ComputeVertexNormal(const int vertex)
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
boost::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
void GetInteriorMap(Array< OneD, unsigned int > &outarray)
Defines a specification for a set of points.
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals(void)
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
virtual void v_FwdTrans_BndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
virtual void v_LaplacianMatrixOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
T Ddot(int n, const Array< OneD, const T > &w, const int incw, const Array< OneD, const T > &x, const int incx, const Array< OneD, const int > &y, const int incy)
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
virtual DNekScalMatSharedPtr v_GetLocMatrix(const MatrixKey &mkey)
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.
virtual void v_GetVertexPhysVals(const int vertex, const Array< OneD, const NekDouble > &inarray, NekDouble &outarray)
boost::shared_ptr< DNekBlkMat > DNekBlkMatSharedPtr
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
SpatialDomains::GeometrySharedPtr GetGeom() const
virtual void v_SetCoeffsToOrientation(Array< OneD, NekDouble > &coeffs, StdRegions::Orientation dir)
boost::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
virtual int v_NumDGBndryCoeffs() const
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...
void MultiplyByElmtInvMass(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Geometry is straight-sided with constant geometric factors.
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 ...
virtual 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)
Evaluate the derivative at the physical quadrature points given by inarray and return in outarray...
boost::shared_ptr< Geometry1D > Geometry1DSharedPtr
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function performs the Backward transformation from coefficient space to physical space...
int GetLeftAdjacentElementEdge() const
virtual int v_GetCoordim()
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Integrate the physical point list inarray over region and return the value.
virtual DNekMatSharedPtr v_CreateStdMatrix(const StdRegions::StdMatrixKey &mkey)
GeomType
Indicates the type of element geometry.
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayofArray
boost::shared_ptr< Basis > BasisSharedPtr
void Zero(int n, T *x, const int incx)
Zero vector.
virtual void v_IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
LibUtilities::PointsType GetPointsType(const int dir) const
This function returns the type of quadrature points used in the dir direction.
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Array< OneD, LibUtilities::BasisSharedPtr > m_base
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Geometry is curved or has non-constant factors.
void GetBoundaryMap(Array< OneD, unsigned int > &outarray)
Describes the specification for a Basis.
virtual NekDouble v_PhysEvaluate(const Array< OneD, const NekDouble > &coord, const Array< OneD, const NekDouble > &physvals)
virtual int v_NumBndryCoeffs() const
1D Gauss-Lobatto-Legendre quadrature points
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 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.
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager