Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
ExpList3D.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ExpList3D.cpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Expansion list 3D definition
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #include <iomanip>
37 #include <MultiRegions/ExpList3D.h>
38 
39 #include <LocalRegions/HexExp.h>
40 #include <LocalRegions/PrismExp.h>
41 #include <LocalRegions/PyrExp.h>
42 #include <LocalRegions/TetExp.h>
43 
46 
47 namespace Nektar
48 {
49  namespace MultiRegions
50  {
51 
53  {
54  SetExpType(e3D);
55  }
56 
58  {
59  SetExpType(e3D);
60  }
61 
63  const std::vector<unsigned int> &eIDs): ExpList(In, eIDs)
64  {
65  SetExpType(e3D);
66 
67  // Setup Default optimisation information.
68  int nel = GetExpSize();
71 
72  SetCoeffPhys();
73 
76  }
77 
79  {
80  }
81 
82 
84  const LibUtilities::BasisKey &TBa,
85  const LibUtilities::BasisKey &TBb,
86  const LibUtilities::BasisKey &TBc,
87  const LibUtilities::BasisKey &HBa,
88  const LibUtilities::BasisKey &HBb,
89  const LibUtilities::BasisKey &HBc,
91  const LibUtilities::PointsType TetNb):
92  ExpList(pSession,graph3D)
93  {
94  SetExpType(e3D);
95 
100 
101  const SpatialDomains::ExpansionMap &expansions = graph3D->GetExpansions();
102 
103  SpatialDomains::ExpansionMap::const_iterator expIt;
104  for (expIt = expansions.begin(); expIt != expansions.end(); ++expIt)
105  {
110 
111  if((TetGeom = boost::dynamic_pointer_cast<SpatialDomains::TetGeom>(expIt->second->m_geomShPtr)))
112  {
114  {
115 // Ntet = MemoryManager<LocalRegions::NodalTetExp>::AllocateSharedPtr(TetBa,TetBb,TetBc,TetNb,TetGeom);
116 // (*m_exp).push_back(Ntet);
117  }
118  else
119  {
121  (*m_exp).push_back(tet);
122  }
123 
125 
126  m_npoints += TBa.GetNumPoints()*TBb.GetNumPoints()*TBc.GetNumPoints();
127  }
128 /*
129  else if((PrismGeom = boost::dynamic_pointer_cast<SpatialDomains::PrismGeom>(expansions[i]->m_geomShPtr)))
130  {
131  prism = MemoryManager<LocalRegions::PrismExp>::AllocateSharedPtr(Ba,Bb,Bc,PrismGeom);
132  (*m_exp).push_back(prism);
133 
134  m_ncoeffs += StdRegions::StdPrismData::getNumberOfCoefficients(Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes());
135  m_npoints += Ba.GetNumPoints()*Bb.GetNumPoints()*Bc.GetNumPoints();
136 
137  }
138  else if((PyrGeom = boost::dynamic_pointer_cast<SpatialDomains::PyrGeom>(expansions[i]->m_geomShPtr)))
139  {
140  pyramid = MemoryManager<LocalRegions::PyrExp>::AllocateSharedPtr(Ba,Bb,Bc,PyrGeom);
141  (*m_exp).push_back(pyramid);
142 
143  m_ncoeffs += StdRegions::StdPyrData::getNumberOfCoefficients(Ba.GetNumModes(), Bb.GetNumModes(), Bc.GetNumModes());
144  m_npoints += Ba.GetNumPoints()*Bb.GetNumPoints()*Bc.GetNumPoints();
145 
146  }
147 */
148  else if((HexGeom = boost::dynamic_pointer_cast<SpatialDomains::HexGeom>(expIt->second->m_geomShPtr)))
149  {
150  hex = MemoryManager<LocalRegions::HexExp>::AllocateSharedPtr(HBa,HBb,HBc, HexGeom);
151  (*m_exp).push_back(hex);
152 
153  m_ncoeffs += HBa.GetNumModes()*HBb.GetNumModes()*HBc.GetNumModes();
154  m_npoints += HBa.GetNumPoints()*HBb.GetNumPoints()*HBc.GetNumPoints();
155  }
156  else
157  {
158  ASSERTL0(false,"dynamic cast to a proper Geometry3D failed");
159  }
160 
161  }
162 
163  // Setup Default optimisation information.
164  int nel = GetExpSize();
167 
168  SetCoeffPhys();
169 
172  }
173 
174  /**
175  * Given a mesh \a graph3D, containing information about the domain and
176  * the spectral/hp element expansion, this constructor fills the list
177  * of local expansions \texttt{m_exp} with the proper expansions,
178  * calculates the total number of quadrature points
179  * \f$\boldsymbol{x}_i\f$ and the local expansion coefficients
180  * \f$\hat{u}^e_n\f$ and allocates memory for the arrays #m_coeffs and
181  * #m_phys.
182  *
183  * @param graph3D A mesh, containing information about the domain
184  * and the spectral/hp element expansion.
185  */
187  const SpatialDomains::MeshGraphSharedPtr &graph3D,
188  const std::string &variable) :
189  ExpList(pSession,graph3D)
190  {
191  SetExpType(e3D);
192 
193  int elmtid = 0;
198 
199  const SpatialDomains::ExpansionMap &expansions
200  = graph3D->GetExpansions(variable);
201 
202  SpatialDomains::ExpansionMap::const_iterator expIt;
203  for (expIt = expansions.begin(); expIt != expansions.end(); ++expIt)
204  {
209 
210  if((TetGeom = boost::dynamic_pointer_cast<
211  SpatialDomains::TetGeom>(expIt->second->m_geomShPtr)))
212  {
214  = expIt->second->m_basisKeyVector[0];
216  = expIt->second->m_basisKeyVector[1];
218  = expIt->second->m_basisKeyVector[2];
219 
222  {
223  ASSERTL0(false,"LocalRegions::NodalTetExp is not "
224  "implemented yet");
225  }
226  else
227  {
229  ::AllocateSharedPtr(TetBa,TetBb,TetBc,
230  TetGeom);
231  tet->SetElmtId(elmtid++);
232  (*m_exp).push_back(tet);
233  }
234  }
235  else if((PrismGeom = boost::dynamic_pointer_cast<SpatialDomains
236  ::PrismGeom>(expIt->second->m_geomShPtr)))
237  {
238  LibUtilities::BasisKey PrismBa
239  = expIt->second->m_basisKeyVector[0];
240  LibUtilities::BasisKey PrismBb
241  = expIt->second->m_basisKeyVector[1];
242  LibUtilities::BasisKey PrismBc
243  = expIt->second->m_basisKeyVector[2];
244 
246  ::AllocateSharedPtr(PrismBa,PrismBb,
247  PrismBc,PrismGeom);
248  prism->SetElmtId(elmtid++);
249  (*m_exp).push_back(prism);
250  }
251  else if((PyrGeom = boost::dynamic_pointer_cast<
252  SpatialDomains::PyrGeom>(expIt->second->m_geomShPtr)))
253  {
255  = expIt->second->m_basisKeyVector[0];
257  = expIt->second->m_basisKeyVector[1];
259  = expIt->second->m_basisKeyVector[2];
260 
262  ::AllocateSharedPtr(PyrBa,PyrBb,PyrBc,
263  PyrGeom);
264  pyramid->SetElmtId(elmtid++);
265  (*m_exp).push_back(pyramid);
266  }
267  else if((HexGeom = boost::dynamic_pointer_cast<
268  SpatialDomains::HexGeom>(expIt->second->m_geomShPtr)))
269  {
271  = expIt->second->m_basisKeyVector[0];
273  = expIt->second->m_basisKeyVector[1];
275  = expIt->second->m_basisKeyVector[2];
276 
278  ::AllocateSharedPtr(HexBa,HexBb,HexBc,
279  HexGeom);
280  hex->SetElmtId(elmtid++);
281  (*m_exp).push_back(hex);
282  }
283  else
284  {
285  ASSERTL0(false,"dynamic cast to a proper Geometry3D "
286  "failed");
287  }
288 
289  }
290 
291  // Setup Default optimisation information.
292  int nel = GetExpSize();
295 
296  SetCoeffPhys();
299  }
300 
301  /**
302  * Given an expansion vector \a expansions, containing
303  * information about the domain and the spectral/hp element
304  * expansion, this constructor fills the list of local
305  * expansions \texttt{m_exp} with the proper expansions,
306  * calculates the total number of quadrature points
307  * \f$\boldsymbol{x}_i\f$ and the local expansion coefficients
308  * \f$\hat{u}^e_n\f$ and allocates memory for the arrays
309  * #m_coeffs and #m_phys.
310  *
311  * @param expansions An expansion vector, containing
312  * information about the domain and the
313  * spectral/hp element expansion.
314  */
316  ExpList()
317  {
318  SetExpType(e3D);
319 
320  int elmtid = 0;
325 
326 
327  for(int i = 0; i < expansions.size(); ++i)
328  {
333 
334  SpatialDomains::ExpansionMap::const_iterator expmap = expansions.find(i);
335  ASSERTL1(expmap != expansions.end(), "Unable to find expansion.");
336  const SpatialDomains::ExpansionShPtr exp = expmap->second;
337 
338  if((TetGeom = boost::dynamic_pointer_cast<
339  SpatialDomains::TetGeom>(exp->m_geomShPtr)))
340  {
342  = exp->m_basisKeyVector[0];
344  = exp->m_basisKeyVector[1];
346  = exp->m_basisKeyVector[2];
347 
350  {
351  ASSERTL0(false,"LocalRegions::NodalTetExp is not "
352  "implemented yet");
353  }
354  else
355  {
357  ::AllocateSharedPtr(TetBa,TetBb,TetBc,
358  TetGeom);
359  tet->SetElmtId(elmtid++);
360  (*m_exp).push_back(tet);
361  }
362  }
363  else if((PrismGeom = boost::dynamic_pointer_cast<
364  SpatialDomains::PrismGeom>(exp->m_geomShPtr)))
365  {
366  LibUtilities::BasisKey PrismBa
367  = exp->m_basisKeyVector[0];
368  LibUtilities::BasisKey PrismBb
369  = exp->m_basisKeyVector[1];
370  LibUtilities::BasisKey PrismBc
371  = exp->m_basisKeyVector[2];
372 
374  ::AllocateSharedPtr(PrismBa,PrismBb,
375  PrismBc,PrismGeom);
376  prism->SetElmtId(elmtid++);
377  (*m_exp).push_back(prism);
378  }
379  else if((PyrGeom = boost::dynamic_pointer_cast<
380  SpatialDomains::PyrGeom>(exp->m_geomShPtr)))
381  {
383  = exp->m_basisKeyVector[0];
385  = exp->m_basisKeyVector[1];
387  = exp->m_basisKeyVector[2];
388 
390  ::AllocateSharedPtr(PyrBa,PyrBb,PyrBc,
391  PyrGeom);
392  pyramid->SetElmtId(elmtid++);
393  (*m_exp).push_back(pyramid);
394  }
395  else if((HexGeom = boost::dynamic_pointer_cast<
396  SpatialDomains::HexGeom>(exp->m_geomShPtr)))
397  {
399  = exp->m_basisKeyVector[0];
401  = exp->m_basisKeyVector[1];
403  = exp->m_basisKeyVector[2];
404 
406  ::AllocateSharedPtr(HexBa,HexBb,HexBc,
407  HexGeom);
408  hex->SetElmtId(elmtid++);
409  (*m_exp).push_back(hex);
410  }
411  else
412  {
413  ASSERTL0(false,"dynamic cast to a proper Geometry3D "
414  "failed");
415  }
416 
417  }
418 
419  // Setup Default optimisation information.
420  int nel = GetExpSize();
423 
424  SetCoeffPhys();
426  }
427 
428  /**
429  * Set up the storage for the concatenated list of
430  * coefficients and physical evaluations at the quadrature
431  * points. Each expansion (local element) is processed in turn
432  * to determine the number of coefficients and physical data
433  * points it contributes to the domain. Three arrays,
434  * #m_coeff_offset, #m_phys_offset and #m_offset_elmt_id, are
435  * also initialised and updated to store the data offsets of
436  * each element in the #m_coeffs and #m_phys arrays, and the
437  * element id that each consecutive block is associated
438  * respectively.
439  */
441  {
442  int i;
443 
444  // Set up offset information and array sizes
448 
449  m_ncoeffs = m_npoints = 0;
450 
451  for(i = 0; i < m_exp->size(); ++i)
452  {
454  m_phys_offset [i] = m_npoints;
455  m_offset_elmt_id[i] = i;
456  m_ncoeffs += (*m_exp)[i]->GetNcoeffs();
457  m_npoints += (*m_exp)[i]->GetTotPoints();
458  }
459 
462  }
463 
465  {
466  Array<OneD, int> NumShape(4,0);
467 
468  for(int i = 0; i < GetExpSize(); ++i)
469  {
470  switch ((*m_exp)[i]->DetShapeType())
471  {
472  case LibUtilities::eTetrahedron: NumShape[0]++; break;
473  case LibUtilities::ePyramid: NumShape[1]++; break;
474  case LibUtilities::ePrism: NumShape[2]++; break;
475  case LibUtilities::eHexahedron: NumShape[3]++; break;
476  default:
477  ASSERTL0(false, "Unknown expansion type.");
478  break;
479  }
480  }
481 
482  int three = 3;
484  ::AllocateSharedPtr(m_session,three,NumShape);
485  }
486 
487  void ExpList3D::v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
488  {
489  int i,j,k;
490  int nquad0 = (*m_exp)[expansion]->GetNumPoints(0);
491  int nquad1 = (*m_exp)[expansion]->GetNumPoints(1);
492  int nquad2 = (*m_exp)[expansion]->GetNumPoints(2);
493  int ntot = nquad0*nquad1*nquad2;
494  int ntotminus = (nquad0-1)*(nquad1-1)*(nquad2-1);
495 
496  Array<OneD,NekDouble> coords[3];
497  coords[0] = Array<OneD,NekDouble>(ntot);
498  coords[1] = Array<OneD,NekDouble>(ntot);
499  coords[2] = Array<OneD,NekDouble>(ntot);
500  (*m_exp)[expansion]->GetCoords(coords[0],coords[1],coords[2]);
501 
502  outfile << " <Piece NumberOfPoints=\""
503  << ntot << "\" NumberOfCells=\""
504  << ntotminus << "\">" << endl;
505  outfile << " <Points>" << endl;
506  outfile << " <DataArray type=\"Float64\" "
507  << "NumberOfComponents=\"3\" format=\"ascii\">" << endl;
508  outfile << " ";
509  for (i = 0; i < ntot; ++i)
510  {
511  for (j = 0; j < 3; ++j)
512  {
513  outfile << setprecision(8) << scientific
514  << (float)coords[j][i] << " ";
515  }
516  outfile << endl;
517  }
518 
519  outfile << endl;
520  outfile << " </DataArray>" << endl;
521  outfile << " </Points>" << endl;
522  outfile << " <Cells>" << endl;
523  outfile << " <DataArray type=\"Int32\" "
524  << "Name=\"connectivity\" format=\"ascii\">" << endl;
525  for (i = 0; i < nquad0-1; ++i)
526  {
527  for (j = 0; j < nquad1-1; ++j)
528  {
529  for (k = 0; k < nquad2-1; ++k)
530  {
531  outfile << k*nquad0*nquad1 + j*nquad0 + i << " "
532  << k*nquad0*nquad1 + j*nquad0 + i + 1 << " "
533  << k*nquad0*nquad1 + (j+1)*nquad0 + i + 1 << " "
534  << k*nquad0*nquad1 + (j+1)*nquad0 + i << " "
535  << (k+1)*nquad0*nquad1 + j*nquad0 + i << " "
536  << (k+1)*nquad0*nquad1 + j*nquad0 + i + 1 << " "
537  << (k+1)*nquad0*nquad1 + (j+1)*nquad0 + i + 1 << " "
538  << (k+1)*nquad0*nquad1 + (j+1)*nquad0 + i << " " << endl;
539  }
540  }
541  }
542  outfile << endl;
543  outfile << " </DataArray>" << endl;
544  outfile << " <DataArray type=\"Int32\" "
545  << "Name=\"offsets\" format=\"ascii\">" << endl;
546  for (i = 0; i < ntotminus; ++i)
547  {
548  outfile << i*8+8 << " ";
549  }
550  outfile << endl;
551  outfile << " </DataArray>" << endl;
552  outfile << " <DataArray type=\"UInt8\" "
553  << "Name=\"types\" format=\"ascii\">" << endl;
554  for (i = 0; i < ntotminus; ++i)
555  {
556  outfile << "12 ";
557  }
558  outfile << endl;
559  outfile << " </DataArray>" << endl;
560  outfile << " </Cells>" << endl;
561  outfile << " <PointData>" << endl;
562  }
563 
565  {
566  int i, j;
567  for (i = 0; i < m_exp->size(); ++i)
568  {
569  for (j = 0; j < (*m_exp)[i]->GetNfaces(); ++j)
570  {
571  (*m_exp)[i]->ComputeFaceNormal(j);
572  }
573  }
574  }
575 
577  const Array<OneD, NekDouble> &inarray,
578  Array<OneD, NekDouble> &outarray)
579  {
580  int cnt,cnt1;
581 
582  cnt = cnt1 = 0;
583  for(int i = 0; i < GetExpSize(); ++i)
584  {
585  // get new points key
586  int pt0 = (*m_exp)[i]->GetNumPoints(0);
587  int pt1 = (*m_exp)[i]->GetNumPoints(1);
588  int pt2 = (*m_exp)[i]->GetNumPoints(2);
589  int npt0 = (int) pt0*scale;
590  int npt1 = (int) pt1*scale;
591  int npt2 = (int) pt2*scale;
592 
593  LibUtilities::PointsKey newPointsKey0(npt0,(*m_exp)[i]->GetPointsType(0));
594  LibUtilities::PointsKey newPointsKey1(npt1,(*m_exp)[i]->GetPointsType(1));
595  LibUtilities::PointsKey newPointsKey2(npt2,(*m_exp)[i]->GetPointsType(2));
596 
597  // Interpolate points;
598  LibUtilities::Interp3D((*m_exp)[i]->GetBasis(0)->GetPointsKey(),
599  (*m_exp)[i]->GetBasis(1)->GetPointsKey(),
600  (*m_exp)[i]->GetBasis(2)->GetPointsKey(),
601  &inarray[cnt], newPointsKey0,
602  newPointsKey1, newPointsKey2,
603  &outarray[cnt1]);
604 
605  cnt += pt0*pt1*pt2;
606  cnt1 += npt0*npt1*npt2;
607  }
608  }
609 
611  const Array<OneD, NekDouble> &inarray,
612  Array<OneD, NekDouble> &outarray)
613  {
614  int cnt,cnt1;
615 
616  cnt = cnt1 = 0;
617  for(int i = 0; i < GetExpSize(); ++i)
618  {
619  // get new points key
620  int pt0 = (*m_exp)[i]->GetNumPoints(0);
621  int pt1 = (*m_exp)[i]->GetNumPoints(1);
622  int pt2 = (*m_exp)[i]->GetNumPoints(2);
623  int npt0 = (int) pt0*scale;
624  int npt1 = (int) pt1*scale;
625  int npt2 = (int) pt2*scale;
626 
627  LibUtilities::PointsKey newPointsKey0(npt0,(*m_exp)[i]->GetPointsType(0));
628  LibUtilities::PointsKey newPointsKey1(npt1,(*m_exp)[i]->GetPointsType(1));
629  LibUtilities::PointsKey newPointsKey2(npt2,(*m_exp)[i]->GetPointsType(2));
630 
631  // Project points;
633  newPointsKey1,
634  newPointsKey2,
635  &inarray[cnt],
636  (*m_exp)[i]->GetBasis(0)->GetPointsKey(),
637  (*m_exp)[i]->GetBasis(1)->GetPointsKey(),
638  (*m_exp)[i]->GetBasis(2)->GetPointsKey(),
639  &outarray[cnt1]);
640 
641  cnt += npt0*npt1*npt2;
642  cnt1 += pt0*pt1*pt2;
643  }
644 
645  }
646  } //end of namespace
647 } //end of namespace
648 
boost::shared_ptr< PyrGeom > PyrGeomSharedPtr
Definition: PyrGeom.h:84
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
virtual void v_SetUpPhysNormals()
Set up the normals on each expansion.
Definition: ExpList3D.cpp:564
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1001
BasisType GetBasisType() const
Return type of expansion basis.
Definition: Basis.h:139
Lagrange Polynomials using the Gauss points .
Definition: BasisType.h:54
boost::shared_ptr< HexExp > HexExpSharedPtr
Definition: HexExp.h:57
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:956
boost::shared_ptr< HexGeom > HexGeomSharedPtr
Definition: HexGeom.h:110
int getNumberOfCoefficients(int Na, int Nb, int Nc)
Definition: ShapeType.hpp:186
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:939
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1896
boost::shared_ptr< TetExp > TetExpSharedPtr
Definition: TetExp.h:223
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual void v_PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList3D.cpp:610
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:988
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:101
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:977
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:991
boost::shared_ptr< PyrExp > PyrExpSharedPtr
Definition: PyrExp.h:161
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:917
Array< OneD, int > m_offset_elmt_id
Array containing the element id m_offset_elmt_id[n] that the n^th consecutive block of data in m_coef...
Definition: ExpList.h:999
int GetNumPoints() const
Return points order at which basis is defined.
Definition: Basis.h:128
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:910
void PhysGalerkinProject3D(const BasisKey &fbasis0, const BasisKey &fbasis1, const BasisKey &fbasis2, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, const BasisKey &tbasis2, Array< OneD, NekDouble > &to)
Defines a specification for a set of points.
Definition: Points.h:58
double NekDouble
boost::shared_ptr< PrismExp > PrismExpSharedPtr
Definition: PrismExp.h:217
virtual void v_ReadGlobalOptimizationParameters()
Definition: ExpList3D.cpp:464
virtual ~ExpList3D()
Destructor.
Definition: ExpList3D.cpp:78
void Interp3D(const BasisKey &fbasis0, const BasisKey &fbasis1, const BasisKey &fbasis2, const Array< OneD, const NekDouble > &from, const BasisKey &tbasis0, const BasisKey &tbasis1, const BasisKey &tbasis2, Array< OneD, NekDouble > &to)
this function interpolates a 3D function evaluated at the quadrature points of the 3D basis...
Definition: Interp.cpp:186
boost::shared_ptr< Expansion > ExpansionShPtr
Definition: MeshGraph.h:173
boost::shared_ptr< PrismGeom > PrismGeomSharedPtr
Definition: PrismGeom.h:109
boost::shared_ptr< TetGeom > TetGeomSharedPtr
Definition: TetGeom.h:106
Abstraction of a three-dimensional multi-elemental expansion which is merely a collection of local ex...
Definition: ExpList3D.h:49
ExpList3D()
Default constructor.
Definition: ExpList3D.cpp:52
void CreateCollections(Collections::ImplementationType ImpType=Collections::eNoImpType)
Construct collections of elements containing a single element type and polynomial order from the list...
Definition: ExpList.cpp:2950
int GetNumModes() const
Returns the order of the basis.
Definition: Basis.h:84
Lagrange for SEM basis .
Definition: BasisType.h:53
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:251
void SetCoeffPhys(void)
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
Definition: ExpList3D.cpp:440
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:191
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
Describes the specification for a Basis.
Definition: Basis.h:50
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
Definition: ExpList3D.cpp:487
std::map< int, ExpansionShPtr > ExpansionMap
Definition: MeshGraph.h:174
virtual void v_PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList3D.cpp:576