43 namespace LocalRegions
51 StdNodalTriExp(Ba,Bb,Ntype),
56 std::string(
"NodalTriExpMatrix")),
57 m_staticCondMatrixManager(
58 boost::bind(&
NodalTriExp::CreateStaticCondMatrix, this, _1),
59 std::string(
"NodalTriExpStaticCondMatrix"))
66 StdRegions::StdNodalTriExp(T),
69 m_matrixManager(T.m_matrixManager),
70 m_staticCondMatrixManager(T.m_staticCondMatrixManager)
101 int nquad0 =
m_base[0]->GetNumPoints();
102 int nquad1 =
m_base[1]->GetNumPoints();
110 Vmath::Vmul(nquad0*nquad1, jac, 1, inarray, 1,tmp, 1);
114 Vmath::Smul(nquad0*nquad1, jac[0], inarray, 1, tmp, 1);
118 ival = StdNodalTriExp::v_Integral(tmp);
125 bool multiplybyweights)
127 int nquad0 =
m_base[0]->GetNumPoints();
128 int nquad1 =
m_base[1]->GetNumPoints();
129 int order1 =
m_base[1]->GetNumModes();
131 if(multiplybyweights)
137 StdTriExp::IProductWRTBase_SumFacKernel(
m_base[0]->GetBdata(),
m_base[1]->GetBdata(),tmp,outarray,wsp);
144 StdTriExp::IProductWRTBase_SumFacKernel(
m_base[0]->GetBdata(),
m_base[1]->GetBdata(),inarray,outarray,wsp);
156 Blas::Dgemv(
'N',
m_ncoeffs,nq,iprodmat->Scale(),(iprodmat->GetOwnedMatrix())->GetPtr().get(),
157 m_ncoeffs, inarray.get(), 1, 0.0, outarray.get(), 1);
164 ASSERTL1((dir==0)||(dir==1)||(dir==2),
"Invalid direction.");
165 ASSERTL1((dir==2)?(
m_geom->GetCoordim()==3):
true,
"Invalid direction.");
168 int nquad0 =
m_base[0]->GetNumPoints();
169 int nquad1 =
m_base[1]->GetNumPoints();
170 int nqtot = nquad0*nquad1;
187 for(i = 0; i < nquad1; ++i)
189 gfac0[i] = 2.0/(1-z1[i]);
191 for(i = 0; i < nquad0; ++i)
193 gfac1[i] = 0.5*(1+z0[i]);
196 for(i = 0; i < nquad1; ++i)
198 Vmath::Smul(nquad0,gfac0[i],&inarray[0]+i*nquad0,1,&tmp0[0]+i*nquad0,1);
201 for(i = 0; i < nquad1; ++i)
203 Vmath::Vmul(nquad0,&gfac1[0],1,&tmp0[0]+i*nquad0,1,&tmp1[0]+i*nquad0,1);
208 Vmath::Vmul(nqtot,&df[2*dir][0], 1,&tmp0[0], 1,&tmp0[0],1);
209 Vmath::Vmul(nqtot,&df[2*dir+1][0],1,&tmp1[0], 1,&tmp1[0],1);
210 Vmath::Vmul(nqtot,&df[2*dir+1][0],1,&inarray[0],1,&tmp2[0],1);
214 Vmath::Smul(nqtot, df[2*dir][0], tmp0, 1, tmp0, 1);
215 Vmath::Smul(nqtot, df[2*dir+1][0], tmp1, 1, tmp1, 1);
216 Vmath::Smul(nqtot, df[2*dir+1][0], inarray, 1, tmp2, 1);
256 ASSERTL1(
false,
"input dir is out of range");
264 Blas::Dgemv(
'N',
m_ncoeffs,nq,iprodmat->Scale(),(iprodmat->GetOwnedMatrix())->GetPtr().get(),
265 m_ncoeffs, inarray.get(), 1, 0.0, outarray.get(), 1);
281 int nquad0 =
m_base[0]->GetNumPoints();
282 int nquad1 =
m_base[1]->GetNumPoints();
283 int nqtot = nquad0*nquad1;
290 StdNodalTriExp::v_PhysDeriv(inarray, diff0, diff1);
294 if(out_d0.num_elements())
297 Vmath::Vvtvp (nqtot,df[1],1,diff1,1, out_d0, 1, out_d0,1);
300 if(out_d1.num_elements())
303 Vmath::Vvtvp (nqtot,df[3],1,diff1,1, out_d1, 1, out_d1,1);
306 if(out_d2.num_elements())
309 Vmath::Vvtvp (nqtot,df[5],1,diff1,1, out_d2, 1, out_d2,1);
314 if(out_d0.num_elements())
316 Vmath::Smul (nqtot, df[0][0], diff0, 1, out_d0, 1);
317 Blas::Daxpy (nqtot, df[1][0], diff1, 1, out_d0, 1);
320 if(out_d1.num_elements())
322 Vmath::Smul (nqtot, df[2][0], diff0, 1, out_d1, 1);
323 Blas::Daxpy (nqtot, df[3][0], diff1, 1, out_d1, 1);
326 if(out_d2.num_elements())
328 Vmath::Smul (nqtot, df[4][0], diff0, 1, out_d2, 1);
329 Blas::Daxpy (nqtot, df[5][0], diff1, 1, out_d2, 1);
369 if(inarray.get() == outarray.get())
374 Blas::Dgemv(
'N',
m_ncoeffs,
m_ncoeffs,mat->Scale(),(mat->GetOwnedMatrix())->GetPtr().get(),
375 m_ncoeffs, tmp.get(), 1, 0.0, outarray.get(), 1);
379 Blas::Dgemv(
'N',
m_ncoeffs,
m_ncoeffs,mat->Scale(),(mat->GetOwnedMatrix())->GetPtr().get(),
380 m_ncoeffs, inarray.get(), 1, 0.0, outarray.get(), 1);
397 ASSERTL1(Lcoords[0] >= -1.0 && Lcoords[1] <= 1.0 &&
398 Lcoords[1] >= -1.0 && Lcoords[1] <=1.0,
399 "Local coordinates are not in region [-1,1]");
403 for(i = 0; i <
m_geom->GetCoordim(); ++i)
405 coords[i] =
m_geom->GetCoord(i,Lcoords);
417 return tmp->GetStdMatrix(mkey);
428 m_geom->GetLocCoords(coord,Lcoord);
430 return StdNodalTriExp::v_PhysEvaluate(Lcoord, physvals);
491 ASSERTL1(
m_geom->GetCoordim() == 2,
"Standard Region Laplacian is only set up for Quads in two-dimensional");
507 int rows = lap00.GetRows();
508 int cols = lap00.GetColumns();
512 (*lap) = gmat[0][0] * lap00 +
513 gmat[1][0] * (lap01 +
Transpose(lap01)) +
530 int rows = LapMat.GetRows();
531 int cols = LapMat.GetColumns();
536 (*helm) = LapMat + factor*MassMat;
557 unsigned int nint = (
unsigned int)(
m_ncoeffs - nbdry);
558 unsigned int exp_size[] = {nbdry,nint};
559 unsigned int nblks = 2;
570 goto UseLocRegionsMatrix;
576 goto UseLocRegionsMatrix;
581 factor = mat->Scale();
582 goto UseStdRegionsMatrix;
615 for(i = 0; i < nbdry; ++i)
617 for(j = 0; j < nbdry; ++j)
619 (*A)(i,j) = mat(bmap[i],bmap[j]);
622 for(j = 0; j < nint; ++j)
624 (*B)(i,j) = mat(bmap[i],imap[j]);
628 for(i = 0; i < nint; ++i)
630 for(j = 0; j < nbdry; ++j)
632 (*C)(i,j) = mat(imap[i],bmap[j]);
635 for(j = 0; j < nint; ++j)
637 (*D)(i,j) = mat(imap[i],imap[j]);
646 (*A) = (*A) - (*B)*(*C);
675 2,
m_base[0]->GetPointsKey());
677 2,
m_base[1]->GetPointsKey());
698 returnval = StdNodalTriExp::v_GenMatrix(mkey);
713 int nqe =
m_base[0]->GetNumPoints();
718 for (i = 0; i < dim; ++i)
739 Vmath::Fill(nqe,df[2*i+1][0] + df[2*i][0],normal[i],1);
749 ASSERTL0(
false,
"Edge is out of range (edge < 3)");
756 fac += normal[i][0]*normal[i][0];
768 int nquad0 = ptsKeys[0].GetNumPoints();
769 int nquad1 = ptsKeys[1].GetNumPoints();
782 for(j = 0; j < nquad0; ++j)
787 normals[i*nquad0+j] = -df[2*i+1][j]*edgejac[j];
790 from_key = ptsKeys[0];
793 for(j = 0; j < nquad1; ++j)
795 edgejac[j] = jac[nquad0*j+nquad0-1];
798 normals[i*nquad1+j] = (df[2*i][nquad0*j + nquad0-1] + df[2*i+1][nquad0*j + nquad0-1])*edgejac[j];
801 from_key = ptsKeys[1];
804 for(j = 0; j < nquad1; ++j)
806 edgejac[j] = jac[nquad0*j];
809 normals[i*nquad1+j] = -df[2*i][nquad0*j]*edgejac[j];
812 from_key = ptsKeys[1];
815 ASSERTL0(
false,
"edge is out of range (edge < 3)");
837 Vmath::Vvtvp(nqe,normal[i],1, normal[i],1,work,1,work,1);
const LibUtilities::PointsKeyVector GetPointsKeys() const
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)
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mod...
std::vector< PointsKey > PointsKeyVector
MatrixType GetMatrixType() const
NodalTriExp(const LibUtilities::BasisKey &Ba, const LibUtilities::BasisKey &Bb, const LibUtilities::PointsType Ntype, const SpatialDomains::TriGeomSharedPtr &geom)
Constructor using BasisKey class for quadrature points and order definition.
void v_ComputeEdgeNormal(const int edge)
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
void MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void GetCoord(const Array< OneD, const NekDouble > &Lcoords, Array< OneD, NekDouble > &coords)
void IProductWRTBase_SumFac(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool multiplybyweights=true)
void IProductWRTDerivBase_SumFac(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
General purpose memory allocation routines with the ability to allocate from thread specific memory p...
virtual StdRegions::StdExpansionSharedPtr v_GetStdExp(void) const
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
void IProductWRTBase_MatOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
void PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2=NullNekDouble1DArray)
Differentiation Methods.
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
std::map< int, StdRegions::NormalVector > m_edgeNormals
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
void Sdiv(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha/y.
LibUtilities::ShapeType GetShapeType() const
SpatialDomains::GeometrySharedPtr m_geom
boost::shared_ptr< DNekMat > DNekMatSharedPtr
void GeneralMatrixOp_MatOp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const StdRegions::StdMatrixKey &mkey)
boost::shared_ptr< StdNodalTriExp > StdNodalTriExpSharedPtr
NekDouble PhysEvaluate(const Array< OneD, const NekDouble > &coord, const Array< OneD, const NekDouble > &physvals)
DNekMatSharedPtr GetStdMatrix(const StdMatrixKey &mkey)
boost::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
void NodalToModalTranspose(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
LibUtilities::NekManager< MatrixKey, DNekScalMat, MatrixKey::opLess > m_matrixManager
int GetTotPoints() const
This function returns the total number of quadrature points used in the element.
int NumBndryCoeffs(void) const
DNekBlkMatSharedPtr GetStdStaticCondMatrix(const StdMatrixKey &mkey)
LibUtilities::PointsKey m_nodalPointsKey
void Reverse(int n, const T *x, const int incx, T *y, const int incy)
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
virtual StdRegions::StdExpansionSharedPtr v_GetLinStdExp(void) const
int getNumberOfCoefficients(int Na)
void IProductWRTDerivBase_MatOp(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
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)
NekMatrix< InnerMatrixType, BlockMatrixTag > Transpose(NekMatrix< InnerMatrixType, BlockMatrixTag > &rhs)
Defines a specification for a set of points.
void 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...
~NodalTriExp()
Destructor.
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
virtual DNekMatSharedPtr v_GenMatrix(const StdRegions::StdMatrixKey &mkey)
LibUtilities::BasisType GetBasisType(const int dir) const
This function returns the type of basis used in the dir direction.
NekDouble Integral(const Array< OneD, const NekDouble > &inarray)
Integrate the physical point list inarray over region.
boost::shared_ptr< DNekBlkMat > DNekBlkMatSharedPtr
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
SpatialDomains::GeometrySharedPtr GetGeom() const
boost::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
#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 ...
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
unsigned int GetNumPoints() const
boost::shared_ptr< TriGeom > TriGeomSharedPtr
void GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3=NullNekDouble1DArray)
void IProductWRTBase_SumFacKernel(const Array< OneD, const NekDouble > &base0, const Array< OneD, const NekDouble > &base1, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wsp, bool doCheckCollDir0=true, bool doCheckCollDir1=true)
GeomType
Indicates the type of element geometry.
void Zero(int n, T *x, const int incx)
Zero vector.
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
DNekScalMatSharedPtr CreateMatrix(const MatrixKey &mkey)
#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 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 ...
void GetBoundaryMap(Array< OneD, unsigned int > &outarray)
Describes the specification for a Basis.
PointsType GetPointsType() const
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.
DNekMatSharedPtr CreateStdMatrix(const StdRegions::StdMatrixKey &mkey)
static VarCoeffMap NullVarCoeffMap
static ConstFactorMap NullConstFactorMap