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

#include <NodalTetElec.h>

Inheritance diagram for Nektar::LibUtilities::NodalTetElec:
Inheritance graph
[legend]
Collaboration diagram for Nektar::LibUtilities::NodalTetElec:
Collaboration graph
[legend]

Public Member Functions

virtual ~NodalTetElec ()
 
const MatrixSharedPtrType GetI (const PointsKey &pkey)
 
const MatrixSharedPtrType GetI (const Array< OneD, const NekDouble > &x, const Array< OneD, const NekDouble > &y, const Array< OneD, const NekDouble > &z)
 
 NodalTetElec (const PointsKey &key)
 
- Public Member Functions inherited from Nektar::LibUtilities::Points< NekDouble >
virtual ~Points ()
 
virtual void Initialize (void)
 
unsigned int GetPointsDim () const
 
unsigned int GetNumPoints () const
 
unsigned int GetTotNumPoints () const
 
PointsType GetPointsType () const
 
const Array< OneD, const
DataType > & 
GetZ () const
 
const Array< OneD, const
DataType > & 
GetW () const
 
void GetZW (Array< OneD, const DataType > &z, Array< OneD, const DataType > &w) const
 
void GetPoints (Array< OneD, const DataType > &x) const
 
void GetPoints (Array< OneD, const DataType > &x, Array< OneD, const DataType > &y) const
 
void GetPoints (Array< OneD, const DataType > &x, Array< OneD, const DataType > &y, Array< OneD, const DataType > &z) const
 
const MatrixSharedPtrTypeGetD (Direction dir=xDir) const
 
virtual const MatrixSharedPtrType GetI (const Array< OneD, const DataType > &x)
 
virtual const MatrixSharedPtrType GetI (unsigned int numpoints, const Array< OneD, const DataType > &x)
 
virtual const MatrixSharedPtrType GetI (const Array< OneD, const DataType > &x, const Array< OneD, const DataType > &y)
 
virtual const MatrixSharedPtrType GetGalerkinProjection (const PointsKey &pkey)
 

Static Public Member Functions

static boost::shared_ptr
< PointsBaseType
Create (const PointsKey &key)
 

Private Member Functions

 NodalTetElec ()
 
void CalculatePoints ()
 
void CalculateWeights ()
 
void CalculateDerivMatrix ()
 
void NodalPointReorder3d ()
 
void CalculateInterpMatrix (const Array< OneD, const NekDouble > &xia, const Array< OneD, const NekDouble > &yia, const Array< OneD, const NekDouble > &zia, Array< OneD, NekDouble > &interp)
 

Private Attributes

boost::shared_ptr
< NodalUtilTetrahedron
m_util
 

Additional Inherited Members

- Public Types inherited from Nektar::LibUtilities::Points< NekDouble >
typedef NekDouble DataType
 
typedef boost::shared_ptr
< NekMatrix< DataType > > 
MatrixSharedPtrType
 
- Protected Member Functions inherited from Nektar::LibUtilities::Points< NekDouble >
 Points (const PointsKey &key)
 
- Protected Attributes inherited from Nektar::LibUtilities::Points< NekDouble >
PointsKey m_pointsKey
 
Array< OneD, DataTypem_points [3]
 
Array< OneD, DataTypem_weights
 
MatrixSharedPtrType m_derivmatrix [3]
 
NekManager< PointsKey,
NekMatrix< DataType >
, PointsKey::opLess
m_InterpManager
 
NekManager< PointsKey,
NekMatrix< DataType >
, PointsKey::opLess
m_GalerkinProjectionManager
 

Detailed Description

Definition at line 52 of file NodalTetElec.h.

Constructor & Destructor Documentation

virtual Nektar::LibUtilities::NodalTetElec::~NodalTetElec ( )
inlinevirtual

Definition at line 55 of file NodalTetElec.h.

56  {
57  }
Nektar::LibUtilities::NodalTetElec::NodalTetElec ( const PointsKey key)
inline

Definition at line 88 of file NodalTetElec.h.

