Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes | List of all members
Nektar::SpatialDomains::QuadGeom Class Reference

#include <QuadGeom.h>

Inheritance diagram for Nektar::SpatialDomains::QuadGeom:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SpatialDomains::QuadGeom:
Collaboration graph
[legend]

Public Member Functions

 QuadGeom ()
 
 QuadGeom (const int id, const PointGeomSharedPtr verts[], const SegGeomSharedPtr edges[], const StdRegions::Orientation eorient[])
 
 QuadGeom (const int id, const SegGeomSharedPtr edges[], const StdRegions::Orientation eorient[])
 
 QuadGeom (const int id, const SegGeomSharedPtr edges[], const StdRegions::Orientation eorient[], const CurveSharedPtr &curve)
 
 QuadGeom (const QuadGeom &in)
 
 ~QuadGeom ()
 
NekDouble GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry2D
 Geometry2D ()
 
 Geometry2D (const int coordim)
 
virtual ~Geometry2D ()
 
int GetFid () const
 
const Geometry1DSharedPtr GetEdge (int i) const
 
const Geometry2DSharedPtr GetFace (int i) const
 
int WhichEdge (SegGeomSharedPtr edge)
 
int WhichFace (Geometry2DSharedPtr face)
 
const LibUtilities::BasisSharedPtr GetEdgeBasis (const int i)
 
StdRegions::Orientation GetFaceOrient (const int i) const
 
StdRegions::Orientation GetCartesianEorient (const int i) const
 
- Public Member Functions inherited from Nektar::SpatialDomains::Geometry
 Geometry ()
 
 Geometry (int coordim)
 
virtual ~Geometry ()
 
bool IsElmtConnected (int gvo_id, int locid) const
 
void AddElmtConnected (int gvo_id, int locid)
 
int NumElmtConnected () const
 
int GetCoordim () const
 
void SetCoordim (int coordim)
 
GeomFactorsSharedPtr GetGeomFactors ()
 
GeomFactorsSharedPtr GetRefGeomFactors (const Array< OneD, const LibUtilities::BasisSharedPtr > &tbasis)
 
GeomFactorsSharedPtr GetMetricInfo ()
 
LibUtilities::ShapeType GetShapeType (void)
 
int GetGlobalID (void)
 
void SetGlobalID (int globalid)
 
int GetVid (int i) const
 
int GetEid (int i) const
 
int GetFid (int i) const
 
int GetTid (int i) const
 
int GetNumVerts () const
 
PointGeomSharedPtr GetVertex (int i) const
 
StdRegions::Orientation GetEorient (const int i) const
 
StdRegions::Orientation GetPorient (const int i) const
 
StdRegions::Orientation GetForient (const int i) const
 
int GetNumEdges () const
 
int GetNumFaces () const
 
int GetShapeDim () const
 
StdRegions::StdExpansionSharedPtr GetXmap () const
 
const Array< OneD, const
NekDouble > & 
GetCoeffs (const int i) const
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
 
bool ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &resid)
 
int GetVertexEdgeMap (int i, int j) const
 
int GetVertexFaceMap (int i, int j) const
 return the id of the $j^{th}$ face attached to the $ i^{th}$ vertex More...
 
int GetEdgeFaceMap (int i, int j) const
 
void FillGeom ()
 Put all quadrature information into face/edge structure and backward transform. More...
 
NekDouble GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 
NekDouble GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 Given local collapsed coordinate Lcoord return the value of physical coordinate in direction i. More...
 
void SetOwnData ()
 
const LibUtilities::BasisSharedPtr GetBasis (const int i)
 Return the j-th basis of the i-th co-ordinate dimension. More...
 
const LibUtilities::PointsKeyVector GetPointsKeys ()
 
void Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 

Static Public Member Functions

static StdRegions::Orientation GetFaceOrientation (const QuadGeom &face1, const QuadGeom &face2)
 Get the orientation of face1. More...
 
static StdRegions::Orientation GetFaceOrientation (const PointGeomVector &face1, const PointGeomVector &face2)
 

Static Public Attributes

static const int kNverts = 4
 
static const int kNedges = 4
 
static const std::string XMLElementType
 

Protected Member Functions

virtual void v_AddElmtConnected (int gvo_id, int locid)
 
virtual int v_NumElmtConnected () const
 
virtual bool v_IsElmtConnected (int gvo_id, int locid) const
 
virtual int v_GetFid () const
 
virtual int v_GetCoordim () const
 
virtual const
LibUtilities::BasisSharedPtr 
v_GetBasis (const int i)
 
virtual const
LibUtilities::BasisSharedPtr 
v_GetEdgeBasis (const int i)
 
virtual NekDouble v_GetCoord (const int i, const Array< OneD, const NekDouble > &Lcoord)
 
void v_GenGeomFactors ()
 
virtual void v_SetOwnData ()
 
virtual void v_FillGeom ()
 Put all quadrature information into edge structure. More...
 
virtual NekDouble v_GetLocCoords (const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
 
virtual int v_GetEid (int i) const
 
virtual int v_GetVid (int i) const
 
virtual PointGeomSharedPtr v_GetVertex (int i) const
 
virtual const Geometry1DSharedPtr v_GetEdge (int i) const
 
virtual StdRegions::Orientation v_GetEorient (const int i) const
 
virtual StdRegions::Orientation v_GetCartesianEorient (const int i) const
 
virtual int v_WhichEdge (SegGeomSharedPtr edge)
 Return the edge number of the given edge. More...
 
virtual int v_GetNumVerts () const
 
virtual int v_GetNumEdges () const
 
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
 
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol)
 
