Nektar++
Public Member Functions | Static Public Member Functions | Private Member Functions | Private Attributes | Static Private Attributes | List of all members
Nektar::LibUtilities::NodalTetElec Class Reference

#include <NodalTetElec.h>

Inheritance diagram for Nektar::LibUtilities::NodalTetElec:
[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, 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 std::shared_ptr< PointsBaseTypeCreate (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

std::shared_ptr< NodalUtilTetrahedronm_util
 

Static Private Attributes

static bool initPointsManager []
 

Additional Inherited Members

- Public Types inherited from Nektar::LibUtilities::Points< NekDouble >
typedef NekDouble DataType
 
typedef std::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::opLessm_InterpManager
 
NekManager< PointsKey, NekMatrix< DataType >, PointsKey::opLessm_GalerkinProjectionManager
 

Detailed Description

Definition at line 51 of file NodalTetElec.h.

Constructor & Destructor Documentation

◆ ~NodalTetElec()

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

Definition at line 54 of file NodalTetElec.h.

References Create(), and LIB_UTILITIES_EXPORT.

55  {
56  }

◆ NodalTetElec() [1/2]

Nektar::LibUtilities::NodalTetElec::NodalTetElec ( const PointsKey key)
inline

Definition at line 87 of file NodalTetElec.h.

87  : PointsBaseType(key)
88  {
89 
90  }
Points< NekDouble > PointsBaseType

◆ NodalTetElec() [2/2]

Nektar::LibUtilities::NodalTetElec::NodalTetElec ( )
inlineprivate

Definition at line 97 of file NodalTetElec.h.

References CalculateDerivMatrix(), CalculateInterpMatrix(), CalculatePoints(), CalculateWeights(), and NodalPointReorder3d().

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

Member Function Documentation

◆ CalculateDerivMatrix()

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

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

Definition at line 232 of file NodalTetElec.cpp.

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

Referenced by NodalTetElec().

233  {
234  // Allocate the derivative matrix.
236 
237  m_derivmatrix[0] = m_util->GetDerivMatrix(0);
238  m_derivmatrix[1] = m_util->GetDerivMatrix(1);
239  m_derivmatrix[2] = m_util->GetDerivMatrix(2);
240  }
std::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:95
MatrixSharedPtrType m_derivmatrix[3]
Definition: Points.h:383

◆ CalculateInterpMatrix()

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 217 of file NodalTetElec.cpp.

References m_util, and Vmath::Vcopy().

Referenced by GetI(), and NodalTetElec().

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

◆ CalculatePoints()

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

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

Definition at line 53 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.

Referenced by NodalTetElec().

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

◆ CalculateWeights()

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

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

Definition at line 203 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.

Referenced by NodalTetElec().

204  {
205  // Allocate the storage for points
207 
208  typedef DataType T;
209 
210  // Solve the Vandermonde system of integrals for the weight vector
211  NekVector<T> w = m_util->GetWeights();
212  m_weights = Array<OneD,T>( w.GetRows(), w.GetPtr() );
213  }
Array< OneD, DataType > m_weights
Definition: Points.h:382
std::shared_ptr< NodalUtilTetrahedron > m_util
Definition: NodalTetElec.h:95

◆ Create()

std::shared_ptr< PointsBaseType > Nektar::LibUtilities::NodalTetElec::Create ( const PointsKey key)
static

Definition at line 242 of file NodalTetElec.cpp.

Referenced by ~NodalTetElec().

243  {
244  std::shared_ptr<PointsBaseType> returnval(MemoryManager<NodalTetElec>::AllocateSharedPtr(key));
245  returnval->Initialize();
246  return returnval;
247  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.

◆ GetI() [1/2]

const MatrixSharedPtrType Nektar::LibUtilities::NodalTetElec::GetI ( const PointsKey pkey)
inlinevirtual

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

Definition at line 61 of file NodalTetElec.h.

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

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

◆ GetI() [2/2]

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 71 of file NodalTetElec.h.

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

75  {
76  size_t numpoints = x.num_elements();
77  unsigned int np = GetTotNumPoints();
78 
79  Array<OneD, NekDouble> interp(np*numpoints);
80  CalculateInterpMatrix(x, y, z, interp);
81 
82  NekDouble* d = interp.data();
83  return MemoryManager<NekMatrix<NekDouble> >
84  ::AllocateSharedPtr(numpoints, np, d);
85  }
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:277

◆ NodalPointReorder3d()

void Nektar::LibUtilities::NodalTetElec::NodalPointReorder3d ( )
private

Definition at line 249 of file NodalTetElec.cpp.

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

Referenced by CalculatePoints(), and NodalTetElec().

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

Member Data Documentation

◆ initPointsManager

bool Nektar::LibUtilities::NodalTetElec::initPointsManager
staticprivate
Initial value:

Definition at line 93 of file NodalTetElec.h.

◆ m_util

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