88  : PointsBaseType(key)
89  {
90 
91  }
Points< NekDouble > PointsBaseType
Nektar::LibUtilities::NodalTetElec::NodalTetElec ( )
inlineprivate

Definition at line 96 of file NodalTetElec.h.

97  {
98 
99  }
static const PointsKey NullPointsKey(0, eNoPointsType)
Points< NekDouble > PointsBaseType

Member Function Documentation

void Nektar::LibUtilities::NodalTetElec::CalculateDerivMatrix ( )
privatevirtual

Reimplemented from Nektar::LibUtilities::Points< NekDouble >.

Definition at line 230 of file NodalTetElec.cpp.

References Nektar::LibUtilities::Points< NekDouble >::CalculateDerivMatrix(), Nektar::LibUtilities::Points< NekDouble >::m_derivmatrix, and m_util.

231  {
232  // Allocate the derivative matrix.
234 
235  m_derivmatrix[0] = m_util->GetDerivMatrix(0);
236  m_derivmatrix[1] = m_util->GetDerivMatrix(1);
237  m_derivmatrix[2] = m_util->GetDerivMatrix(2);
238  }
MatrixSharedPtrType m_derivmatrix[3]
Definition: Points.h:373
boost::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:94
void Nektar::LibUtilities::NodalTetElec::CalculateInterpMatrix ( const Array< OneD, const NekDouble > &  xia,
const Array< OneD, const NekDouble > &  yia,
const Array< OneD, const NekDouble > &  zia,
Array< OneD, NekDouble > &  interp 
)
private

Definition at line 215 of file NodalTetElec.cpp.

References m_util, and Vmath::Vcopy().

Referenced by GetI().

218  {
219  Array<OneD, Array<OneD, NekDouble> > xi(3);
220  xi[0] = xia;
221  xi[1] = yia;
222  xi[2] = zia;
223 
224  boost::shared_ptr<NekMatrix<NekDouble> > mat =
225  m_util->GetInterpolationMatrix(xi);
226  Vmath::Vcopy(mat->GetRows() * mat->GetColumns(), mat->GetRawPtr(),
227  1, &interp[0], 1);
228  }
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
boost::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:94
void Nektar::LibUtilities::NodalTetElec::CalculatePoints ( )
privatevirtual

Reimplemented from Nektar::LibUtilities::Points< NekDouble >.

Definition at line 51 of file NodalTetElec.cpp.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), ASSERTL1, Nektar::LibUtilities::Points< DataT >::CalculatePoints(), Nektar::LibUtilities::Points< NekDouble >::GetNumPoints(), Nektar::LibUtilities::PointsKey::GetTotNumPoints(), Nektar::LibUtilities::Points< NekDouble >::m_points, Nektar::LibUtilities::Points< NekDouble >::m_pointsKey, m_util, NodalPointReorder3d(), Nektar::LibUtilities::NodalTetElecData, Nektar::LibUtilities::NodalTetElecNPTS, Nektar::LibUtilities::perm12A_3d, Nektar::LibUtilities::perm12B_3d, Nektar::LibUtilities::perm12C_3d, Nektar::LibUtilities::perm24_3d, Nektar::LibUtilities::perm4_3d, and Nektar::LibUtilities::perm6_3d.