virtual bool v_ContainsPoint (const Array< OneD, const NekDouble > &gloCoord, Array< OneD, NekDouble > &locCoord, NekDouble tol, NekDouble &resid)
 
virtual void v_Reset (CurveMap &curvedEdges, CurveMap &curvedFaces)
 Reset this geometry object: unset the current state and remove allocated GeomFactors. More...
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry2D
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)
 
- Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
void GenGeomFactors ()
 
virtual int v_GetFid (int i) const
 
virtual StdRegions::Orientation v_GetPorient (const int i) const
 
virtual StdRegions::Orientation v_GetForient (const int i) const
 
virtual int v_GetNumFaces () const
 
virtual
StdRegions::StdExpansionSharedPtr 
v_GetXmap () const
 
virtual int v_GetVertexEdgeMap (int i, int j) const
 
virtual int v_GetVertexFaceMap (int i, int j) const
 
virtual int v_GetEdgeFaceMap (int i, int j) const
 
void SetUpCoeffs (const int nCoeffs)
 Initialise the m_coeffs array. More...
 

Protected Attributes

PointGeomVector m_verts
 
SegGeomVector m_edges
 
StdRegions::Orientation m_eorient [kNedges]
 
int m_fid
 
bool m_ownVerts
 
std::list< CompToElmtm_elmtMap
 
CurveSharedPtr m_curve
 
- Protected Attributes inherited from Nektar::SpatialDomains::Geometry
int m_coordim
 coordinate dimension More...
 
GeomFactorsSharedPtr m_geomFactors
 
GeomState m_geomFactorsState
 
StdRegions::StdExpansionSharedPtr m_xmap
 
GeomState m_state
 enum identifier to determine if quad points are filled More...
 
GeomType m_geomType
 
LibUtilities::ShapeType m_shapeType
 
int m_globalID
 
Array< OneD, Array< OneD,
NekDouble > > 
m_coeffs
 

Private Member Functions

void SetUpXmap ()
 

Private Attributes

bool m_ownData
 Boolean indicating whether object owns the data. More...
 

Additional Inherited Members

- Static Protected Member Functions inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsSharedPtr ValidateRegGeomFactor (GeomFactorsSharedPtr geomFactor)
 
- Static Protected Attributes inherited from Nektar::SpatialDomains::Geometry
static GeomFactorsVector m_regGeomFactorsManager
 

Detailed Description

Definition at line 60 of file QuadGeom.h.

Constructor & Destructor Documentation

Nektar::SpatialDomains::QuadGeom::QuadGeom ( )
Nektar::SpatialDomains::QuadGeom::QuadGeom ( const int  id,
const PointGeomSharedPtr  verts[],
const SegGeomSharedPtr  edges[],
const StdRegions::Orientation  eorient[] 
)

Copy the vert shared pointers.

Copy the edge shared pointers.

Definition at line 63 of file QuadGeom.cpp.

References ASSERTL0, Nektar::LibUtilities::eQuadrilateral, kNedges, kNverts, Nektar::SpatialDomains::Geometry::m_coordim, m_edges, m_eorient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), and SetUpXmap().

66  :
67  Geometry2D(verts[0]->GetCoordim()),
68  m_fid(id)
69  {
70  m_globalID = id;
72 
73  /// Copy the vert shared pointers.
74  m_verts.insert(m_verts.begin(), verts, verts+QuadGeom::kNverts);
75 
76  /// Copy the edge shared pointers.
77  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
78 
79 
80  for (int j=0; j<kNedges; ++j)
81  {
82  m_eorient[j] = eorient[j];
83  }
84 
85  m_coordim = verts[0]->GetCoordim();
86  ASSERTL0(m_coordim > 1,
87  "Cannot call function with dim == 1");
88 
89  SetUpXmap();
90  SetUpCoeffs(m_xmap->GetNcoeffs());
91  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
int m_coordim
coordinate dimension
Definition: Geometry.h:169
Nektar::SpatialDomains::QuadGeom::QuadGeom ( const int  id,
const SegGeomSharedPtr  edges[],
const StdRegions::Orientation  eorient[] 
)

Copy the edge shared pointers.

Definition at line 143 of file QuadGeom.cpp.

References ASSERTL0, Nektar::StdRegions::eForwards, Nektar::LibUtilities::eQuadrilateral, Nektar::SpatialDomains::Geometry::GetVertex(), kNedges, Nektar::SpatialDomains::Geometry::m_coordim, m_edges, m_eorient, m_fid, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), and SetUpXmap().

145  :
146  Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
147  m_fid(id)
148  {
149  int j;
150 
151  m_globalID = m_fid;
153 
154  /// Copy the edge shared pointers.
155  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
156 
157  for(j=0; j <kNedges; ++j)
158  {
159  if(eorient[j] == StdRegions::eForwards)
160  {
161  m_verts.push_back(edges[j]->GetVertex(0));
162  }
163  else
164  {
165  m_verts.push_back(edges[j]->GetVertex(1));
166  }
167  }
168 
169  for (j=0; j<kNedges; ++j)
170  {
171  m_eorient[j] = eorient[j];
172  }
173 
174  m_coordim = edges[0]->GetVertex(0)->GetCoordim();
175  ASSERTL0(m_coordim > 1,
176  "Cannot call function with dim == 1");
177 
178  SetUpXmap();
179  SetUpCoeffs(m_xmap->GetNcoeffs());
180  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
int m_coordim
coordinate dimension
Definition: Geometry.h:169
Nektar::SpatialDomains::QuadGeom::QuadGeom ( const int  id,
const SegGeomSharedPtr  edges[],
const StdRegions::Orientation  eorient[],
const CurveSharedPtr curve 
)

Copy the edge shared pointers.

Definition at line 97 of file QuadGeom.cpp.

References ASSERTL0, Nektar::StdRegions::eForwards, Nektar::LibUtilities::eQuadrilateral, Nektar::SpatialDomains::Geometry::GetVertex(), kNedges, Nektar::SpatialDomains::Geometry::m_coordim, m_edges, m_eorient, m_fid, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_shapeType, m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), and SetUpXmap().

