49 namespace SpatialDomains
62 TriGeom::TriGeom(
int id,
const int coordim):
96 "Cannot call function with dim == 1");
108 Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
134 m_coordim = edges[0]->GetVertex(0)->GetCoordim();
150 :
Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
177 m_coordim = edges[0]->GetVertex(0)->GetCoordim();
196 std::list<CompToElmt>::const_iterator def;
205 for (
int i = 0; i <
kNedges; i++)
229 "Geometry is not in physical space");
234 return m_xmap->PhysEvaluate(Lcoord, tmp);
248 int i, j, vmap[3] = {-1,-1,-1};
249 NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
255 for (i = 0; i < 3; ++i)
257 cx += (*face2[i])(0) - (*face1[i])(0);
258 cy += (*face2[i])(1) - (*face1[i])(1);
259 cz += (*face2[i])(2) - (*face1[i])(2);
268 for (i = 0; i < 3; ++i)
273 for (j = 0; j < 3; ++j)
275 x1 = (*face2[j])(0)-cx;
276 y1 = (*face2[j])(1)-cy;
277 z1 = (*face2[j])(2)-cz;
278 if (sqrt((x1-x)*(x1-x)+(y1-y)*(y1-y)+(z1-z)*(z1-z)) < 1e-8)
286 if (vmap[1] == (vmap[0]+1) % 3)
317 ASSERTL0(
false,
"Unable to determine triangle orientation");
346 std::list<CompToElmt>::const_iterator def;
379 return m_xmap->GetBasis(i);
388 ASSERTL1(i <= 2,
"edge is out of range");
392 return m_xmap->GetBasis(0);
396 return m_xmap->GetBasis(1);
424 if(
m_xmap->GetBasisNumModes(0) != 2 ||
425 m_xmap->GetBasisNumModes(1) != 2)
458 int nEdgeCoeffs =
m_xmap->GetEdgeNcoeffs(0);
471 int N =
m_curve->m_points.size();
473 -1+(int)sqrt(static_cast<NekDouble>(8*N+1)))/2;
475 ASSERTL0(nEdgePts*(nEdgePts+1)/2 == N,
476 "NUMPOINTS should be a triangle number for"
484 "Number of edge points does not correspond "
485 "to number of face points in triangle "
498 max(nEdgePts*nEdgePts,
m_xmap->GetTotPoints()));
508 for (j = 0; j < N; ++j)
510 phys[j] = (
m_curve->m_points[j]->GetPtr())[i];
513 t->BwdTrans(phys, tmp);
518 m_xmap->GetBasis(0)->GetPointsKey(),
519 m_xmap->GetBasis(1)->GetPointsKey(),
529 int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
538 ASSERTL0(nEdgePts * nEdgePts == npts,
539 "NUMPOINTS should be a square number for"
547 "Number of edge points does not correspond to "
548 "number of face points in triangle "
554 for (j = 0; j <
npts; ++j)
556 tmp[j] = (
m_curve->m_points[j]->GetPtr())[i];
562 curveKey, curveKey, tmp,
563 m_xmap->GetBasis(0)->GetPointsKey(),
564 m_xmap->GetBasis(1)->GetPointsKey(),
573 ASSERTL0(
false,
"Only 1D/2D points distributions "
586 nEdgeCoeffs =
m_edges[i]->GetXmap()->GetNcoeffs();
590 for(k = 0; k < nEdgeCoeffs; k++)
593 signArray[k] *
m_edges[i]->GetCoeffs(j)[k];
627 orth1.
Mult(norm,dv1);
628 orth2.
Mult(norm,dv2);
635 Lcoords[0] = xin.
dot(orth2)/dv1.
dot(orth2);
636 Lcoords[0] = 2*Lcoords[0]-1;
637 Lcoords[1] = xin.
dot(orth1)/dv2.
dot(orth1);
638 Lcoords[1] = 2*Lcoords[1]-1;
661 Lcoords[0] = za[min_i%za.num_elements()];
662 Lcoords[1] = zb[min_i/za.num_elements()];
665 Lcoords[0] = (1.0+Lcoords[0])*(1.0-Lcoords[1])/2 -1.0;
679 ASSERTL2((i >=0) && (i <= 2),
"Edge id must be between 0 and 2");
689 ASSERTL2((i >=0) && (i <= 2),
"Vertex id must be between 0 and 2");
699 ASSERTL2((i >=0) && (i <= 2),
"Vertex id must be between 0 and 2");
709 ASSERTL2((i >=0) && (i <= 2),
"Edge id must be between 0 and 3");
719 ASSERTL2((i >=0) && (i <= 2),
"Edge id must be between 0 and 2");
729 ASSERTL2((i >=0) && (i <= 3),
"Edge id must be between 0 and 3");
758 for (i=0,edgeIter =
m_edges.begin(); edgeIter !=
m_edges.end(); ++edgeIter,++i)
760 if (*edgeIter == edge)
817 ASSERTL1(gloCoord.num_elements() >= 2,
818 "Two dimensional geometry expects at least two coordinates.");
821 if (stdCoord[0] >= -(1+tol) && stdCoord[1] >= -(1+tol)
822 && stdCoord[0] + stdCoord[1] <= tol)
836 if (it != curvedFaces.end())
841 for (
int i = 0; i < 3; ++i)
843 m_edges[i]->Reset(curvedEdges, curvedFaces);
852 int order0 =
m_edges[0]->GetBasis(0)->GetNumModes();
853 int order1 = max(order0,
StdRegions::StdExpansionSharedPtr m_xmap
virtual int v_NumElmtConnected() const
#define ASSERTL0(condition, msg)
virtual int v_GetFid() const
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
std::vector< PointGeomSharedPtr > PointGeomVector
GeomFactorsSharedPtr m_geomFactors
Structure holding graphvertexobject id and local element facet id.
NekDouble dot(PointGeom &a)
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
virtual void v_GenGeomFactors()
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 Sub(PointGeom &a, PointGeom &b)
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.
StdRegions::StdExpansionSharedPtr GetXmap() const
GeomState m_geomFactorsState
boost::shared_ptr< StdNodalTriExp > StdNodalTriExpSharedPtr
virtual void v_FillGeom()
Put all quadrature information into edge structure.
boost::shared_ptr< Curve > CurveSharedPtr
Gauss Radau pinned at x=-1, .
StdRegions::Orientation m_eorient[kNedges]
const LibUtilities::BasisSharedPtr GetBasis(const int i)
Return the j-th basis of the i-th co-ordinate dimension.
Principle Orthogonal Functions .
void Interp2D(const BasisKey &fbasis0, const BasisKey &fbasis1, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, Array< OneD, NekDouble > &to)
this function interpolates a 2D function evaluated at the quadrature points of the 2D basis...
virtual void v_SetOwnData()
virtual PointGeomSharedPtr v_GetVertex(int i) const
virtual StdRegions::Orientation v_GetEorient(const int i) const
boost::shared_ptr< SegGeom > SegGeomSharedPtr
Principle Modified Functions .
std::list< CompToElmt > m_elmtMap
virtual const LibUtilities::BasisSharedPtr v_GetEdgeBasis(const int i)
void NewtonIterationForLocCoord(const Array< OneD, const NekDouble > &coords, const Array< OneD, const NekDouble > &ptsx, const Array< OneD, const NekDouble > &ptsy, Array< OneD, NekDouble > &Lcoords, NekDouble &resid)
virtual NekDouble v_GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
virtual void v_AddElmtConnected(int gvo_id, int locid)
PointsManagerT & PointsManager(void)
Principle Orthogonal Functions .
Defines a specification for a set of points.
virtual int v_GetVid(int i) const
virtual int v_GetCoordim() const
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
virtual NekDouble v_GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
PointGeomSharedPtr GetVertex(int i) const
static const int kNedges
Get the orientation of face1.
virtual const LibUtilities::BasisSharedPtr v_GetBasis(const int i)
void Mult(PointGeom &a, PointGeom &b)
NekDouble GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
virtual StdRegions::Orientation v_GetCartesianEorient(const int i) const
virtual int v_WhichEdge(SegGeomSharedPtr edge)
Return the edge number of the given edge.
Array< OneD, Array< OneD, NekDouble > > m_coeffs
#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.
virtual int v_GetEid(int i) const
boost::shared_ptr< Geometry1D > Geometry1DSharedPtr
NekDouble GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
LibUtilities::ShapeType m_shapeType
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
static StdRegions::Orientation GetFaceOrientation(const TriGeom &face1, const TriGeom &face2)
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Determines if a point specified in global coordinates is located within this tetrahedral geometry...
boost::unordered_map< int, CurveSharedPtr > CurveMap
virtual int v_GetNumEdges() const
Geometric information has been generated.
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
GeomFactorsSharedPtr GetMetricInfo()
virtual int v_GetNumVerts() const
GeomType
Indicates the type of element geometry.
boost::shared_ptr< Basis > BasisSharedPtr
virtual bool v_IsElmtConnected(int gvo_id, int locid) const
GeomState m_state
enum identifier to determine if quad points are filled
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Geometry is curved or has non-constant factors.
int m_coordim
coordinate dimension
Describes the specification for a Basis.
virtual const Geometry1DSharedPtr v_GetEdge(int i) const
boost::shared_ptr< PointGeom > PointGeomSharedPtr
1D Gauss-Lobatto-Legendre quadrature points
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.
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.