Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ExpList3DHomogeneous1D.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ExpList3DHomogeneous1D.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: An ExpList which is homogeneous in 1 direction and so
33 // uses much of the functionality from a ExpList2D and its daughters
34 //
35 ///////////////////////////////////////////////////////////////////////////////
36 
38 #include <MultiRegions/ExpList2D.h>
39 
40 using namespace std;
41 
42 namespace Nektar
43 {
44  namespace MultiRegions
45  {
46  // Forward declaration for typedefs
47  ExpList3DHomogeneous1D::ExpList3DHomogeneous1D():
49  {
51  }
52 
55  const LibUtilities::BasisKey &HomoBasis,
56  const NekDouble lhom,
57  const bool useFFT,
58  const bool dealiasing):
59  ExpListHomogeneous1D(pSession,HomoBasis,lhom,useFFT,dealiasing)
60  {
62  }
63 
64  // Constructor for ExpList3DHomogeneous1D to act as a Explist2D field
67  const LibUtilities::BasisKey &HomoBasis,
68  const NekDouble lhom,
69  const bool useFFT,
70  const bool dealiasing,
72  const std::string &var,
73  const Collections::ImplementationType ImpType):
74  ExpListHomogeneous1D(pSession,HomoBasis,lhom,useFFT,dealiasing)
75  {
77 
78  GenExpList3DHomogeneous1D(graph2D->GetExpansions(var),ImpType);
79  }
80 
81  // Constructor for ExpList3DHomogeneous1D to act as a Explist2D field
84  const LibUtilities::BasisKey &HomoBasis,
85  const NekDouble lhom,
86  const bool useFFT,
87  const bool dealiasing,
88  const SpatialDomains::ExpansionMap &expansions,
89  const Collections::ImplementationType ImpType):
90  ExpListHomogeneous1D(pSession,HomoBasis,lhom,useFFT,dealiasing)
91  {
93 
94  GenExpList3DHomogeneous1D(expansions,ImpType);
95  }
96 
98  const SpatialDomains::ExpansionMap &expansions,
99  const Collections::ImplementationType ImpType)
100  {
101  int n,j,nel;
102  bool False = false;
103  ExpList2DSharedPtr plane_zero;
104 
105  // note that nzplanes can be larger than nzmodes
106  m_planes[0] = plane_zero = MemoryManager<ExpList2D>::AllocateSharedPtr(m_session, expansions, False,ImpType);
107 
109  nel = m_planes[0]->GetExpSize();
110 
111  for(j = 0; j < nel; ++j)
112  {
113  (*m_exp).push_back(m_planes[0]->GetExp(j));
114  }
115 
116  for(n = 1; n < m_planes.num_elements(); ++n)
117  {
119  for(j = 0; j < nel; ++j)
120  {
121  (*m_exp).push_back((*m_exp)[j]);
122  }
123  }
124 
125  // Setup Default optimisation information.
126  nel = GetExpSize();
129 
130  SetCoeffPhys();
131  }
132 
133 
134  /**
135  * @param In ExpList3DHomogeneous1D object to copy.
136  */
137  ExpList3DHomogeneous1D::ExpList3DHomogeneous1D(const ExpList3DHomogeneous1D &In, bool DeclarePlanesSetCoeffPhys):
139  {
141 
142  if(DeclarePlanesSetCoeffPhys)
143  {
144  bool False = false;
145  ExpList2DSharedPtr zero_plane = boost::dynamic_pointer_cast<ExpList2D> (In.m_planes[0]);
146 
147  for(int n = 0; n < m_planes.num_elements(); ++n)
148  {
150  }
151 
152  SetCoeffPhys();
153  }
154  }
155 
156  /**
157  * @param In ExpList3DHomogeneous1D object to copy.
158  * @param eIDs Id of elements that should be copied.
159  */
161  const ExpList3DHomogeneous1D &In,
162  const std::vector<unsigned int> &eIDs,
163  bool DeclarePlanesSetCoeffPhys,
164  const Collections::ImplementationType ImpType):
165  ExpListHomogeneous1D(In, eIDs)
166  {
168 
169  if(DeclarePlanesSetCoeffPhys)
170  {
171  bool False = false;
172  std::vector<unsigned int> eIDsPlane;
173  int nel = eIDs.size()/m_planes.num_elements();
174 
175  for (int i = 0; i < nel; ++i)
176  {
177  eIDsPlane.push_back(eIDs[i]);
178  }
179 
180  ExpList2DSharedPtr zero_plane_old =
181  boost::dynamic_pointer_cast<ExpList2D> (In.m_planes[0]);
182 
183  ExpList2DSharedPtr zero_plane =
184  MemoryManager<ExpList2D>::AllocateSharedPtr(*(zero_plane_old), eIDsPlane, ImpType);
185 
186  for(int n = 0; n < m_planes.num_elements(); ++n)
187  {
189  }
190 
191  SetCoeffPhys();
192  }
193  }
194 
195  /**
196  * Destructor
197  */
199  {
200  }
201 
203  {
204  int i,n,cnt;
205  int ncoeffs_per_plane = m_planes[0]->GetNcoeffs();
206  int npoints_per_plane = m_planes[0]->GetTotPoints();
207 
208  int nzplanes = m_planes.num_elements();
209 
210  // Set total coefficients and points
211  m_ncoeffs = ncoeffs_per_plane*nzplanes;
212  m_npoints = npoints_per_plane*nzplanes;
213 
216 
217  int nel = m_planes[0]->GetExpSize();
218  m_coeff_offset = Array<OneD,int>(nel*nzplanes);
219  m_phys_offset = Array<OneD,int>(nel*nzplanes);
220  Array<OneD, NekDouble> tmparray;
221 
222  for(cnt = n = 0; n < nzplanes; ++n)
223  {
224  m_planes[n]->SetCoeffsArray(tmparray= m_coeffs + ncoeffs_per_plane*n);
225  m_planes[n]->SetPhysArray(tmparray = m_phys + npoints_per_plane*n);
226 
227  for(i = 0; i < nel; ++i)
228  {
229  m_coeff_offset[cnt] = m_planes[n]->GetCoeff_Offset(i) + n*ncoeffs_per_plane;
230  m_phys_offset[cnt++] = m_planes[n]->GetPhys_Offset(i) + n*npoints_per_plane;
231  }
232  }
233  }
234 
239  {
240  int n;
241  Array<OneD, NekDouble> tmp_xc;
242  int nzplanes = m_planes.num_elements();
243  int npoints = GetTotPoints(eid);
244 
245  (*m_exp)[eid]->GetCoords(xc0,xc1);
246 
247  // Fill z-direction
249  Array<OneD, NekDouble> local_pts(m_planes.num_elements());
250 
251  for(n = 0; n < m_planes.num_elements(); n++)
252  {
253  local_pts[n] = pts[m_transposition->GetPlaneID(n)];
254  }
255 
256  Array<OneD, NekDouble> z(nzplanes);
257 
258  Vmath::Smul(nzplanes,m_lhom/2.0,local_pts,1,z,1);
259  Vmath::Sadd(nzplanes,m_lhom/2.0,z,1,z,1);
260 
261  for(n = 0; n < nzplanes; ++n)
262  {
263  Vmath::Fill(npoints,z[n],tmp_xc = xc2 + npoints*n,1);
264  if(n)
265  {
266  Vmath::Vcopy(npoints,xc0,1,tmp_xc = xc0+npoints*n,1);
267  Vmath::Vcopy(npoints,xc1,1,tmp_xc = xc1+npoints*n,1);
268  }
269  }
270  }
271 
272  /**
273  * The operation calls the 2D plane coordinates through the
274  * function ExpList#GetCoords and then evaluated the third
275  * coordinate using the member \a m_lhom
276  *
277  * @param coord_0 After calculation, the \f$x_1\f$ coordinate
278  * will be stored in this array.
279  *
280  * @param coord_1 After calculation, the \f$x_2\f$ coordinate
281  * will be stored in this array.
282  *
283  * @param coord_2 After calculation, the \f$x_3\f$ coordinate
284  * will be stored in this array. This
285  * coordinate is evaluated using the
286  * predefined value \a m_lhom
287  */
291  {
292  int n;
293  Array<OneD, NekDouble> tmp_xc;
294  int nzplanes = m_planes.num_elements();
295  int npoints = m_planes[0]->GetTotPoints();
296 
297  m_planes[0]->GetCoords(xc0,xc1);
298 
299  // Fill z-direction
301 
302  Array<OneD, NekDouble> local_pts(m_planes.num_elements());
303 
304  for(n = 0; n < m_planes.num_elements(); n++)
305  {
306  local_pts[n] = pts[m_transposition->GetPlaneID(n)];
307  }
308 
309  Array<OneD, NekDouble> z(nzplanes);
310 
311  Vmath::Smul(nzplanes,m_lhom/2.0,local_pts,1,z,1);
312  Vmath::Sadd(nzplanes,m_lhom/2.0,z,1,z,1);
313 
314  for(n = 0; n < nzplanes; ++n)
315  {
316  Vmath::Fill(npoints,z[n],tmp_xc = xc2 + npoints*n,1);
317  if(n)
318  {
319  Vmath::Vcopy(npoints,xc0,1,tmp_xc = xc0+npoints*n,1);
320  Vmath::Vcopy(npoints,xc1,1,tmp_xc = xc1+npoints*n,1);
321  }
322  }
323  }
324 
325  void ExpList3DHomogeneous1D::v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
326  {
327  ASSERTL0(expansion == -1, "Multi-zone output not supported for homogeneous expansions.");
328 
329  const int nPtsPlane = m_planes[0]->GetNpoints();
330  const int nElmt = m_planes[0]->GetExpSize();
331  const int nPlanes = m_planes.num_elements();
332 
333  int cnt = 0;
334  int cnt2 = 0;
335  for (int i = 0; i < nElmt; ++i)
336  {
337  const int np0 = (*m_exp)[i]->GetNumPoints(0);
338  const int np1 = (*m_exp)[i]->GetNumPoints(1);
339 
340  for (int n = 1; n < nPlanes; ++n)
341  {
342  const int o1 = (n-1) * nPtsPlane;
343  const int o2 = n * nPtsPlane;
344  for (int j = 1; j < np1; ++j)
345  {
346  for(int k = 1; k < np0; ++k)
347  {
348  outfile << cnt + (j-1)*np0 + (k-1) + o1 + 1 << " ";
349  outfile << cnt + (j-1)*np0 + (k-1) + o2 + 1 << " ";
350  outfile << cnt + (j-1)*np0 + k + o2 + 1 << " ";
351  outfile << cnt + (j-1)*np0 + k + o1 + 1 << " ";
352  outfile << cnt + j *np0 + (k-1) + o1 + 1 << " ";
353  outfile << cnt + j *np0 + (k-1) + o2 + 1 << " ";
354  outfile << cnt + j *np0 + k + o2 + 1 << " ";
355  outfile << cnt + j *np0 + k + o1 + 1 << endl;
356  cnt2++;
357  }
358  }
359  }
360 
361  cnt += np0*np1;
362  }
363  }
364 
365 
367  std::ostream &outfile,
368  int expansion,
369  int istrip)
370  {
371  // If there is only one plane (e.g. HalfMode), we write a 2D plane.
372  if (m_planes.num_elements() == 1)
373  {
374  m_planes[0]->WriteVtkPieceHeader(outfile, expansion);
375  return;
376  }
377 
378  // If we are using Fourier points, output extra plane to fill domain
379  int outputExtraPlane = 0;
380  if ( m_homogeneousBasis->GetBasisType() == LibUtilities::eFourier
381  && m_homogeneousBasis->GetPointsType() ==
383  {
384  outputExtraPlane = 1;
385  }
386 
387  int i,j,k;
388  int nq0 = (*m_exp)[expansion]->GetNumPoints(0);
389  int nq1 = (*m_exp)[expansion]->GetNumPoints(1);
390  int nq2 = m_planes.num_elements() + outputExtraPlane;
391  int ntot = nq0*nq1*nq2;
392  int ntotminus = (nq0-1)*(nq1-1)*(nq2-1);
393 
394  Array<OneD,NekDouble> coords[3];
395  coords[0] = Array<OneD,NekDouble>(ntot);
396  coords[1] = Array<OneD,NekDouble>(ntot);
397  coords[2] = Array<OneD,NekDouble>(ntot);
398  GetCoords(expansion,coords[0],coords[1],coords[2]);
399 
400  if (outputExtraPlane)
401  {
402  // Copy coords[0] and coords[1] to extra plane
404  Vmath::Vcopy (nq0*nq1, coords[0], 1,
405  tmp = coords[0] + (nq2-1)*nq0*nq1, 1);
406  Vmath::Vcopy (nq0*nq1, coords[1], 1,
407  tmp = coords[1] + (nq2-1)*nq0*nq1, 1);
408  // Fill coords[2] for extra plane
409  NekDouble z = coords[2][nq0*nq1*m_planes.num_elements()-1] +
410  (coords[2][nq0*nq1] - coords[2][0]);
411  Vmath::Fill(nq0*nq1, z, tmp = coords[2] + (nq2-1)*nq0*nq1, 1);
412  }
413 
414  NekDouble DistStrip;
415  m_session->LoadParameter("DistStrip", DistStrip, 0);
416  // Reset the z-coords for homostrips
417  for(int i = 0; i < ntot; i++)
418  {
419  coords[2][i] += istrip*DistStrip;
420  }
421 
422  outfile << " <Piece NumberOfPoints=\""
423  << ntot << "\" NumberOfCells=\""
424  << ntotminus << "\">" << endl;
425  outfile << " <Points>" << endl;
426  outfile << " <DataArray type=\"Float64\" "
427  << "NumberOfComponents=\"3\" format=\"ascii\">" << endl;
428  outfile << " ";
429  for (i = 0; i < ntot; ++i)
430  {
431  for (j = 0; j < 3; ++j)
432  {
433  outfile << coords[j][i] << " ";
434  }
435  outfile << endl;
436  }
437  outfile << endl;
438  outfile << " </DataArray>" << endl;
439  outfile << " </Points>" << endl;
440  outfile << " <Cells>" << endl;
441  outfile << " <DataArray type=\"Int32\" "
442  << "Name=\"connectivity\" format=\"ascii\">" << endl;
443  for (i = 0; i < nq0-1; ++i)
444  {
445  for (j = 0; j < nq1-1; ++j)
446  {
447  for (k = 0; k < nq2-1; ++k)
448  {
449  outfile << k*nq0*nq1 + j*nq0 + i << " "
450  << k*nq0*nq1 + j*nq0 + i + 1 << " "
451  << k*nq0*nq1 + (j+1)*nq0 + i + 1 << " "
452  << k*nq0*nq1 + (j+1)*nq0 + i << " "
453  << (k+1)*nq0*nq1 + j*nq0 + i << " "
454  << (k+1)*nq0*nq1 + j*nq0 + i + 1 << " "
455  << (k+1)*nq0*nq1 + (j+1)*nq0 + i + 1 << " "
456  << (k+1)*nq0*nq1 + (j+1)*nq0 + i << endl;
457  }
458  }
459  }
460  outfile << endl;
461  outfile << " </DataArray>" << endl;
462  outfile << " <DataArray type=\"Int32\" "
463  << "Name=\"offsets\" format=\"ascii\">" << endl;
464  for (i = 0; i < ntotminus; ++i)
465  {
466  outfile << i*8+8 << " ";
467  }
468  outfile << endl;
469  outfile << " </DataArray>" << endl;
470  outfile << " <DataArray type=\"UInt8\" "
471  << "Name=\"types\" format=\"ascii\">" << endl;
472  for (i = 0; i < ntotminus; ++i)
473  {
474  outfile << "12 ";
475  }
476  outfile << endl;
477  outfile << " </DataArray>" << endl;
478  outfile << " </Cells>" << endl;
479  outfile << " <PointData>" << endl;
480  }
481 
483  const Array<OneD, const NekDouble> &inarray,
484  const Array<OneD, const NekDouble> &soln)
485  {
486  int cnt = 0;
487  NekDouble errL2, err = 0.0;
489  Array<OneD, NekDouble> local_w(m_planes.num_elements());
490 
491  for(int n = 0; n < m_planes.num_elements(); n++)
492  {
493  local_w[n] = w[m_transposition->GetPlaneID(n)];
494  }
495 
496  if (soln == NullNekDouble1DArray)
497  {
498  for(int n = 0; n < m_planes.num_elements(); ++n)
499  {
500  errL2 = m_planes[n]->L2(inarray + cnt);
501  cnt += m_planes[n]->GetTotPoints();
502  err += errL2*errL2*local_w[n]*m_lhom*0.5;
503  }
504  }
505  else
506  {
507  for(int n = 0; n < m_planes.num_elements(); ++n)
508  {
509  errL2 = m_planes[n]->L2(inarray + cnt, soln + cnt);
510  cnt += m_planes[n]->GetTotPoints();
511  err += errL2*errL2*local_w[n]*m_lhom*0.5;
512  }
513  }
514  m_comm->GetColumnComm()->AllReduce(err, LibUtilities::ReduceSum);
515 
516  return sqrt(err);
517  }
518 
520  {
521  Array<OneD, NekDouble> energy(m_planes.num_elements()/2);
522  NekDouble area = 0.0;
523 
524  // Calculate total area of elements.
525  for (int n = 0; n < m_planes[0]->GetExpSize(); ++n)
526  {
527  Array<OneD, NekDouble> inarray(m_planes[0]->GetExp(n)->GetTotPoints(), 1.0);
528  area += m_planes[0]->GetExp(n)->Integral(inarray);
529  }
530 
531  m_comm->GetRowComm()->AllReduce(area, LibUtilities::ReduceSum);
532 
533  // Calculate L2 norm of real/imaginary planes.
534  for (int n = 0; n < m_planes.num_elements(); n += 2)
535  {
536  NekDouble err;
537 
538  energy[n/2] = 0;
539 
540  for(int i = 0; i < m_planes[n]->GetExpSize(); ++i)
541  {
542  StdRegions::StdExpansionSharedPtr exp = m_planes[n]->GetExp(i);
543  Array<OneD, NekDouble> phys(exp->GetTotPoints());
544  exp->BwdTrans(m_planes[n]->GetCoeffs()+m_planes[n]->GetCoeff_Offset(i),
545  phys);
546  err = exp->L2(phys);
547  energy[n/2] += err*err;
548 
549  exp = m_planes[n+1]->GetExp(i);
550  exp->BwdTrans(m_planes[n+1]->GetCoeffs()+m_planes[n+1]->GetCoeff_Offset(i),
551  phys);
552  err = exp->L2(phys);
553  energy[n/2] += err*err;
554  }
555 
556  m_comm->GetRowComm()->AllReduce(energy[n/2], LibUtilities::ReduceSum);
557  energy[n/2] /= 2.0*area;
558  }
559 
560  return energy;
561  }
562  } //end of namespace
563 } //end of namespace
564 
565 
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
const Array< OneD, const NekDouble > & GetCoeffs() const
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:1938
virtual void v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int GetCoeff_Offset(int n) const
Get the start offset position for a global list of m_coeffs correspoinding to element n...
Definition: ExpList.h:2076
virtual Array< OneD, const NekDouble > v_HomogeneousEnergy(void)
static Array< OneD, NekDouble > NullNekDouble1DArray
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1052
LibUtilities::TranspositionSharedPtr m_transposition
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:46
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:2067
STL namespace.
NekDouble m_lhom
Width of homogeneous direction.
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:1015
Fourier Expansion .
Definition: BasisType.h:52
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:998
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:2046
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual void v_GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2)
void SetCoeffPhys(void)
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:1047
LibUtilities::BasisSharedPtr m_homogeneousBasis
Definition of the total number of degrees of freedom and quadrature points. Sets up the storage for m...
1D Evenly-spaced points using Fourier Fit
Definition: PointsType.h:66
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:1036
int GetTotPoints(void) const
Returns the total number of quadrature points m_npoints .
Definition: ExpList.h:1535
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:1050
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:976
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:969
Array< OneD, ExpListSharedPtr > m_planes
double NekDouble
boost::shared_ptr< ExpList2D > ExpList2DSharedPtr
Shared pointer to an ExpList2D object.
Definition: ExpList2D.h:49
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:315
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: ExpList.h:966
void GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1=NullNekDouble1DArray, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
This function calculates the coordinates of all the elemental quadrature points . ...
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
Definition: ExpList2D.h:60
void GenExpList3DHomogeneous1D(const SpatialDomains::ExpansionMap &expansions, const Collections::ImplementationType ImpType)
virtual NekDouble v_L2(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:277
boost::shared_ptr< StdExpansion > StdExpansionSharedPtr
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
Abstraction of a two-dimensional multi-elemental expansion which is merely a collection of local expa...
Describes the specification for a Basis.
Definition: Basis.h:50
std::map< int, ExpansionShPtr > ExpansionMap
Definition: MeshGraph.h:174