52  {
53  // Allocate the storage for points
55 
56  int index=0,isum=0;
57  const int offset = 5; //offset to match Datafile
58  NekDouble b,c,d;
59  unsigned int numPoints = GetNumPoints();
60 
61  // initialize values
62  for(unsigned int i=0; i < numPoints-2; ++i)
63  {
64  index += NodalTetElecNPTS[i];
65  }
66 
67  for(unsigned int i=0; i < NodalTetElecNPTS[numPoints-2]; ++i, ++index)
68  {
69  // 1 Point Symmetry: aaaa
70  if(int(NodalTetElecData[index][0]))
71  {
72  b = NodalTetElecData[index][6];
73  c = NodalTetElecData[index][7];
74  d = NodalTetElecData[index][8];
75 
76  m_points[0][isum] = 2.0*b - 1.0;
77  m_points[1][isum] = 2.0*c - 1.0;
78  m_points[2][isum] = 2.0*d - 1.0;
79  isum++;
80  continue;
81  }//end symmetry 1
82 
83 
84  // 4 Point symmetry: aaab or abbb
85  if(int(NodalTetElecData[index][1]))
86  {
87  for(unsigned int j=0; j < 4; ++j)
88  {
89  b = NodalTetElecData[index][offset + perm4_3d[j][1]];
90  c = NodalTetElecData[index][offset + perm4_3d[j][2]];
91  d = NodalTetElecData[index][offset + perm4_3d[j][3]];
92 
93  m_points[0][isum] = 2.0*b - 1.0;
94  m_points[1][isum] = 2.0*c - 1.0;
95  m_points[2][isum] = 2.0*d - 1.0;
96  isum++;
97  }//end j
98  continue;
99  }//end symmetry 4
100 
101 
102  // 6 Point symmetry: aabb
103  if(int(NodalTetElecData[index][2]))
104  {
105  for(unsigned int j=0; j < 6; ++j)
106  {
107  b = NodalTetElecData[index][offset + perm6_3d[j][1]];
108  c = NodalTetElecData[index][offset + perm6_3d[j][2]];
109  d = NodalTetElecData[index][offset + perm6_3d[j][3]];
110 
111  m_points[0][isum] = 2.0*b - 1.0;
112  m_points[1][isum] = 2.0*c - 1.0;
113  m_points[2][isum] = 2.0*d - 1.0;
114  isum++;
115  }//end j
116  continue;
117  }//end symmetry6
118 
119 
120  // 12 Point symmetry: case aabc
121  if(int(NodalTetElecData[index][3]) == 1)
122  {
123  for(unsigned int j=0; j < 12; ++j)
124  {
125  b = NodalTetElecData[index][offset + perm12A_3d[j][1]];
126  c = NodalTetElecData[index][offset + perm12A_3d[j][2]];
127  d = NodalTetElecData[index][offset + perm12A_3d[j][3]];
128 
129  m_points[0][isum] = 2.0*b - 1.0;
130  m_points[1][isum] = 2.0*c - 1.0;
131  m_points[2][isum] = 2.0*d - 1.0;
132  isum++;
133  }//end j
134  continue;
135  }//end symmetry 12 aabc
136 
137 
138  // 12 Point symmetry: case abcc
139  if(int(NodalTetElecData[index][3]) == 2)
140  {
141  for(unsigned int j=0; j < 12; ++j)
142  {
143  b = NodalTetElecData[index][offset + perm12B_3d[j][1]];
144  c = NodalTetElecData[index][offset + perm12B_3d[j][2]];
145  d = NodalTetElecData[index][offset + perm12B_3d[j][3]];
146 
147  m_points[0][isum] = 2.0*b - 1.0;
148  m_points[1][isum] = 2.0*c - 1.0;
149  m_points[2][isum] = 2.0*d - 1.0;
150  isum++;
151  }//end j
152  continue;
153  }//end symmetry 12 abcc
154 
155 
156  // 12 Point symmetry: case abbc
157  if(int(NodalTetElecData[index][3]) == 3)
158  {
159  for(unsigned int j=0; j < 12; ++j)
160  {
161  b = NodalTetElecData[index][offset + perm12C_3d[j][1]];
162  c = NodalTetElecData[index][offset + perm12C_3d[j][2]];
163  d = NodalTetElecData[index][offset + perm12C_3d[j][3]];
164 
165  m_points[0][isum] = 2.0*b - 1.0;
166  m_points[1][isum] = 2.0*c - 1.0;
167  m_points[2][isum] = 2.0*d - 1.0;
168  isum++;
169  }//end j
170  continue;
171  }//end symmetry 12 abbc
172 
173  // 24 Point symmetry: case abcd
174  if(int(NodalTetElecData[index][4]))
175  {
176  for(unsigned int j=0; j < 24; ++j)
177  {
178  b = NodalTetElecData[index][offset + perm24_3d[j][1]];
179  c = NodalTetElecData[index][offset + perm24_3d[j][2]];
180  d = NodalTetElecData[index][offset + perm24_3d[j][3]];
181 
182  m_points[0][isum] = 2.0*b - 1.0;
183  m_points[1][isum] = 2.0*c - 1.0;
184  m_points[2][isum] = 2.0*d - 1.0;
185  isum++;
186  }//end j
187  continue;
188  }//end symmetry24abcd
189 
190 
191  }//end npts
192 
194 
195  ASSERTL1((static_cast<unsigned int>(isum)==m_pointsKey.GetTotNumPoints()),"sum not equal to npts");
196 
198  numPoints - 1, m_points[0], m_points[1], m_points[2]);
199  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
virtual void CalculatePoints()
Definition: Points.h:377
Array< OneD, DataType > m_points[3]
Definition: Points.h:371
static const unsigned int perm12C_3d[12][4]
static const unsigned int perm6_3d[6][4]
unsigned int GetNumPoints() const
Definition: Points.h:268
unsigned int GetTotNumPoints() const
Definition: Points.h:179
double NekDouble
static const unsigned int perm4_3d[4][4]
static const unsigned int NodalTetElecNPTS[NodalTetElecAvailable]
static const unsigned int perm24_3d[24][4]
static const NekDouble NodalTetElecData[][9]
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
static const unsigned int perm12A_3d[12][4]
boost::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:94
static const unsigned int perm12B_3d[12][4]
void Nektar::LibUtilities::NodalTetElec::CalculateWeights ( )
privatevirtual