100  :
101  Geometry2D(edges[0]->GetVertex(0)->GetCoordim()),
102  m_fid(id),
103  m_curve(curve)
104  {
105  int j;
106 
107  m_globalID = m_fid;
108 
110 
111  /// Copy the edge shared pointers.
112  m_edges.insert(m_edges.begin(), edges, edges+QuadGeom::kNedges);
113 
114  for(j=0; j <kNedges; ++j)
115  {
116  if(eorient[j] == StdRegions::eForwards)
117  {
118  m_verts.push_back(edges[j]->GetVertex(0));
119  }
120  else
121  {
122  m_verts.push_back(edges[j]->GetVertex(1));
123  }
124  }
125 
126  for (j=0; j<kNedges; ++j)
127  {
128  m_eorient[j] = eorient[j];
129  }
130 
131  m_coordim = edges[0]->GetVertex(0)->GetCoordim();
132  ASSERTL0(m_coordim > 1,
133  "Cannot call function with dim == 1");
134 
135  SetUpXmap();
136  SetUpCoeffs(m_xmap->GetNcoeffs());
137  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
PointGeomSharedPtr GetVertex(int i) const
Definition: Geometry.h:348
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
int m_coordim
coordinate dimension
Definition: Geometry.h:169
Nektar::SpatialDomains::QuadGeom::QuadGeom ( const QuadGeom in)

Definition at line 186 of file QuadGeom.cpp.

References kNedges, m_edges, m_elmtMap, m_eorient, m_fid, Nektar::SpatialDomains::Geometry::m_globalID, m_ownData, m_ownVerts, Nektar::SpatialDomains::Geometry::m_shapeType, and m_verts.

187  {
188  // From Geometry
189  m_shapeType = in.m_shapeType;
190 
191  // From QuadFaceComponent
192  m_fid = in.m_fid;
193  m_globalID = m_fid;
194  m_ownVerts = in.m_ownVerts;
195  std::list<CompToElmt>::const_iterator def;
196  for(def = in.m_elmtMap.begin(); def != in.m_elmtMap.end(); def++)
197  {
198  m_elmtMap.push_back(*def);
199  }
200 
201  // From QuadGeom
202  m_verts = in.m_verts;
203  m_edges = in.m_edges;
204  for (int i = 0; i < kNedges; i++)
205  {
206  m_eorient[i] = in.m_eorient[i];
207  }
208  m_ownData = in.m_ownData;
209  }
std::list< CompToElmt > m_elmtMap
Definition: QuadGeom.h:110
LibUtilities::ShapeType m_shapeType
Definition: Geometry.h:177
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
bool m_ownData
Boolean indicating whether object owns the data.
Definition: QuadGeom.h:193
Nektar::SpatialDomains::QuadGeom::~QuadGeom ( )

Definition at line 215 of file QuadGeom.cpp.

216  {
217  }

Member Function Documentation

NekDouble Nektar::SpatialDomains::QuadGeom::GetCoord ( const int  i,
const Array< OneD, const NekDouble > &  Lcoord 
)

Definition at line 242 of file QuadGeom.cpp.

References ASSERTL1, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_state, and Nektar::SpatialDomains::Geometry::m_xmap.

Referenced by v_GetCoord().

244  {
246  "Geometry is not in physical space");
247 
248  Array<OneD, NekDouble> tmp(m_xmap->GetTotPoints());
249  m_xmap->BwdTrans(m_coeffs[i], tmp);
250 
251  return m_xmap->PhysEvaluate(Lcoord, tmp);
252  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometric information has been generated.
GeomState m_state
enum identifier to determine if quad points are filled
Definition: Geometry.h:175
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::GetFaceOrientation ( const QuadGeom face1,
const QuadGeom face2 
)
static

Get the orientation of face1.

Definition at line 254 of file QuadGeom.cpp.

References m_verts.

Referenced by Nektar::MultiRegions::DisContField3D::FindPeriodicFaces().

257  {
258  return GetFaceOrientation(face1.m_verts, face2.m_verts);
259  }
static StdRegions::Orientation GetFaceOrientation(const QuadGeom &face1, const QuadGeom &face2)
Get the orientation of face1.
Definition: QuadGeom.cpp:254
StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::GetFaceOrientation ( const PointGeomVector face1,
const PointGeomVector face2 
)
static

Calculate the orientation of face2 to face1 (note this is not face1 to face2!).

Definition at line 265 of file QuadGeom.cpp.

References ASSERTL0, Nektar::StdRegions::eDir1BwdDir1_Dir2BwdDir2, Nektar::StdRegions::eDir1BwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1BwdDir2_Dir2BwdDir1, Nektar::StdRegions::eDir1BwdDir2_Dir2FwdDir1, Nektar::StdRegions::eDir1FwdDir1_Dir2BwdDir2, Nektar::StdRegions::eDir1FwdDir1_Dir2FwdDir2, Nektar::StdRegions::eDir1FwdDir2_Dir2BwdDir1, and Nektar::StdRegions::eDir1FwdDir2_Dir2FwdDir1.

268  {
269  int i, j, vmap[4] = {-1,-1,-1,-1};
270  NekDouble x, y, z, x1, y1, z1, cx = 0.0, cy = 0.0, cz = 0.0;
271 
272  // For periodic faces, we calculate the vector between the centre
273  // points of the two faces. (For connected faces this will be
274  // zero). We can then use this to determine alignment later in the
275  // algorithm.
276  for (i = 0; i < 4; ++i)
277  {
278  cx += (*face2[i])(0) - (*face1[i])(0);
279  cy += (*face2[i])(1) - (*face1[i])(1);
280  cz += (*face2[i])(2) - (*face1[i])(2);
281  }
282  cx /= 4;
283  cy /= 4;
284  cz /= 4;
285 
286  // Now construct a mapping which takes us from the vertices of one
287  // face to the other. That is, vertex j of face2 corresponds to
288  // vertex vmap[j] of face1.
289  for (i = 0; i < 4; ++i)
290  {
291  x = (*face1[i])(0);
292  y = (*face1[i])(1);
293  z = (*face1[i])(2);
294  for (j = 0; j < 4; ++j)
295  {
296  x1 = (*face2[j])(0)-cx;
297  y1 = (*face2[j])(1)-cy;
298  z1 = (*face2[j])(2)-cz;
299  if (sqrt((x1-x)*(x1-x)+(y1-y)*(y1-y)+(z1-z)*(z1-z)) < 1e-8)
300  {
301  vmap[j] = i;
302  break;
303  }
304  }
305  }
306 
307  // Use the mapping to determine the eight alignment options between
308  // faces.
309  if (vmap[1] == (vmap[0]+1) % 4)
310  {
311  switch (vmap[0])
312  {
313  case 0:
315  break;
316  case 1:
318  break;
319  case 2:
321  break;
322  case 3:
324  break;
325  }
326  }
327  else
328  {
329  switch (vmap[0])
330  {
331  case 0:
333  break;
334  case 1:
336  break;
337  case 2:
339  break;
340  case 3:
342  break;
343  }
344  }
345 
346  ASSERTL0(false, "unable to determine face orientation");
348  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
double NekDouble
void Nektar::SpatialDomains::QuadGeom::SetUpXmap ( )
private

Definition at line 219 of file QuadGeom.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::LibUtilities::eGaussLobattoLegendre, Nektar::LibUtilities::eModified_A, Nektar::SpatialDomains::Geometry::GetBasis(), m_edges, and Nektar::SpatialDomains::Geometry::m_xmap.

Referenced by QuadGeom(), and v_Reset().

220  {
221  int order0 = max(m_edges[0]->GetBasis(0)->GetNumModes(),
222  m_edges[2]->GetBasis(0)->GetNumModes());
223  int order1 = max(m_edges[1]->GetBasis(0)->GetNumModes(),
224  m_edges[3]->GetBasis(0)->GetNumModes());
225 
226  const LibUtilities::BasisKey B0(
228  LibUtilities::PointsKey(
230  const LibUtilities::BasisKey B1(
232  LibUtilities::PointsKey(
234 
236  ::AllocateSharedPtr(B0,B1);
237  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Principle Modified Functions .
Definition: BasisType.h:49
const LibUtilities::BasisSharedPtr GetBasis(const int i)
Return the j-th basis of the i-th co-ordinate dimension.
Definition: Geometry.h:475
1D Gauss-Lobatto-Legendre quadrature points
Definition: PointsType.h:50
void Nektar::SpatialDomains::QuadGeom::v_AddElmtConnected ( int  gvo_id,
int  locid 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 354 of file QuadGeom.cpp.

References m_elmtMap.

355  {
356  CompToElmt ee(gvo_id,locid);
357  m_elmtMap.push_back(ee);
358  }
std::list< CompToElmt > m_elmtMap
Definition: QuadGeom.h:110
bool Nektar::SpatialDomains::QuadGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
NekDouble  tol = 0.0 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 789 of file QuadGeom.cpp.

References Nektar::SpatialDomains::Geometry::GetCoordim().

Referenced by v_ContainsPoint().

791  {
792  Array<OneD,NekDouble> locCoord(GetCoordim(),0.0);
793  return v_ContainsPoint(gloCoord,locCoord,tol);
794 
795  }
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Definition: QuadGeom.cpp:789
bool Nektar::SpatialDomains::QuadGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 799 of file QuadGeom.cpp.

References v_ContainsPoint().

802  {
803  NekDouble resid;
804  return v_ContainsPoint(gloCoord,stdCoord,tol,resid);
805  }
virtual bool v_ContainsPoint(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0)
Definition: QuadGeom.cpp:789
double NekDouble
bool Nektar::SpatialDomains::QuadGeom::v_ContainsPoint ( const Array< OneD, const NekDouble > &  gloCoord,
Array< OneD, NekDouble > &  locCoord,
NekDouble  tol,
NekDouble resid 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 807 of file QuadGeom.cpp.

References ASSERTL1, and Nektar::SpatialDomains::Geometry::GetLocCoords().

811  {
812  ASSERTL1(gloCoord.num_elements() >= 2,
813  "Two dimensional geometry expects at least two coordinates.");
814 
815  resid = GetLocCoords(gloCoord, stdCoord);
816  if (stdCoord[0] >= -(1+tol) && stdCoord[1] >= -(1+tol)
817  && stdCoord[0] <= (1+tol) && stdCoord[1] <= (1+tol))
818  {
819  return true;
820  }
821  return false;
822  }
NekDouble GetLocCoords(const Array< OneD, const NekDouble > &coords, Array< OneD, NekDouble > &Lcoords)
Definition: Geometry.h:450
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
void Nektar::SpatialDomains::QuadGeom::v_FillGeom ( )
protectedvirtual

Put all quadrature information into edge structure.

Note verts and edges are listed according to anticlockwise convention but points in _coeffs have to be in array format from left to right.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 528 of file QuadGeom.cpp.

References ASSERTL0, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::Geometry::GetXmap(), Nektar::LibUtilities::Interp2D(), kNedges, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, m_curve, m_edges, m_eorient, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_state, Nektar::SpatialDomains::Geometry::m_xmap, and npts.

Referenced by v_GenGeomFactors(), and v_GetLocCoords().

529  {
530  // check to see if geometry structure is already filled
531  if(m_state != ePtsFilled)
532  {
533  int i,j,k;
534  int nEdgeCoeffs;
535 
536  if (m_curve)
537  {
538  int npts = m_curve->m_points.size();
539  int nEdgePts = (int)sqrt(static_cast<NekDouble>(npts));
540  Array<OneD, NekDouble> tmp(npts);
541  Array<OneD, NekDouble> tmp2(m_xmap->GetTotPoints());
542  LibUtilities::PointsKey curveKey(
543  nEdgePts, m_curve->m_ptype);
544 
545  // Sanity checks:
546  // - Curved faces should have square number of points;
547  // - Each edge should have sqrt(npts) points.
548  ASSERTL0(nEdgePts * nEdgePts == npts,
549  "NUMPOINTS should be a square number in"
550  " quadrilteral "
551  + boost::lexical_cast<string>(m_globalID));
552 
553  for (i = 0; i < kNedges; ++i)
554  {
555  ASSERTL0(
556  m_edges[i]->GetXmap()->GetNcoeffs() == nEdgePts,
557  "Number of edge points does not correspond to "
558  "number of face points in quadrilateral "
559  + boost::lexical_cast<string>(m_globalID));
560  }
561 
562  for (i = 0; i < m_coordim; ++i)
563  {
564  for (j = 0; j < npts; ++j)
565  {
566  tmp[j] = (m_curve->m_points[j]->GetPtr())[i];
567  }
568 
569  // Interpolate m_curve points to GLL points
571  curveKey, curveKey, tmp,
572  m_xmap->GetBasis(0)->GetPointsKey(),
573  m_xmap->GetBasis(1)->GetPointsKey(),
574  tmp2);
575 
576  // Forwards transform to get coefficient space.
577  m_xmap->FwdTrans(tmp2, m_coeffs[i]);
578  }
579  }
580 
581  // Now fill in edges.
582  Array<OneD, unsigned int> mapArray;
583  Array<OneD, int> signArray;
584 
585  for(i = 0; i < kNedges; i++)
586  {
587  m_edges[i]->FillGeom();
588  m_xmap->GetEdgeToElementMap(i,m_eorient[i],
589  mapArray,signArray);
590 
591  nEdgeCoeffs = m_edges[i]->GetXmap()->GetNcoeffs();
592 
593  for(j = 0; j < m_coordim; j++)
594  {
595  for(k = 0; k < nEdgeCoeffs; k++)
596  {
597  m_coeffs[j][mapArray[k]]
598  = signArray[k]*(m_edges[i]->GetCoeffs(j))[k];
599  }
600  }
601  }
602 
604  }
605  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
StdRegions::StdExpansionSharedPtr GetXmap() const
Definition: Geometry.h:383
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...
Definition: Interp.cpp:116
static std::string npts
Definition: InputFld.cpp:43
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometric information has been generated.
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
GeomState m_state
enum identifier to determine if quad points are filled
Definition: Geometry.h:175
int m_coordim
coordinate dimension
Definition: Geometry.h:169
void Nektar::SpatialDomains::QuadGeom::v_GenGeomFactors ( )
protectedvirtual

Set up GeoFac for this geometry using Coord quadrature distribution

Implements Nektar::SpatialDomains::Geometry.

Definition at line 446 of file QuadGeom.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), Nektar::SpatialDomains::eDeformed, Nektar::SpatialDomains::ePtsFilled, Nektar::SpatialDomains::eRegular, Nektar::NekConstants::kNekZeroTol, Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, Nektar::SpatialDomains::Geometry::m_geomFactors, Nektar::SpatialDomains::Geometry::m_geomFactorsState, m_verts, Nektar::SpatialDomains::Geometry::m_xmap, and v_FillGeom().

447  {
449  {
450  int i;
451  GeomType Gtype = eRegular;
452 
454 
455  // We will first check whether we have a regular or deformed
456  // geometry. We will define regular as those cases where the
457  // Jacobian and the metric terms of the derivative are constants
458  // (i.e. not coordinate dependent)
459 
460  // Check to see if expansions are linear
461  // If not linear => deformed geometry
462  for(i = 0; i < m_coordim; ++i)
463  {
464  if((m_xmap->GetBasisNumModes(0) != 2)||
465  (m_xmap->GetBasisNumModes(1) != 2))
466  {
467  Gtype = eDeformed;
468  }
469  }
470 
471  // For linear expansions, the mapping from standard to local
472  // element is given by the relation:
473  // x_i = 0.25 * [ ( x_i^A + x_i^B + x_i^C + x_i^D) +
474  // (-x_i^A + x_i^B + x_i^C - x_i^D)*xi_1 +
475  // (-x_i^A - x_i^B + x_i^C + x_i^D)*xi_2 +
476  // ( x_i^A - x_i^B + x_i^C - x_i^D)*xi_1*xi_2 ]
477  //
478  // The jacobian of the transformation and the metric terms
479  // dxi_i/dx_j, involve only terms of the form dx_i/dxi_j (both
480  // for coordim == 2 or 3). Inspecting the formula above, it can
481  // be appreciated that the derivatives dx_i/dxi_j will be
482  // constant, if the coefficient of the non-linear term is zero.
483  //
484  // That is why for regular geometry, we require
485  //
486  // x_i^A - x_i^B + x_i^C - x_i^D = 0
487  //
488  // or equivalently
489  //
490  // x_i^A - x_i^B = x_i^D - x_i^C
491  //
492  // This corresponds to quadrilaterals which are paralellograms.
493  if(Gtype == eRegular)
494  {
495  for(i = 0; i < m_coordim; i++)
496  {
497  if( fabs( (*m_verts[0])(i) - (*m_verts[1])(i) +
498  (*m_verts[2])(i) - (*m_verts[3])(i) )
500  {
501  Gtype = eDeformed;
502  break;
503  }
504  }
505  }
506 
508  Gtype, m_coordim, m_xmap, m_coeffs);
510  }
511  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
GeomFactorsSharedPtr m_geomFactors
Definition: Geometry.h:170
static const NekDouble kNekZeroTol
virtual void v_FillGeom()
Put all quadrature information into edge structure.
Definition: QuadGeom.cpp:528
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometry is straight-sided with constant geometric factors.
Geometric information has been generated.
GeomType
Indicates the type of element geometry.
Geometry is curved or has non-constant factors.
int m_coordim
coordinate dimension
Definition: Geometry.h:169
const LibUtilities::BasisSharedPtr Nektar::SpatialDomains::QuadGeom::v_GetBasis ( const int  i)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 411 of file QuadGeom.cpp.

References Nektar::SpatialDomains::Geometry::m_xmap.

412  {
413  return m_xmap->GetBasis(i);
414  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::v_GetCartesianEorient ( const int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 727 of file QuadGeom.cpp.

References ASSERTL2, Nektar::StdRegions::eBackwards, Nektar::StdRegions::eForwards, and m_eorient.

728  {
729  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
730  if(i < 2)
731  {
732  return m_eorient[i];
733  }
734  else
735  {
737  {
738  return StdRegions::eBackwards;
739  }
740  else
741  {
742  return StdRegions::eForwards;
743  }
744  }
745  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
NekDouble Nektar::SpatialDomains::QuadGeom::v_GetCoord ( const int  i,
const Array< OneD, const NekDouble > &  Lcoord 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 437 of file QuadGeom.cpp.

References GetCoord().

438  {
439  return GetCoord(i,Lcoord);
440  }
NekDouble GetCoord(const int i, const Array< OneD, const NekDouble > &Lcoord)
Definition: QuadGeom.cpp:242
int Nektar::SpatialDomains::QuadGeom::v_GetCoordim ( void  ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 402 of file QuadGeom.cpp.

References Nektar::SpatialDomains::Geometry::m_coordim.

403  {
404  return m_coordim;
405  }
int m_coordim
coordinate dimension
Definition: Geometry.h:169
const Geometry1DSharedPtr Nektar::SpatialDomains::QuadGeom::v_GetEdge ( int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 707 of file QuadGeom.cpp.

References ASSERTL2, and m_edges.

708  {
709  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
710  return m_edges[i];
711  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
const LibUtilities::BasisSharedPtr Nektar::SpatialDomains::QuadGeom::v_GetEdgeBasis ( const int  i)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 420 of file QuadGeom.cpp.

References ASSERTL1, and Nektar::SpatialDomains::Geometry::m_xmap.

421  {
422  ASSERTL1(i <= 3,"edge is out of range");
423 
424  if (i == 0 || i == 2)
425  {
426  return m_xmap->GetBasis(0);
427  }
428  else
429  {
430  return m_xmap->GetBasis(1);
431  }
432  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:218
int Nektar::SpatialDomains::QuadGeom::v_GetEid ( int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 677 of file QuadGeom.cpp.

References ASSERTL2, and m_edges.

678  {
679  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
680  return m_edges[i]->GetEid();
681  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::v_GetEorient ( const int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 717 of file QuadGeom.cpp.

References ASSERTL2, and m_eorient.

718  {
719  ASSERTL2((i >=0) && (i <= 3),"Edge id must be between 0 and 3");
720  return m_eorient[i];
721  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
StdRegions::Orientation m_eorient[kNedges]
Definition: QuadGeom.h:107
int Nektar::SpatialDomains::QuadGeom::v_GetFid ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 393 of file QuadGeom.cpp.

References m_fid.

394  {
395  return m_fid;
396  }
NekDouble Nektar::SpatialDomains::QuadGeom::v_GetLocCoords ( const Array< OneD, const NekDouble > &  coords,
Array< OneD, NekDouble > &  Lcoords 
)
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 610 of file QuadGeom.cpp.

References Nektar::SpatialDomains::PointGeom::dot(), Nektar::SpatialDomains::eRegular, Nektar::SpatialDomains::Geometry::GetMetricInfo(), Vmath::Imin(), Nektar::SpatialDomains::Geometry::m_coeffs, Nektar::SpatialDomains::Geometry::m_coordim, m_verts, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::PointGeom::Mult(), Nektar::SpatialDomains::Geometry2D::NewtonIterationForLocCoord(), npts, Vmath::Sadd(), Nektar::SpatialDomains::PointGeom::Sub(), v_FillGeom(), Vmath::Vmul(), and Vmath::Vvtvp().

612  {
613  NekDouble resid = 0.0;
614  if(GetMetricInfo()->GetGtype() == eRegular)
615  {
616  NekDouble coords2 = (m_coordim == 3)? coords[2]: 0.0;
617  PointGeom dv1, dv2, norm, orth1, orth2;
618  PointGeom xin(m_coordim,0,coords[0],coords[1],coords2);
619 
620  // Calculate edge vectors from 0-1 and 0-3 edges.
621  dv1.Sub(*m_verts[1],*m_verts[0]);
622  dv2.Sub(*m_verts[3],*m_verts[0]);
623 
624  // Obtain normal to plane in which dv1 and dv2 lie
625  norm.Mult(dv1,dv2);
626 
627  // Obtain vector which are normal to dv1 and dv2.
628  orth1.Mult(norm,dv1);
629  orth2.Mult(norm,dv2);
630 
631  // Start with vector of desired points minus vertex_0
632  xin -= *m_verts[0];
633 
634  // Calculate length using L/|dv1| = (x-v0).n1/(dv1.n1) for coordiante 1
635  // Then rescale to [-1,1].
636  Lcoords[0] = xin.dot(orth2)/dv1.dot(orth2);
637  Lcoords[0] = 2*Lcoords[0]-1;
638  Lcoords[1] = xin.dot(orth1)/dv2.dot(orth1);
639  Lcoords[1] = 2*Lcoords[1]-1;
640  }
641  else
642  {
644 
645  // Determine nearest point of coords to values in m_xmap
646  int npts = m_xmap->GetTotPoints();
647  Array<OneD, NekDouble> ptsx(npts), ptsy(npts);
648  Array<OneD, NekDouble> tmpx(npts), tmpy(npts);
649 
650  m_xmap->BwdTrans(m_coeffs[0], ptsx);
651  m_xmap->BwdTrans(m_coeffs[1], ptsy);
652 
653  const Array<OneD, const NekDouble> za = m_xmap->GetPoints(0);
654  const Array<OneD, const NekDouble> zb = m_xmap->GetPoints(1);
655 
656  //guess the first local coords based on nearest point
657  Vmath::Sadd(npts, -coords[0], ptsx,1,tmpx,1);
658  Vmath::Sadd(npts, -coords[1], ptsy,1,tmpy,1);
659  Vmath::Vmul (npts, tmpx,1,tmpx,1,tmpx,1);
660  Vmath::Vvtvp(npts, tmpy,1,tmpy,1,tmpx,1,tmpx,1);
661 
662  int min_i = Vmath::Imin(npts,tmpx,1);
663 
664  Lcoords[0] = za[min_i%za.num_elements()];
665  Lcoords[1] = zb[min_i/za.num_elements()];
666 
667  // Perform newton iteration to find local coordinates
668  NewtonIterationForLocCoord(coords, ptsx, ptsy, Lcoords,resid);
669  }
670  return resid;
671  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
int Imin(int n, const T *x, const int incx)
Return the index of the minimum element in x.
Definition: Vmath.cpp:833
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
Definition: Vmath.cpp:428
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)
Definition: Geometry2D.cpp:103
static std::string npts
Definition: InputFld.cpp:43
double NekDouble
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.cpp:301
virtual void v_FillGeom()
Put all quadrature information into edge structure.
Definition: QuadGeom.cpp:528
Array< OneD, Array< OneD, NekDouble > > m_coeffs
Definition: Geometry.h:180
Geometry is straight-sided with constant geometric factors.
GeomFactorsSharedPtr GetMetricInfo()
Definition: Geometry.h:299
int m_coordim
coordinate dimension
Definition: Geometry.h:169
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.
Definition: Vmath.cpp:169
int Nektar::SpatialDomains::QuadGeom::v_GetNumEdges ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 783 of file QuadGeom.cpp.

References kNedges.

784  {
785  return kNedges;
786  }
int Nektar::SpatialDomains::QuadGeom::v_GetNumVerts ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 774 of file QuadGeom.cpp.

References kNverts.

775  {
776  return kNverts;
777  }
PointGeomSharedPtr Nektar::SpatialDomains::QuadGeom::v_GetVertex ( int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 697 of file QuadGeom.cpp.

References ASSERTL2, and m_verts.

698  {
699  ASSERTL2((i >=0) && (i <= 3),"Vertex id must be between 0 and 3");
700  return m_verts[i];
701  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
int Nektar::SpatialDomains::QuadGeom::v_GetVid ( int  i) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 687 of file QuadGeom.cpp.

References ASSERTL2, and m_verts.

688  {
689  ASSERTL2((i >=0) && (i <= 3),"Verted id must be between 0 and 3");
690  return m_verts[i]->GetVid();
691  }
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:240
bool Nektar::SpatialDomains::QuadGeom::v_IsElmtConnected ( int  gvo_id,
int  locid 
) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 373 of file QuadGeom.cpp.

References Nektar::StdRegions::find(), and m_elmtMap.

374  {
375  std::list<CompToElmt>::const_iterator def;
376  CompToElmt ee(gvo_id,locid);
377 
378  def = find(m_elmtMap.begin(),m_elmtMap.end(),ee);
379 
380  // Found the element connectivity object in the list
381  if(def != m_elmtMap.end())
382  {
383  return(true);
384  }
385 
386  return(false);
387  }
std::list< CompToElmt > m_elmtMap
Definition: QuadGeom.h:110
InputIterator find(InputIterator first, InputIterator last, InputIterator startingpoint, const EqualityComparable &value)
Definition: StdRegions.hpp:315
int Nektar::SpatialDomains::QuadGeom::v_NumElmtConnected ( ) const
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 364 of file QuadGeom.cpp.

References m_elmtMap.

365  {
366  return int(m_elmtMap.size());
367  }
std::list< CompToElmt > m_elmtMap
Definition: QuadGeom.h:110
void Nektar::SpatialDomains::QuadGeom::v_Reset ( CurveMap curvedEdges,
CurveMap curvedFaces 
)
protectedvirtual

Reset this geometry object: unset the current state and remove allocated GeomFactors.

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 824 of file QuadGeom.cpp.

References Nektar::iterator, m_curve, m_edges, Nektar::SpatialDomains::Geometry::m_globalID, Nektar::SpatialDomains::Geometry::m_xmap, Nektar::SpatialDomains::Geometry::SetUpCoeffs(), SetUpXmap(), and Nektar::SpatialDomains::Geometry::v_Reset().

826  {
827  Geometry::v_Reset(curvedEdges, curvedFaces);
828  CurveMap::iterator it = curvedFaces.find(m_globalID);
829 
830  if (it != curvedFaces.end())
831  {
832  m_curve = it->second;
833  }
834 
835  for (int i = 0; i < 4; ++i)
836  {
837  m_edges[i]->Reset(curvedEdges, curvedFaces);
838  }
839 
840  SetUpXmap();
841  SetUpCoeffs(m_xmap->GetNcoeffs());
842  }
StdRegions::StdExpansionSharedPtr m_xmap
Definition: Geometry.h:172
virtual void v_Reset(CurveMap &curvedEdges, CurveMap &curvedFaces)
Reset this geometry object: unset the current state and remove allocated GeomFactors.
Definition: Geometry.cpp:307
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
void SetUpCoeffs(const int nCoeffs)
Initialise the m_coeffs array.
Definition: Geometry.h:484
void Nektar::SpatialDomains::QuadGeom::v_SetOwnData ( )
protectedvirtual

Reimplemented from Nektar::SpatialDomains::Geometry.

Definition at line 517 of file QuadGeom.cpp.

References m_ownData.

518  {
519  m_ownData = true;
520  }
bool m_ownData
Boolean indicating whether object owns the data.
Definition: QuadGeom.h:193
int Nektar::SpatialDomains::QuadGeom::v_WhichEdge ( SegGeomSharedPtr  edge)
protectedvirtual

Return the edge number of the given edge.

Reimplemented from Nektar::SpatialDomains::Geometry2D.

Definition at line 751 of file QuadGeom.cpp.

References Nektar::iterator, and m_edges.

752  {
753  int returnval = -1;
754 
755  SegGeomVector::iterator edgeIter;
756  int i;
757 
758  for (i=0,edgeIter = m_edges.begin(); edgeIter != m_edges.end(); ++edgeIter,++i)
759  {
760  if (*edgeIter == edge)
761  {
762  returnval = i;
763  break;
764  }
765  }
766 
767  return returnval;
768  }
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator

Member Data Documentation

const int Nektar::SpatialDomains::QuadGeom::kNedges = 4
static
const int Nektar::SpatialDomains::QuadGeom::kNverts = 4
static
CurveSharedPtr Nektar::SpatialDomains::QuadGeom::m_curve
protected

Definition at line 111 of file QuadGeom.h.

Referenced by v_FillGeom(), and v_Reset().

SegGeomVector Nektar::SpatialDomains::QuadGeom::m_edges
protected

Definition at line 106 of file QuadGeom.h.

Referenced by QuadGeom(), SetUpXmap(), v_FillGeom(), v_GetEdge(), v_GetEid(), v_Reset(), and v_WhichEdge().

std::list<CompToElmt> Nektar::SpatialDomains::QuadGeom::m_elmtMap
protected

Definition at line 110 of file QuadGeom.h.

Referenced by QuadGeom(), v_AddElmtConnected(), v_IsElmtConnected(), and v_NumElmtConnected().

StdRegions::Orientation Nektar::SpatialDomains::QuadGeom::m_eorient[kNedges]
protected

Definition at line 107 of file QuadGeom.h.

Referenced by QuadGeom(), v_FillGeom(), v_GetCartesianEorient(), and v_GetEorient().

int Nektar::SpatialDomains::QuadGeom::m_fid
protected

Definition at line 108 of file QuadGeom.h.

Referenced by QuadGeom(), and v_GetFid().

bool Nektar::SpatialDomains::QuadGeom::m_ownData
private

Boolean indicating whether object owns the data.

Definition at line 193 of file QuadGeom.h.

Referenced by QuadGeom(), and v_SetOwnData().

bool Nektar::SpatialDomains::QuadGeom::m_ownVerts
protected

Definition at line 109 of file QuadGeom.h.

Referenced by QuadGeom().

PointGeomVector Nektar::SpatialDomains::QuadGeom::m_verts
protected
const std::string Nektar::SpatialDomains::QuadGeom::XMLElementType
static

Definition at line 102 of file QuadGeom.h.