Reimplemented from Nektar::LibUtilities::Points< NekDouble >.

Definition at line 201 of file NodalTetElec.cpp.

References Nektar::LibUtilities::Points< NekDouble >::CalculateWeights(), Nektar::NekVector< DataType >::GetPtr(), Nektar::NekVector< DataType >::GetRows(), m_util, and Nektar::LibUtilities::Points< NekDouble >::m_weights.

202  {
203  // Allocate the storage for points
205 
206  typedef DataType T;
207 
208  // Solve the Vandermonde system of integrals for the weight vector
209  NekVector<T> w = m_util->GetWeights();
210  m_weights = Array<OneD,T>( w.GetRows(), w.GetPtr() );
211  }
Array< OneD, DataType > m_weights
Definition: Points.h:372
boost::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:94
boost::shared_ptr< PointsBaseType > Nektar::LibUtilities::NodalTetElec::Create ( const PointsKey key)
static

Definition at line 240 of file NodalTetElec.cpp.

241  {
242  boost::shared_ptr<PointsBaseType> returnval(MemoryManager<NodalTetElec>::AllocateSharedPtr(key));
243  returnval->Initialize();
244  return returnval;
245  }
const MatrixSharedPtrType Nektar::LibUtilities::NodalTetElec::GetI ( const PointsKey pkey)
inlinevirtual

Reimplemented from Nektar::LibUtilities::Points< NekDouble >.

Definition at line 62 of file NodalTetElec.h.

References ASSERTL0, Nektar::LibUtilities::PointsKey::GetPointsDim(), and Nektar::LibUtilities::PointsManager().

63  {
64  ASSERTL0(pkey.GetPointsDim() == 3,
65  "NodalTetElec Points can only interp to other 3d "
66  "point distributions");
67  Array<OneD, const NekDouble> x, y, z;
68  PointsManager()[pkey]->GetPoints(x, y, z);
69  return GetI(x, y, z);
70  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
const MatrixSharedPtrType GetI(const PointsKey &pkey)
Definition: NodalTetElec.h:62
PointsManagerT & PointsManager(void)
const MatrixSharedPtrType Nektar::LibUtilities::NodalTetElec::GetI ( const Array< OneD, const NekDouble > &  x,
const Array< OneD, const NekDouble > &  y,
const Array< OneD, const NekDouble > &  z 
)
inlinevirtual

Reimplemented from Nektar::LibUtilities::Points< NekDouble >.

Definition at line 72 of file NodalTetElec.h.

References CalculateInterpMatrix(), and Nektar::LibUtilities::Points< NekDouble >::GetTotNumPoints().

76  {
77  int numpoints = x.num_elements();
78  unsigned int np = GetTotNumPoints();
79 
80  Array<OneD, NekDouble> interp(np*numpoints);
81  CalculateInterpMatrix(x, y, z, interp);
82 
83  NekDouble* d = interp.data();
84  return MemoryManager<NekMatrix<NekDouble> >
85  ::AllocateSharedPtr(numpoints, np, d);
86  }
void CalculateInterpMatrix(const Array< OneD, const NekDouble > &xia, const Array< OneD, const NekDouble > &yia, const Array< OneD, const NekDouble > &zia, Array< OneD, NekDouble > &interp)
double NekDouble
unsigned int GetTotNumPoints() const
Definition: Points.h:273
void Nektar::LibUtilities::NodalTetElec::NodalPointReorder3d ( )
private

Definition at line 247 of file NodalTetElec.cpp.

References Nektar::LibUtilities::Points< NekDouble >::GetNumPoints(), Nektar::NekConstants::kNekZeroTol, and Nektar::LibUtilities::Points< NekDouble >::m_points.

Referenced by CalculatePoints().

248  {
249  int cnt;
250  int istart,iend;
251 
252  const int nVerts = 4;
253  const int nEdgeInteriorPoints = GetNumPoints()-2;
254  const int nFaceInteriorPoints = (GetNumPoints()-3)*(GetNumPoints()-2)/2;
255  //const int nBoundaryPoints = 4 + 6*nEdgeInteriorPoints + 4*nFaceInteriorPoints;
256  const int nAllPoints = GetNumPoints()*(GetNumPoints()+1)*(GetNumPoints()+2)/6;
257  if(nEdgeInteriorPoints==0)
258  {
259  return;
260  }
261 
262  //group all edge 1 points
263  istart = nVerts;
264  for(int i = cnt = istart; i < nAllPoints; i++)
265  {
266  if( fabs(m_points[1][i] + 1.0) < NekConstants::kNekZeroTol &&
267  fabs(m_points[2][i] + 1.0) < NekConstants::kNekZeroTol)
268  {
269  std::swap(m_points[0][cnt], m_points[0][i]);
270  std::swap(m_points[1][cnt], m_points[1][i]);
271  std::swap(m_points[2][cnt], m_points[2][i]);
272  cnt++;
273  }
274  }
275 
276  // bubble sort edge 1 (counterclockwise numbering)
277  iend = istart + nEdgeInteriorPoints;
278  for(int i = istart; i < iend; i++)
279  {
280  for(int j = istart+1; j < iend; j++)
281  {
282  if(m_points[0][j] < m_points[0][j-1])
283  {
284  std::swap(m_points[0][j], m_points[0][j-1]);
285  std::swap(m_points[1][j], m_points[1][j-1]);
286  std::swap(m_points[2][j], m_points[2][j-1]);
287  }
288  }
289  }
290 
291  // group the points of edge 2 together;
292  istart = iend;
293  for(int i = cnt = istart; i < nAllPoints; i++)
294  {
295  if( fabs(m_points[1][i]+m_points[0][i]) < NekConstants::kNekZeroTol &&
296  fabs(m_points[2][i] + 1.0) < NekConstants::kNekZeroTol)
297  {
298  std::swap(m_points[0][cnt], m_points[0][i]);
299  std::swap(m_points[1][cnt], m_points[1][i]);
300  std::swap(m_points[2][cnt], m_points[2][i]);
301  cnt++;
302  }
303  }
304 
305  // bubble sort edge 2 (counterclockwise numbering)
306  iend = istart + nEdgeInteriorPoints;
307  for(int i = istart; i < iend; i++)
308  {
309  for(int j = istart+1; j < iend; j++)
310  {
311  if(m_points[1][j] < m_points[1][j-1])
312  {
313  std::swap(m_points[0][j], m_points[0][j-1]);
314  std::swap(m_points[1][j], m_points[1][j-1]);
315  std::swap(m_points[2][j], m_points[2][j-1]);
316  }
317  }
318  }
319 
320  // group the points of edge 3 together;
321  istart = iend;
322  for(int i = cnt = istart; i < nAllPoints; i++)
323  {
324  if( fabs(m_points[0][i] + 1.0) < NekConstants::kNekZeroTol &&
325  fabs(m_points[2][i] + 1.0) < NekConstants::kNekZeroTol)
326  {
327  std::swap(m_points[0][cnt], m_points[0][i]);
328  std::swap(m_points[1][cnt], m_points[1][i]);
329  std::swap(m_points[2][cnt], m_points[2][i]);
330  cnt++;
331  }
332  }
333 
334  // bubble sort edge 3 (counterclockwise numbering)
335  iend = istart + nEdgeInteriorPoints;
336  for(int i = istart; i < iend; i++)
337  {
338  for(int j = istart+1; j < iend; j++)
339  {
340  if(m_points[1][j] > m_points[1][j-1])
341  {
342  std::swap(m_points[0][j], m_points[0][j-1]);
343  std::swap(m_points[1][j], m_points[1][j-1]);
344  std::swap(m_points[2][j], m_points[2][j-1]);
345  }
346  }
347  }
348 
349  // group the points of edge 4 together;
350  istart = iend;
351  for(int i = cnt = istart; i < nAllPoints; i++)
352  {
353  if( fabs(m_points[0][i] + 1.0) < NekConstants::kNekZeroTol &&
354  fabs(m_points[1][i] + 1.0) < NekConstants::kNekZeroTol)
355  {
356  std::swap(m_points[0][cnt], m_points[0][i]);
357  std::swap(m_points[1][cnt], m_points[1][i]);
358  std::swap(m_points[2][cnt], m_points[2][i]);
359  cnt++;
360  }
361  }
362 
363  // bubble sort edge 3 (counterclockwise numbering)
364  iend = istart + nEdgeInteriorPoints;
365  for(int i = istart; i < iend; i++)
366  {
367  for(int j = istart+1; j < iend; j++)
368  {
369  if(m_points[2][j] < m_points[2][j-1])
370  {
371  std::swap(m_points[0][j], m_points[0][j-1]);
372  std::swap(m_points[1][j], m_points[1][j-1]);
373  std::swap(m_points[2][j], m_points[2][j-1]);
374  }
375  }
376  }
377 
378  // group the points of edge 5 together;
379  istart = iend;
380  for(int i = cnt = istart; i < nAllPoints; i++)
381  {
382  if( fabs(m_points[0][i]+m_points[2][i]) < NekConstants::kNekZeroTol &&
383  fabs(m_points[1][i] + 1.0) < NekConstants::kNekZeroTol)
384  {
385  std::swap(m_points[0][cnt], m_points[0][i]);
386  std::swap(m_points[1][cnt], m_points[1][i]);
387  std::swap(m_points[2][cnt], m_points[2][i]);
388  cnt++;
389  }
390  }
391 
392  // bubble sort edge 5 (counterclockwise numbering)
393  iend = istart + nEdgeInteriorPoints;
394  for(int i = istart; i < iend; i++)
395  {
396  for(int j = istart+1; j < iend; j++)
397  {
398  if(m_points[2][j] < m_points[2][j-1])
399  {
400  std::swap(m_points[0][j], m_points[0][j-1]);
401  std::swap(m_points[1][j], m_points[1][j-1]);
402  std::swap(m_points[2][j], m_points[2][j-1]);
403  }
404  }
405  }
406 
407  // group the points of edge 6 together;
408  istart = iend;
409  for(int i = cnt = istart; i < nAllPoints; i++)
410  {
411  if( fabs(m_points[1][i]+m_points[2][i]) < NekConstants::kNekZeroTol &&
412  fabs(m_points[0][i] + 1.0) < NekConstants::kNekZeroTol)
413  {
414  std::swap(m_points[0][cnt], m_points[0][i]);
415  std::swap(m_points[1][cnt], m_points[1][i]);
416  std::swap(m_points[2][cnt], m_points[2][i]);
417  cnt++;
418  }
419  }
420 
421  // bubble sort edge 6 (counterclockwise numbering)
422  iend = istart + nEdgeInteriorPoints;
423  for(int i = istart; i < iend; i++)
424  {
425  for(int j = istart+1; j < iend; j++)
426  {
427  if(m_points[2][j] < m_points[2][j-1])
428  {
429  std::swap(m_points[0][j], m_points[0][j-1]);
430  std::swap(m_points[1][j], m_points[1][j-1]);
431  std::swap(m_points[2][j], m_points[2][j-1]);
432  }
433  }
434  }
435 
436  if(GetNumPoints() < 4)
437  {
438  //no face points
439  return;
440  }
441 
442  // group the points of face 1 together;
443  istart = iend;
444  for(int i = cnt = istart; i < nAllPoints; i++)
445  {
446  if(fabs(m_points[2][i] + 1.0) < NekConstants::kNekZeroTol)
447  {
448  std::swap(m_points[0][cnt], m_points[0][i]);
449  std::swap(m_points[1][cnt], m_points[1][i]);
450  std::swap(m_points[2][cnt], m_points[2][i]);
451  cnt++;
452  }
453  }
454 
455  // bubble sort face1 (tensor numbering)
456  iend = istart + nFaceInteriorPoints;
457  bool repeat = true;
458  while(repeat)
459  {
460  repeat = false;
461  for(int i = istart; i < iend - 1; i++)
462  {
463  if(m_points[1][i] > m_points[1][i+1])
464  {
465  std::swap(m_points[0][i+1], m_points[0][i]);
466  std::swap(m_points[1][i+1], m_points[1][i]);
467  std::swap(m_points[2][i+1], m_points[2][i]);
468  repeat = true;
469  }
470  }
471  }
472  int offset = 0;
473  int npl = GetNumPoints() - 3;
474  while(npl > 1)
475  {
476  repeat = true;
477  while(repeat)
478  {
479  repeat = false;
480  for(int i = offset+istart; i < offset+istart + npl - 1; i++)
481  {
482  if(m_points[0][i] > m_points[0][i+1])
483  {
484  std::swap(m_points[0][i+1], m_points[0][i]);
485  std::swap(m_points[1][i+1], m_points[1][i]);
486  std::swap(m_points[2][i+1], m_points[2][i]);
487  repeat = true;
488  }
489  }
490  }
491  offset += npl;
492  npl--;
493  }
494 
495  // group the points of face 2 together;
496  istart = iend;
497  for(int i = cnt = istart; i < nAllPoints; i++)
498  {
499  if(fabs(m_points[1][i] + 1.0) < NekConstants::kNekZeroTol)
500  {
501  std::swap(m_points[0][cnt], m_points[0][i]);
502  std::swap(m_points[1][cnt], m_points[1][i]);
503  std::swap(m_points[2][cnt], m_points[2][i]);
504  cnt++;
505  }
506  }
507 
508  // bubble sort face2 (tensor numbering)
509  iend = istart + nFaceInteriorPoints;
510  repeat = true;
511  while(repeat)
512  {
513  repeat = false;
514  for(int i = istart; i < iend - 1; i++)
515  {
516  if(m_points[2][i] > m_points[2][i+1])
517  {
518  std::swap(m_points[0][i+1], m_points[0][i]);
519  std::swap(m_points[1][i+1], m_points[1][i]);
520  std::swap(m_points[2][i+1], m_points[2][i]);
521  repeat = true;
522  }
523  }
524  }
525  offset = 0;
526  npl = GetNumPoints() - 3;
527  while(npl > 1)
528  {
529  repeat = true;
530  while(repeat)
531  {
532  repeat = false;
533  for(int i = offset+istart; i < offset+istart + npl - 1; i++)
534  {
535  if(m_points[0][i] > m_points[0][i+1])
536  {
537  std::swap(m_points[0][i+1], m_points[0][i]);
538  std::swap(m_points[1][i+1], m_points[1][i]);
539  std::swap(m_points[2][i+1], m_points[2][i]);
540  repeat = true;
541  }
542  }
543  }
544  offset += npl;
545  npl--;
546  }
547 
548  // group the points of face 3 together;
549  istart = iend;
550  for(int i = cnt = istart; i < nAllPoints; i++)
551  {
552  if(fabs(m_points[1][i]+m_points[0][i]+m_points[2][i]+1.0) < 1E-9) //nek zero tol too small
553  {
554  std::swap(m_points[0][cnt], m_points[0][i]);
555  std::swap(m_points[1][cnt], m_points[1][i]);
556  std::swap(m_points[2][cnt], m_points[2][i]);
557  cnt++;
558  }
559  }
560 
561  // bubble sort face3 (tensor numbering)
562  iend = istart + nFaceInteriorPoints;
563  repeat = true;
564  while(repeat)
565  {
566  repeat = false;
567  for(int i = istart; i < iend - 1; i++)
568  {
569  if(m_points[2][i] > m_points[2][i+1])
570  {
571  std::swap(m_points[0][i+1], m_points[0][i]);
572  std::swap(m_points[1][i+1], m_points[1][i]);
573  std::swap(m_points[2][i+1], m_points[2][i]);
574  repeat = true;
575  }
576  }
577  }
578  offset = 0;
579  npl = GetNumPoints() - 3;
580  while(npl > 1)
581  {
582  repeat = true;
583  while(repeat)
584  {
585  repeat = false;
586  for(int i = offset+istart; i < offset+istart + npl - 1; i++)
587  {
588  if(m_points[1][i] > m_points[1][i+1])
589  {
590  std::swap(m_points[0][i+1], m_points[0][i]);
591  std::swap(m_points[1][i+1], m_points[1][i]);
592  std::swap(m_points[2][i+1], m_points[2][i]);
593  repeat = true;
594  }
595  }
596  }
597  offset += npl;
598  npl--;
599  }
600 
601  // group the points of face 4 together;
602  istart = iend;
603  for(int i = cnt = istart; i < nAllPoints; i++)
604  {
605  if(fabs(m_points[0][i] + 1.0) < NekConstants::kNekZeroTol)
606  {
607  std::swap(m_points[0][cnt], m_points[0][i]);
608  std::swap(m_points[1][cnt], m_points[1][i]);
609  std::swap(m_points[2][cnt], m_points[2][i]);
610  cnt++;
611  }
612  }
613 
614  // bubble sort face4 (tensor numbering)
615  iend = istart + nFaceInteriorPoints;
616  repeat = true;
617  while(repeat)
618  {
619  repeat = false;
620  for(int i = istart; i < iend - 1; i++)
621  {
622  if(m_points[2][i] > m_points[2][i+1])
623  {
624  std::swap(m_points[0][i+1], m_points[0][i]);
625  std::swap(m_points[1][i+1], m_points[1][i]);
626  std::swap(m_points[2][i+1], m_points[2][i]);
627  repeat = true;
628  }
629  }
630  }
631  offset = 0;
632  npl = GetNumPoints() - 3;
633  while(npl > 1)
634  {
635  repeat = true;
636  while(repeat)
637  {
638  repeat = false;
639  for(int i = offset+istart; i < offset+istart + npl - 1; i++)
640  {
641  if(m_points[1][i] > m_points[1][i+1])
642  {
643  std::swap(m_points[0][i+1], m_points[0][i]);
644  std::swap(m_points[1][i+1], m_points[1][i]);
645  std::swap(m_points[2][i+1], m_points[2][i]);
646  repeat = true;
647  }
648  }
649  }
650  offset += npl;
651  npl--;
652  }
653 
654  }
Array< OneD, DataType > m_points[3]
Definition: Points.h:371
static const NekDouble kNekZeroTol
unsigned int GetNumPoints() const
Definition: Points.h:268

Member Data Documentation

boost::shared_ptr<NodalUtilTetrahedron> Nektar::LibUtilities::NodalTetElec::m_util
private