Nektar++
ExpList.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ExpList.h
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 top class definition
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef NEKTAR_LIBS_MULTIREGIONS_EXPLIST_H
37 #define NEKTAR_LIBS_MULTIREGIONS_EXPLIST_H
38 
43 #include <LocalRegions/Expansion.h>
44 #include <Collections/Collection.h>
51 #include <boost/enable_shared_from_this.hpp>
52 #include <tinyxml.h>
53 
54 namespace Nektar
55 {
56  namespace MultiRegions
57  {
58  // Forward declarations
59  class GlobalLinSys;
60  class AssemblyMapDG;
61 
62  class AssemblyMapCG;
63  class GlobalLinSysKey;
64  class GlobalMatrix;
65 
66  enum Direction
67  {
68  eX,
69  eY,
70  eZ,
71  eS,
73  };
74 
76  {
77  e0D,
78  e1D,
79  e2D,
82  e3D,
84  };
85 
87  {
88  eX,
89  eY,
90  eZ
91  };
92 
93  /// A map between global matrix keys and their associated block
94  /// matrices.
95  typedef map<GlobalMatrixKey,DNekScalBlkMatSharedPtr> BlockMatrixMap;
96  /// A shared pointer to a BlockMatrixMap.
97  typedef boost::shared_ptr<BlockMatrixMap> BlockMatrixMapShPtr;
98 
99 
100  /// Base class for all multi-elemental spectral/hp expansions.
101  class ExpList: public boost::enable_shared_from_this<ExpList>
102  {
103  public:
104  /// The default constructor.
106 
107  /// The default constructor.
109  const LibUtilities::SessionReaderSharedPtr &pSession);
110 
111  /// The default constructor.
113  const LibUtilities::SessionReaderSharedPtr &pSession,
114  const SpatialDomains::MeshGraphSharedPtr &pGraph);
115 
116  /// The copy constructor.
118  const ExpList &in,
119  const bool DeclareCoeffPhysArrays = true);
120 
121  /// The default destructor.
122  MULTI_REGIONS_EXPORT virtual ~ExpList();
123 
124  /// Returns the total number of local degrees of freedom
125  /// \f$N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_m\f$.
126  inline int GetNcoeffs(void) const;
127 
128  /// Returns the total number of local degrees of freedom
129  /// for element eid
130  MULTI_REGIONS_EXPORT int GetNcoeffs(const int eid) const;
131 
132  /// Returns the type of the expansion
134 
135  /// Returns the type of the expansion
137 
138  /// Evaulates the maximum number of modes in the elemental basis
139  /// order over all elements
140  inline int EvalBasisNumModesMax(void) const;
141 
142  /// Returns the vector of the number of modes in the elemental
143  /// basis order over all elements.
145  EvalBasisNumModesMaxPerExp(void) const;
146 
147  /// Returns the total number of quadrature points #m_npoints
148  /// \f$=Q_{\mathrm{tot}}\f$.
149  inline int GetTotPoints(void) const;
150 
151  /// Returns the total number of quadrature points for eid's element
152  /// \f$=Q_{\mathrm{tot}}\f$.
153  inline int GetTotPoints(const int eid) const;
154 
155  /// Returns the total number of quadrature points #m_npoints
156  /// \f$=Q_{\mathrm{tot}}\f$.
157  inline int GetNpoints(void) const;
158 
159 
160  /// Returns the total number of qudature points scaled by
161  /// the factor scale on each 1D direction
162  inline int Get1DScaledTotPoints(const NekDouble scale) const;
163 
164  /// Sets the wave space to the one of the possible configuration
165  /// true or false
166  inline void SetWaveSpace(const bool wavespace);
167 
168 
169  ///Set Modified Basis for the stability analysis
170  inline void SetModifiedBasis(const bool modbasis);
171 
172  /// Set the \a i th value of \a m_phys to value \a val
173  inline void SetPhys(int i, NekDouble val);
174 
175  /// This function returns the third direction expansion condition,
176  /// which can be in wave space (coefficient) or not
177  /// It is stored in the variable m_WaveSpace.
178  inline bool GetWaveSpace(void) const;
179 
180  /// Fills the array #m_phys
181  inline void SetPhys(const Array<OneD, const NekDouble> &inarray);
182 
183  /// Sets the array #m_phys
184  inline void SetPhysArray(Array<OneD, NekDouble> &inarray);
185 
186  /// This function manually sets whether the array of physical
187  /// values \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) is
188  /// filled or not.
189  inline void SetPhysState(const bool physState);
190 
191  /// This function indicates whether the array of physical values
192  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) is filled or
193  /// not.
194  inline bool GetPhysState(void) const;
195 
196  /// This function integrates a function \f$f(\boldsymbol{x})\f$
197  /// over the domain consisting of all the elements of the expansion.
199 
200  /// This function integrates a function \f$f(\boldsymbol{x})\f$
201  /// over the domain consisting of all the elements of the expansion.
203  const Array<OneD,
204  const NekDouble> &inarray);
205 
206  /// This function calculates the inner product of a function
207  /// \f$f(\boldsymbol{x})\f$ with respect to all \emph{local}
208  /// expansion modes \f$\phi_n^e(\boldsymbol{x})\f$.
209  inline void IProductWRTBase_IterPerExp(
210  const Array<OneD, const NekDouble> &inarray,
211  Array<OneD, NekDouble> &outarray);
212 
213  ///
214  inline void IProductWRTBase(
215  const Array<OneD, const NekDouble> &inarray,
216  Array<OneD, NekDouble> &outarray,
217  CoeffState coeffstate = eLocal);
218 
219  /// This function calculates the inner product of a function
220  /// \f$f(\boldsymbol{x})\f$ with respect to the derivative (in
221  /// direction \param dir) of all \emph{local} expansion modes
222  /// \f$\phi_n^e(\boldsymbol{x})\f$.
224  (const int dir,
225  const Array<OneD, const NekDouble> &inarray,
226  Array<OneD, NekDouble> &outarray);
227 
228  /// This function calculates the inner product of a function
229  /// \f$f(\boldsymbol{x})\f$ with respect to the derivative (in
230  /// direction \param dir) of all \emph{local} expansion modes
231  /// \f$\phi_n^e(\boldsymbol{x})\f$.
233  (const Array<OneD, const Array<OneD, NekDouble> > &inarray,
234  Array<OneD, NekDouble> &outarray);
235 
236  /// This function elementally evaluates the forward transformation
237  /// of a function \f$u(\boldsymbol{x})\f$ onto the global
238  /// spectral/hp expansion.
239  inline void FwdTrans_IterPerExp (
240  const Array<OneD,
241  const NekDouble> &inarray,
242  Array<OneD,NekDouble> &outarray);
243 
244  ///
245  inline void FwdTrans(
246  const Array<OneD,
247  const NekDouble> &inarray,
248  Array<OneD, NekDouble> &outarray,
249  CoeffState coeffstate = eLocal);
250 
251  /// This function elementally mulplies the coefficient space of
252  /// Sin my the elemental inverse of the mass matrix.
254  const Array<OneD,
255  const NekDouble> &inarray,
256  Array<OneD, NekDouble> &outarray);
257 
258  ///
259  inline void MultiplyByInvMassMatrix(
260  const Array<OneD,const NekDouble> &inarray,
261  Array<OneD, NekDouble> &outarray,
262  CoeffState coeffstate = eLocal);
263 
264  /// Smooth a field across elements
265  inline void SmoothField(Array<OneD,NekDouble> &field);
266 
267  /// Solve helmholtz problem
268  inline void HelmSolve(
269  const Array<OneD, const NekDouble> &inarray,
270  Array<OneD, NekDouble> &outarray,
271  const FlagList &flags,
272  const StdRegions::ConstFactorMap &factors,
273  const StdRegions::VarCoeffMap &varcoeff =
275  const Array<OneD, const NekDouble> &dirForcing =
277 
278  /// Solve Advection Diffusion Reaction
280  const Array<OneD, Array<OneD, NekDouble> > &velocity,
281  const Array<OneD, const NekDouble> &inarray,
282  Array<OneD, NekDouble> &outarray,
283  const NekDouble lambda,
284  CoeffState coeffstate = eLocal,
286  dirForcing = NullNekDouble1DArray);
287 
288 
289  /// Solve Advection Diffusion Reaction
290  inline void LinearAdvectionReactionSolve(
291  const Array<OneD, Array<OneD, NekDouble> > &velocity,
292  const Array<OneD, const NekDouble> &inarray,
293  Array<OneD, NekDouble> &outarray,
294  const NekDouble lambda,
295  CoeffState coeffstate = eLocal,
297  dirForcing = NullNekDouble1DArray);
298 
299  ///
301  const Array<OneD, const NekDouble> &inarray,
302  Array<OneD, NekDouble> &outarray);
303 
304 
305  /// This function elementally evaluates the backward transformation
306  /// of the global spectral/hp element expansion.
307  inline void BwdTrans_IterPerExp (
308  const Array<OneD, const NekDouble> &inarray,
309  Array<OneD,NekDouble> &outarray);
310 
311  ///
312  inline void BwdTrans (
313  const Array<OneD,
314  const NekDouble> &inarray,
315  Array<OneD,NekDouble> &outarray,
316  CoeffState coeffstate = eLocal);
317 
318  /// This function calculates the coordinates of all the elemental
319  /// quadrature points \f$\boldsymbol{x}_i\f$.
320  inline void GetCoords(
321  Array<OneD, NekDouble> &coord_0,
324 
325  // Homogeneous transforms
326  inline void HomogeneousFwdTrans(
327  const Array<OneD, const NekDouble> &inarray,
328  Array<OneD, NekDouble> &outarray,
329  CoeffState coeffstate = eLocal,
330  bool Shuff = true,
331  bool UnShuff = true);
332 
333  inline void HomogeneousBwdTrans(
334  const Array<OneD, const NekDouble> &inarray,
335  Array<OneD, NekDouble> &outarray,
336  CoeffState coeffstate = eLocal,
337  bool Shuff = true,
338  bool UnShuff = true);
339 
340  inline void DealiasedProd(
341  const Array<OneD, NekDouble> &inarray1,
342  const Array<OneD, NekDouble> &inarray2,
343  Array<OneD, NekDouble> &outarray,
344  CoeffState coeffstate = eLocal);
345 
346  inline void GetBCValues(
347  Array<OneD, NekDouble> &BndVals,
348  const Array<OneD, NekDouble> &TotField,
349  int BndID);
350 
351  inline void NormVectorIProductWRTBase(
354  Array<OneD, NekDouble> &outarray,
355  int BndID);
356 
357  /// Apply geometry information to each expansion.
359 
360  /// Reset geometry information and reset matrices
362  {
363  v_Reset();
364  }
365 
366  void WriteTecplotHeader(std::ostream &outfile,
367  std::string var = "")
368  {
369  v_WriteTecplotHeader(outfile, var);
370  }
371 
373  std::ostream &outfile,
374  int expansion = -1)
375  {
376  v_WriteTecplotZone(outfile, expansion);
377  }
378 
379  void WriteTecplotField(std::ostream &outfile,
380  int expansion = -1)
381  {
382  v_WriteTecplotField(outfile, expansion);
383  }
384 
385  void WriteTecplotConnectivity(std::ostream &outfile,
386  int expansion = -1)
387  {
388  v_WriteTecplotConnectivity(outfile, expansion);
389  }
390 
391  MULTI_REGIONS_EXPORT void WriteVtkHeader(std::ostream &outfile);
392  MULTI_REGIONS_EXPORT void WriteVtkFooter(std::ostream &outfile);
393 
394  void WriteVtkPieceHeader(std::ostream &outfile, int expansion,
395  int istrip = 0)
396  {
397  v_WriteVtkPieceHeader(outfile, expansion, istrip);
398  }
399 
401  std::ostream &outfile,
402  int expansion);
403 
405  std::ostream &outfile,
406  int expansion,
407  std::string var = "v")
408  {
409  v_WriteVtkPieceData(outfile, expansion, var);
410  }
411 
412  /// This function returns the dimension of the coordinates of the
413  /// element \a eid.
414  // inline
415  MULTI_REGIONS_EXPORT int GetCoordim(int eid);
416 
417  /// Set the \a i th coefficiient in \a m_coeffs to value \a val
418  inline void SetCoeff(int i, NekDouble val);
419 
420  /// Set the \a i th coefficiient in #m_coeffs to value \a val
421  inline void SetCoeffs(int i, NekDouble val);
422 
423  /// Set the #m_coeffs array to inarray
424  inline void SetCoeffsArray(Array<OneD, NekDouble> &inarray);
425 
426  /// This function returns (a reference to) the array
427  /// \f$\boldsymbol{\hat{u}}_l\f$ (implemented as #m_coeffs)
428  /// containing all local expansion coefficients.
429  inline const Array<OneD, const NekDouble> &GetCoeffs() const;
430 
431  /// Impose Dirichlet Boundary Conditions onto Array
432  inline void ImposeDirichletConditions(
433  Array<OneD,NekDouble>& outarray);
434 
435  /// Fill Bnd Condition expansion from the values stored in expansion
436  inline void FillBndCondFromField(void);
437 
438  /// Put the coefficients into global ordering using m_coeffs
439  inline void LocalToGlobal(void);
440 
441  /// Put the coefficients into local ordering and place in m_coeffs
442  inline void GlobalToLocal(void);
443 
444  /// Get the \a i th value (coefficient) of #m_coeffs
445  inline NekDouble GetCoeff(int i);
446 
447  /// Get the \a i th value (coefficient) of #m_coeffs
448  inline NekDouble GetCoeffs(int i);
449 
450  /// This function returns (a reference to) the array
451  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) containing the
452  /// function \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the
453  /// quadrature points.
454  // inline
456  &GetPhys() const;
457 
458  /// This function calculates the \f$L_\infty\f$ error of the global
459  /// spectral/hp element approximation.
461  const Array<OneD, const NekDouble> &inarray,
463 
464  /// This function calculates the \f$L_2\f$ error with
465  /// respect to soln of the global
466  /// spectral/hp element approximation.
468  const Array<OneD, const NekDouble> &inarray,
470  {
471  return v_L2(inarray, soln);
472  }
473 
474  /// Calculates the \f$H^1\f$ error of the global spectral/hp
475  /// element approximation.
477  const Array<OneD, const NekDouble> &inarray,
479 
481  {
482  return v_Integral(inarray);
483  }
484 
485  /// This function calculates the energy associated with
486  /// each one of the modesof a 3D homogeneous nD expansion
488  {
489  return v_HomogeneousEnergy();
490  }
491 
492  /// This function sets the Spectral Vanishing Viscosity
493  /// in homogeneous1D expansion.
495  {
497  }
498 
499  /// This function returns a vector containing the wave
500  /// numbers in z-direction associated
501  /// with the 3D homogenous expansion. Required if a
502  /// parellelisation is applied in the Fourier direction
504  {
505  return v_GetZIDs();
506  }
507 
508  /// This function returns the transposition class
509  /// associaed with the homogeneous expansion.
511  {
512  return v_GetTransposition();
513  }
514 
515  /// This function returns the Width of homogeneous direction
516  /// associaed with the homogeneous expansion.
518  {
519  return v_GetHomoLen();
520  }
521 
522  /// This function returns a vector containing the wave
523  /// numbers in y-direction associated
524  /// with the 3D homogenous expansion. Required if a
525  /// parellelisation is applied in the Fourier direction
527  {
528  return v_GetYIDs();
529  }
530 
531  /// This function interpolates the physical space points in
532  /// \a inarray to \a outarray using the same points defined in the
533  /// expansion but where the number of points are rescaled
534  /// by \a 1DScale
536  const NekDouble scale,
537  const Array<OneD, NekDouble> &inarray,
538  Array<OneD, NekDouble> &outarray)
539  {
540  v_PhysInterp1DScaled(scale, inarray,outarray);
541  }
542 
543  /// This function Galerkin projects the physical space points in
544  /// \a inarray to \a outarray where inarray is assumed to
545  /// be defined in the expansion but where the number of
546  /// points are rescaled by \a 1DScale
548  const NekDouble scale,
549  const Array<OneD, NekDouble> &inarray,
550  Array<OneD, NekDouble> &outarray)
551  {
552  v_PhysGalerkinProjection1DScaled(scale, inarray, outarray);
553  }
554 
555  /// This function returns the number of elements in the expansion.
556  inline int GetExpSize(void);
557 
558 
559  /// This function returns the number of elements in the
560  /// expansion which may be different for a homogeoenous extended
561  /// expansionp.
562  inline int GetNumElmts(void)
563  {
564  return v_GetNumElmts();
565  }
566 
567  /// This function returns the vector of elements in the expansion.
568  inline const boost::shared_ptr<LocalRegions::ExpansionVector>
569  GetExp() const;
570 
571  /// This function returns (a shared pointer to) the local elemental
572  /// expansion of the \f$n^{\mathrm{th}}\f$ element.
573  inline LocalRegions::ExpansionSharedPtr& GetExp(int n) const;
574 
575  /// This function returns (a shared pointer to) the local elemental
576  /// expansion containing the arbitrary point given by \a gloCoord.
578  const Array<OneD, const NekDouble> &gloCoord);
579 
580  /** This function returns the index of the local elemental
581  * expansion containing the arbitrary point given by \a gloCoord.
582  **/
584  const Array<OneD, const NekDouble> &gloCoord,
585  NekDouble tol = 0.0,
586  bool returnNearestElmt = false);
587 
588  /** This function returns the index and the Local
589  * Cartesian Coordinates \a locCoords of the local
590  * elemental expansion containing the arbitrary point
591  * given by \a gloCoords.
592  **/
594  const Array<OneD, const NekDouble> &gloCoords,
595  Array<OneD, NekDouble> &locCoords,
596  NekDouble tol = 0.0,
597  bool returnNearestElmt = false);
598 
599  /// Get the start offset position for a global list of #m_coeffs
600  /// correspoinding to element n.
601  inline int GetCoeff_Offset(int n) const;
602 
603  /// Get the start offset position for a global list of m_phys
604  /// correspoinding to element n.
605  inline int GetPhys_Offset(int n) const;
606 
607  /// Get the element id associated with the n th
608  /// consecutive block of data in #m_phys and #m_coeffs
609  inline int GetOffset_Elmt_Id(int n) const;
610 
611  /// This function returns (a reference to) the array
612  /// \f$\boldsymbol{\hat{u}}_l\f$ (implemented as #m_coeffs)
613  /// containing all local expansion coefficients.
615 
616  /// This function returns (a reference to) the array
617  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) containing the
618  /// function \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the
619  /// quadrature points.
621 
622  inline void PhysDeriv(
623  Direction edir,
624  const Array<OneD, const NekDouble> &inarray,
625  Array<OneD, NekDouble> &out_d);
626 
627  /// This function discretely evaluates the derivative of a function
628  /// \f$f(\boldsymbol{x})\f$ on the domain consisting of all
629  /// elements of the expansion.
630  inline void PhysDeriv(
631  const Array<OneD, const NekDouble> &inarray,
632  Array<OneD, NekDouble> &out_d0,
635 
636  inline void PhysDeriv(
637  const int dir,
638  const Array<OneD, const NekDouble> &inarray,
639  Array<OneD, NekDouble> &out_d);
640 
641 
642  // functions associated with DisContField
645 
646  inline boost::shared_ptr<ExpList> &UpdateBndCondExpansion(int i);
647 
648  inline void Upwind(
649  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
650  const Array<OneD, const NekDouble> &Fwd,
651  const Array<OneD, const NekDouble> &Bwd,
653 
654  inline void Upwind(
656  const Array<OneD, const NekDouble> &Fwd,
657  const Array<OneD, const NekDouble> &Bwd,
659 
660  /**
661  * Return a reference to the trace space associated with this
662  * expansion list.
663  */
664  inline boost::shared_ptr<ExpList> &GetTrace();
665 
666  inline boost::shared_ptr<AssemblyMapDG> &GetTraceMap(void);
667 
668  inline const Array<OneD, const int> &GetTraceBndMap(void);
669 
670  inline void GetNormals(Array<OneD, Array<OneD, NekDouble> > &normals);
671 
672  inline void AddTraceIntegral(
675  Array<OneD, NekDouble> &outarray);
676 
677  inline void AddTraceIntegral(
679  Array<OneD, NekDouble> &outarray);
680 
681  inline void AddFwdBwdTraceIntegral(
682  const Array<OneD, const NekDouble> &Fwd,
683  const Array<OneD, const NekDouble> &Bwd,
684  Array<OneD, NekDouble> &outarray);
685 
686  inline void GetFwdBwdTracePhys(
688  Array<OneD,NekDouble> &Bwd);
689 
690  inline void GetFwdBwdTracePhys(
691  const Array<OneD,const NekDouble> &field,
693  Array<OneD,NekDouble> &Bwd);
694 
695  inline void ExtractTracePhys(Array<OneD,NekDouble> &outarray);
696 
697  inline void ExtractTracePhys(
698  const Array<OneD, const NekDouble> &inarray,
699  Array<OneD,NekDouble> &outarray);
700 
701  inline const Array<OneD, const SpatialDomains::
703 
704  inline Array<OneD, SpatialDomains::
706 
707  inline void EvaluateBoundaryConditions(
708  const NekDouble time = 0.0,
709  const std::string varName = "",
712 
713  // Routines for continous matrix solution
714  /// This function calculates the result of the multiplication of a
715  /// matrix of type specified by \a mkey with a vector given by \a
716  /// inarray.
717  inline void GeneralMatrixOp(
718  const GlobalMatrixKey &gkey,
719  const Array<OneD,const NekDouble> &inarray,
720  Array<OneD, NekDouble> &outarray,
721  CoeffState coeffstate = eLocal);
722 
724  const GlobalMatrixKey &gkey,
725  const Array<OneD,const NekDouble> &inarray,
726  Array<OneD, NekDouble> &outarray);
727 
728  inline void SetUpPhysNormals();
729 
730  inline void GetBoundaryToElmtMap(Array<OneD, int> &ElmtID,
731  Array<OneD,int> &EdgeID);
732 
734  std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef,
735  int NumHomoDir = 0,
736  int NumHomoStrip = 1,
739  std::vector<NekDouble> &HomoLen =
741  std::vector<unsigned int> &HomoZIDs =
743  std::vector<unsigned int> &HomoYIDs =
745 
746 
748  {
749  return m_globalOptParam;
750  }
751 
752  map<int, RobinBCInfoSharedPtr> GetRobinBCInfo()
753  {
754  return v_GetRobinBCInfo();
755  }
756 
758  PeriodicMap &periodicVerts,
759  PeriodicMap &periodicEdges,
760  PeriodicMap &periodicFaces = NullPeriodicMap)
761  {
762  v_GetPeriodicEntities(periodicVerts, periodicEdges, periodicFaces);
763  }
764 
765  std::vector<LibUtilities::FieldDefinitionsSharedPtr>
767  {
768  return v_GetFieldDefinitions();
769  }
770 
771 
772  void GetFieldDefinitions(std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef)
773  {
774  v_GetFieldDefinitions(fielddef);
775  }
776 
777 
778 
779  /// Append the element data listed in elements
780  /// fielddef->m_ElementIDs onto fielddata
783  std::vector<NekDouble> &fielddata)
784  {
785  v_AppendFieldData(fielddef,fielddata);
786  }
787 
788 
789  /// Append the data in coeffs listed in elements
790  /// fielddef->m_ElementIDs onto fielddata
793  std::vector<NekDouble> &fielddata,
794  Array<OneD, NekDouble> &coeffs)
795  {
796  v_AppendFieldData(fielddef,fielddata,coeffs);
797  }
798 
799 
800  /** \brief Extract the data in fielddata into the coeffs
801  * using the basic ExpList Elemental expansions rather
802  * than planes in homogeneous case
803  */
806  std::vector<NekDouble> &fielddata,
807  std::string &field,
808  Array<OneD, NekDouble> &coeffs);
809 
810 
811  /** \brief Extract the data from fromField using
812  * fromExpList the coeffs using the basic ExpList
813  * Elemental expansions rather than planes in homogeneous
814  * case
815  */
817  const boost::shared_ptr<ExpList> &fromExpList,
818  const Array<OneD, const NekDouble> &fromCoeffs,
819  Array<OneD, NekDouble> &toCoeffs);
820 
821 
822  //Extract data in fielddata into the m_coeffs_list for the 3D stability analysis (base flow is 2D)
825  std::vector<NekDouble> &fielddata,
826  std::string &field,
827  Array<OneD, NekDouble> &coeffs);
828 
829 
830  /// Returns a shared pointer to the current object.
831  boost::shared_ptr<ExpList> GetSharedThisPtr()
832  {
833  return shared_from_this();
834  }
835 
836  /// Returns the session object
837  boost::shared_ptr<LibUtilities::SessionReader> GetSession()
838  {
839  return m_session;
840  }
841 
842  /// Returns the comm object
843  boost::shared_ptr<LibUtilities::Comm> GetComm()
844  {
845  return m_comm;
846  }
847 
849  {
850  return m_graph;
851  }
852 
853  // Wrapper functions for Homogeneous Expansions
855  {
856  return v_GetHomogeneousBasis();
857  }
858 
859  boost::shared_ptr<ExpList> &GetPlane(int n)
860  {
861  return v_GetPlane(n);
862  }
863 
864  //expansion type
866 
868  Collections::ImplementationType ImpType
869  = Collections::eNoImpType);
870 
871  protected:
872  boost::shared_ptr<DNekMat> GenGlobalMatrixFull(
873  const GlobalLinSysKey &mkey,
874  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
875 
876  /// Communicator
878 
879  /// Session
881 
882  /// Mesh associated with this expansion list.
884 
885  /// The total number of local degrees of freedom. #m_ncoeffs
886  /// \f$=N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_l\f$
888 
889  /** The total number of quadrature points. #m_npoints
890  *\f$=Q_{\mathrm{tot}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_Q\f$
891  **/
893 
894  /**
895  * \brief Concatenation of all local expansion coefficients.
896  *
897  * The array of length #m_ncoeffs\f$=N_{\mathrm{eof}}\f$ which is
898  * the concatenation of the local expansion coefficients
899  * \f$\hat{u}_n^e\f$ over all \f$N_{\mathrm{el}}\f$ elements
900  * \f[\mathrm{\texttt{m\_coeffs}}=\boldsymbol{\hat{u}}_{l} =
901  * \underline{\boldsymbol{\hat{u}}}^e = \left [ \begin{array}{c}
902  * \boldsymbol{\hat{u}}^{1} \ \
903  * \boldsymbol{\hat{u}}^{2} \ \
904  * \vdots \ \
905  * \boldsymbol{\hat{u}}^{{{N_{\mathrm{el}}}}} \end{array} \right ],
906  * \quad
907  * \mathrm{where}\quad \boldsymbol{\hat{u}}^{e}[n]=\hat{u}_n^{e}\f]
908  */
910 
911  /**
912  * \brief The global expansion evaluated at the quadrature points
913  *
914  * The array of length #m_npoints\f$=Q_{\mathrm{tot}}\f$ containing
915  * the evaluation of \f$u^{\delta}(\boldsymbol{x})\f$ at the
916  * quadrature points \f$\boldsymbol{x}_i\f$.
917  * \f[\mathrm{\texttt{m\_phys}}=\boldsymbol{u}_{l} =
918  * \underline{\boldsymbol{u}}^e = \left [ \begin{array}{c}
919  * \boldsymbol{u}^{1} \ \
920  * \boldsymbol{u}^{2} \ \
921  * \vdots \ \
922  * \boldsymbol{u}^{{{N_{\mathrm{el}}}}} \end{array} \right ],\quad
923  * \mathrm{where}\quad
924  * \boldsymbol{u}^{e}[i]=u^{\delta}(\boldsymbol{x}_i)\f]
925  */
927 
928  /**
929  * \brief The state of the array #m_phys.
930  *
931  * Indicates whether the array #m_phys, created to contain the
932  * evaluation of \f$u^{\delta}(\boldsymbol{x})\f$ at the quadrature
933  * points \f$\boldsymbol{x}_i\f$, is filled with these values.
934  */
936 
937  /**
938  * \brief The list of local expansions.
939  *
940  * The (shared pointer to the) vector containing (shared pointers
941  * to) all local expansions. The fact that the local expansions are
942  * all stored as a (pointer to a) #StdExpansion, the abstract base
943  * class for all local expansions, allows a general implementation
944  * where most of the routines for the derived classes are defined
945  * in the #ExpList base class.
946  */
947  boost::shared_ptr<LocalRegions::ExpansionVector> m_exp;
948 
949  Collections::CollectionVector m_collections;
950 
951  /// Offset of elemental data into the array #m_coeffs
952  std::vector<int> m_coll_coeff_offset;
953 
954  /// Offset of elemental data into the array #m_phys
955  std::vector<int> m_coll_phys_offset;
956 
957  /// Offset of elemental data into the array #m_coeffs
959 
960  /// Offset of elemental data into the array #m_phys
962 
963  /// Array containing the element id #m_offset_elmt_id[n]
964  /// that the n^th consecutive block of data in #m_coeffs
965  /// and #m_phys is associated, i.e. for an array of
966  /// constant expansion size and single shape elements
967  /// m_phys[n*m_npoints] is the data related to
968  /// m_exp[m_offset_elmt_id[n]];
970 
972 
973  BlockMatrixMapShPtr m_blockMat;
974 
975  //@todo should this be in ExpList or ExpListHomogeneous1D.cpp
976  // it's a bool which determine if the expansion is in the wave space (coefficient space)
977  // or not
979 
980  /// This function assembles the block diagonal matrix of local
981  /// matrices of the type \a mtype.
983  const GlobalMatrixKey &gkey);
984 
986  const GlobalMatrixKey &gkey);
987 
989  const GlobalMatrixKey &gkey,
990  const Array<OneD,const NekDouble> &inarray,
991  Array<OneD, NekDouble> &outarray);
992 
993  /// Generates a global matrix from the given key and map.
994  boost::shared_ptr<GlobalMatrix> GenGlobalMatrix(
995  const GlobalMatrixKey &mkey,
996  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
997 
998 
999  void GlobalEigenSystem(
1000  const boost::shared_ptr<DNekMat> &Gmat,
1001  Array<OneD, NekDouble> &EigValsReal,
1002  Array<OneD, NekDouble> &EigValsImag,
1003  Array<OneD, NekDouble> &EigVecs
1005 
1006 
1007  /// This operation constructs the global linear system of type \a
1008  /// mkey.
1009  boost::shared_ptr<GlobalLinSys> GenGlobalLinSys(
1010  const GlobalLinSysKey &mkey,
1011  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
1012 
1013  /// Generate a GlobalLinSys from information provided by the key
1014  /// "mkey" and the mapping provided in LocToGloBaseMap.
1015  boost::shared_ptr<GlobalLinSys> GenGlobalBndLinSys(
1016  const GlobalLinSysKey &mkey,
1017  const AssemblyMapSharedPtr &locToGloMap);
1018 
1020  {
1022  }
1023 
1024  // Virtual prototypes
1025 
1026  virtual int v_GetNumElmts(void)
1027  {
1028  return (*m_exp).size();
1029  }
1030 
1032  &v_GetBndCondExpansions(void);
1033 
1034  virtual boost::shared_ptr<ExpList> &v_UpdateBndCondExpansion(int i);
1035 
1036  virtual void v_Upwind(
1037  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
1038  const Array<OneD, const NekDouble> &Fwd,
1039  const Array<OneD, const NekDouble> &Bwd,
1041 
1042  virtual void v_Upwind(
1043  const Array<OneD, const NekDouble> &Vn,
1044  const Array<OneD, const NekDouble> &Fwd,
1045  const Array<OneD, const NekDouble> &Bwd,
1047 
1048  virtual boost::shared_ptr<ExpList> &v_GetTrace();
1049 
1050  virtual boost::shared_ptr<AssemblyMapDG> &v_GetTraceMap();
1051 
1052  virtual const Array<OneD, const int> &v_GetTraceBndMap();
1053 
1054  virtual void v_GetNormals(
1055  Array<OneD, Array<OneD, NekDouble> > &normals);
1056 
1057  virtual void v_AddTraceIntegral(
1058  const Array<OneD, const NekDouble> &Fx,
1059  const Array<OneD, const NekDouble> &Fy,
1060  Array<OneD, NekDouble> &outarray);
1061 
1062  virtual void v_AddTraceIntegral(
1063  const Array<OneD, const NekDouble> &Fn,
1064  Array<OneD, NekDouble> &outarray);
1065 
1066  virtual void v_AddFwdBwdTraceIntegral(
1067  const Array<OneD, const NekDouble> &Fwd,
1068  const Array<OneD, const NekDouble> &Bwd,
1069  Array<OneD, NekDouble> &outarray);
1070 
1071  virtual void v_GetFwdBwdTracePhys(
1072  Array<OneD,NekDouble> &Fwd,
1073  Array<OneD,NekDouble> &Bwd);
1074 
1075  virtual void v_GetFwdBwdTracePhys(
1076  const Array<OneD,const NekDouble> &field,
1077  Array<OneD,NekDouble> &Fwd,
1078  Array<OneD,NekDouble> &Bwd);
1079 
1080  virtual void v_ExtractTracePhys(
1081  Array<OneD,NekDouble> &outarray);
1082 
1083  virtual void v_ExtractTracePhys(
1084  const Array<OneD, const NekDouble> &inarray,
1085  Array<OneD,NekDouble> &outarray);
1086 
1087  virtual void v_MultiplyByInvMassMatrix(
1088  const Array<OneD,const NekDouble> &inarray,
1089  Array<OneD, NekDouble> &outarray,
1090  CoeffState coeffstate);
1091 
1092  virtual void v_HelmSolve(
1093  const Array<OneD, const NekDouble> &inarray,
1094  Array<OneD, NekDouble> &outarray,
1095  const FlagList &flags,
1096  const StdRegions::ConstFactorMap &factors,
1097  const StdRegions::VarCoeffMap &varcoeff,
1098  const Array<OneD, const NekDouble> &dirForcing);
1099 
1101  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1102  const Array<OneD, const NekDouble> &inarray,
1103  Array<OneD, NekDouble> &outarray,
1104  const NekDouble lambda,
1105  CoeffState coeffstate = eLocal,
1107  dirForcing = NullNekDouble1DArray);
1108 
1109  virtual void v_LinearAdvectionReactionSolve(
1110  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1111  const Array<OneD, const NekDouble> &inarray,
1112  Array<OneD, NekDouble> &outarray,
1113  const NekDouble lambda,
1114  CoeffState coeffstate = eLocal,
1116  dirForcing = NullNekDouble1DArray);
1117 
1118  // wrapper functions about virtual functions
1119  virtual void v_ImposeDirichletConditions(Array<OneD,NekDouble>& outarray);
1120 
1121  virtual void v_FillBndCondFromField();
1122 
1123  virtual void v_Reset();
1124 
1125  virtual void v_LocalToGlobal(void);
1126 
1127  virtual void v_GlobalToLocal(void);
1128 
1129  virtual void v_BwdTrans(
1130  const Array<OneD,const NekDouble> &inarray,
1131  Array<OneD, NekDouble> &outarray,
1132  CoeffState coeffstate);
1133 
1134  virtual void v_BwdTrans_IterPerExp(
1135  const Array<OneD,const NekDouble> &inarray,
1136  Array<OneD,NekDouble> &outarray);
1137 
1138  virtual void v_FwdTrans(
1139  const Array<OneD,const NekDouble> &inarray,
1140  Array<OneD, NekDouble> &outarray,
1141  CoeffState coeffstate);
1142 
1143  virtual void v_FwdTrans_IterPerExp(
1144  const Array<OneD,const NekDouble> &inarray,
1145  Array<OneD,NekDouble> &outarray);
1146 
1147  virtual void v_SmoothField(Array<OneD,NekDouble> &field);
1148 
1149  virtual void v_IProductWRTBase(
1150  const Array<OneD, const NekDouble> &inarray,
1151  Array<OneD, NekDouble> &outarray,
1152  CoeffState coeffstate);
1153 
1154  virtual void v_IProductWRTBase_IterPerExp(
1155  const Array<OneD,const NekDouble> &inarray,
1156  Array<OneD, NekDouble> &outarray);
1157 
1158  virtual void v_GeneralMatrixOp(
1159  const GlobalMatrixKey &gkey,
1160  const Array<OneD,const NekDouble> &inarray,
1161  Array<OneD, NekDouble> &outarray,
1162  CoeffState coeffstate);
1163 
1164  virtual void v_GetCoords(
1165  Array<OneD, NekDouble> &coord_0,
1166  Array<OneD, NekDouble> &coord_1,
1168 
1169  virtual void v_PhysDeriv(
1170  const Array<OneD, const NekDouble> &inarray,
1171  Array<OneD, NekDouble> &out_d0,
1172  Array<OneD, NekDouble> &out_d1,
1173  Array<OneD, NekDouble> &out_d2);
1174 
1175  virtual void v_PhysDeriv(
1176  const int dir,
1177  const Array<OneD, const NekDouble> &inarray,
1178  Array<OneD, NekDouble> &out_d);
1179 
1180  virtual void v_PhysDeriv(
1181  Direction edir,
1182  const Array<OneD, const NekDouble> &inarray,
1183  Array<OneD, NekDouble> &out_d);
1184 
1185  virtual void v_HomogeneousFwdTrans(
1186  const Array<OneD, const NekDouble> &inarray,
1187  Array<OneD, NekDouble> &outarray,
1188  CoeffState coeffstate = eLocal,
1189  bool Shuff = true,
1190  bool UnShuff = true);
1191 
1192  virtual void v_HomogeneousBwdTrans(
1193  const Array<OneD, const NekDouble> &inarray,
1194  Array<OneD, NekDouble> &outarray,
1195  CoeffState coeffstate = eLocal,
1196  bool Shuff = true,
1197  bool UnShuff = true);
1198 
1199  virtual void v_DealiasedProd(
1200  const Array<OneD, NekDouble> &inarray1,
1201  const Array<OneD, NekDouble> &inarray2,
1202  Array<OneD, NekDouble> &outarray,
1203  CoeffState coeffstate = eLocal);
1204 
1205  virtual void v_GetBCValues(
1206  Array<OneD, NekDouble> &BndVals,
1207  const Array<OneD, NekDouble> &TotField,
1208  int BndID);
1209 
1210  virtual void v_NormVectorIProductWRTBase(
1213  Array<OneD, NekDouble> &outarray,
1214  int BndID);
1215 
1216  virtual void v_SetUpPhysNormals();
1217 
1218  virtual void v_GetBoundaryToElmtMap(Array<OneD, int> &ElmtID,
1219  Array<OneD,int> &EdgeID);
1220 
1221  virtual void v_ReadGlobalOptimizationParameters();
1222 
1223  virtual std::vector<LibUtilities::FieldDefinitionsSharedPtr>
1224  v_GetFieldDefinitions(void);
1225 
1226  virtual void v_GetFieldDefinitions(
1227  std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef);
1228 
1229  virtual void v_AppendFieldData(
1231  std::vector<NekDouble> &fielddata);
1232 
1233  virtual void v_AppendFieldData(
1235  std::vector<NekDouble> &fielddata,
1236  Array<OneD, NekDouble> &coeffs);
1237 
1238  virtual void v_ExtractDataToCoeffs(
1240  std::vector<NekDouble> &fielddata, std::string &field,
1241  Array<OneD, NekDouble> &coeffs);
1242 
1243  virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr<ExpList> &fromExpList, const Array<OneD, const NekDouble> &fromCoeffs, Array<OneD, NekDouble> &toCoeffs);
1244 
1245  virtual void v_WriteTecplotHeader(std::ostream &outfile,
1246  std::string var = "");
1247  virtual void v_WriteTecplotZone(std::ostream &outfile,
1248  int expansion);
1249  virtual void v_WriteTecplotField(std::ostream &outfile,
1250  int expansion);
1251  virtual void v_WriteTecplotConnectivity(std::ostream &outfile,
1252  int expansion);
1253  virtual void v_WriteVtkPieceHeader(
1254  std::ostream &outfile,
1255  int expansion,
1256  int istrip);
1257 
1258  virtual void v_WriteVtkPieceData(
1259  std::ostream &outfile,
1260  int expansion,
1261  std::string var);
1262 
1263  virtual NekDouble v_L2(
1264  const Array<OneD, const NekDouble> &phys,
1266 
1267  virtual NekDouble v_Integral (
1268  const Array<OneD, const NekDouble> &inarray);
1269 
1272  virtual NekDouble v_GetHomoLen(void);
1275 
1276  // 1D Scaling and projection
1277  virtual void v_PhysInterp1DScaled(
1278  const NekDouble scale, const Array<OneD, NekDouble> &inarray,
1279  Array<OneD, NekDouble> &outarray);
1280 
1281  virtual void v_PhysGalerkinProjection1DScaled(
1282  const NekDouble scale,
1283  const Array<OneD, NekDouble> &inarray,
1284  Array<OneD, NekDouble> &outarray);
1285 
1286  void ExtractFileBCs(const std::string &fileName,
1287  const std::string &varName,
1288  const boost::shared_ptr<ExpList> locExpList);
1289 
1290  // Utility function for a common case of retrieving a
1291  // BoundaryCondition from a boundary condition collection.
1294  GetBoundaryCondition(const SpatialDomains::
1295  BoundaryConditionCollection& collection,
1296  unsigned int index, const std::string& variable);
1297 
1298 
1299  private:
1301 
1304 
1305  virtual void v_EvaluateBoundaryConditions(
1306  const NekDouble time = 0.0,
1307  const std::string varName = "",
1310 
1311  virtual map<int, RobinBCInfoSharedPtr> v_GetRobinBCInfo(void);
1312 
1313 
1314  virtual void v_GetPeriodicEntities(
1315  PeriodicMap &periodicVerts,
1316  PeriodicMap &periodicEdges,
1317  PeriodicMap &periodicFaces);
1318 
1319  // Homogeneous direction wrapper functions.
1321  {
1322  ASSERTL0(false,
1323  "This method is not defined or valid for this class type");
1325  }
1326 
1327  // wrapper function to set viscosity for Homo1D expansion
1329  {
1330  ASSERTL0(false,
1331  "This method is not defined or valid for this class type");
1332  }
1333 
1334 
1335  virtual boost::shared_ptr<ExpList> &v_GetPlane(int n);
1336  };
1337 
1338 
1339  /// Shared pointer to an ExpList object.
1340  typedef boost::shared_ptr<ExpList> ExpListSharedPtr;
1341  /// An empty ExpList object.
1343  static ExpListSharedPtr NullExpListSharedPtr;
1344 
1345  // Inline routines follow.
1346 
1347  /**
1348  * Returns the total number of local degrees of freedom
1349  * \f$N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_m\f$.
1350  */
1351  inline int ExpList::GetNcoeffs() const
1352  {
1353  return m_ncoeffs;
1354  }
1355 
1356  inline int ExpList::GetNcoeffs(const int eid) const
1357  {
1358  return (*m_exp)[eid]->GetNcoeffs();
1359  }
1360 
1361  /**
1362  * Evaulates the maximum number of modes in the elemental basis
1363  * order over all elements
1364  */
1366  {
1367  unsigned int i;
1368  int returnval = 0;
1369 
1370  for(i= 0; i < (*m_exp).size(); ++i)
1371  {
1372  returnval = max(returnval,
1373  (*m_exp)[i]->EvalBasisNumModesMax());
1374  }
1375 
1376  return returnval;
1377  }
1378 
1379  /**
1380  *
1381  */
1383  const
1384  {
1385  unsigned int i;
1386  Array<OneD,int> returnval((*m_exp).size(),0);
1387 
1388  for(i= 0; i < (*m_exp).size(); ++i)
1389  {
1390  returnval[i]
1391  = max(returnval[i],(*m_exp)[i]->EvalBasisNumModesMax());
1392  }
1393 
1394  return returnval;
1395  }
1396 
1397 
1398  /**
1399  *
1400  */
1401  inline int ExpList::GetTotPoints() const
1402  {
1403  return m_npoints;
1404  }
1405 
1406  inline int ExpList::GetTotPoints(const int eid) const
1407  {
1408  return (*m_exp)[eid]->GetTotPoints();
1409  }
1410 
1411 
1412  inline int ExpList::Get1DScaledTotPoints(const NekDouble scale) const
1413  {
1414  int returnval = 0;
1415  int cnt;
1416  int nbase = (*m_exp)[0]->GetNumBases();
1417 
1418  for(int i = 0; i < (*m_exp).size(); ++i)
1419  {
1420  cnt = 1;
1421  for(int j = 0; j < nbase; ++j)
1422  {
1423  cnt *= (int)(scale*((*m_exp)[i]->GetNumPoints(j)));
1424  }
1425  returnval += cnt;
1426  }
1427  return returnval;
1428  }
1429 
1430  /**
1431  *
1432  */
1433  inline int ExpList::GetNpoints() const
1434  {
1435  return m_npoints;
1436  }
1437 
1438 
1439  /**
1440  *
1441  */
1442  inline void ExpList::SetWaveSpace(const bool wavespace)
1443  {
1444  m_WaveSpace = wavespace;
1445  }
1446 
1447  /**
1448  *
1449  */
1450  inline bool ExpList::GetWaveSpace() const
1451  {
1452  return m_WaveSpace;
1453  }
1454 
1455  /// Set the \a i th value of\a m_phys to value \a val
1456  inline void ExpList::SetPhys(int i, NekDouble val)
1457  {
1458  m_phys[i] = val;
1459  }
1460 
1461 
1462  /**
1463  * This function fills the array \f$\boldsymbol{u}_l\f$, the evaluation
1464  * of the expansion at the quadrature points (implemented as #m_phys),
1465  * with the values of the array \a inarray.
1466  *
1467  * @param inarray The array containing the values where
1468  * #m_phys should be filled with.
1469  */
1470  inline void ExpList::SetPhys(
1471  const Array<OneD, const NekDouble> &inarray)
1472  {
1473  ASSERTL0(inarray.num_elements() == m_npoints,
1474  "Input array does not have correct number of elements.");
1475 
1476  Vmath::Vcopy(m_npoints,&inarray[0],1,&m_phys[0],1);
1477  m_physState = true;
1478  }
1479 
1480 
1482  {
1483  m_phys = inarray;
1484  }
1485 
1486 
1487  /**
1488  * @param physState \a true (=filled) or \a false (=not filled).
1489  */
1490  inline void ExpList::SetPhysState(const bool physState)
1491  {
1492  m_physState = physState;
1493  }
1494 
1495 
1496  /**
1497  * @return physState \a true (=filled) or \a false (=not filled).
1498  */
1499  inline bool ExpList::GetPhysState() const
1500  {
1501  return m_physState;
1502  }
1503 
1504  /**
1505  *
1506  */
1508  const Array<OneD, const NekDouble> &inarray,
1509  Array<OneD, NekDouble> &outarray,
1510  CoeffState coeffstate)
1511  {
1512  v_IProductWRTBase(inarray,outarray, coeffstate);
1513  }
1514 
1515  /**
1516  *
1517  */
1519  const Array<OneD, const NekDouble> &inarray,
1520  Array<OneD, NekDouble> &outarray)
1521  {
1522  v_IProductWRTBase_IterPerExp(inarray,outarray);
1523  }
1524 
1525  /**
1526  *
1527  */
1528  inline void ExpList::FwdTrans(
1529  const Array<OneD, const NekDouble> &inarray,
1530  Array<OneD, NekDouble> &outarray,
1531  CoeffState coeffstate)
1532  {
1533  v_FwdTrans(inarray,outarray,coeffstate);
1534  }
1535 
1536  /**
1537  *
1538  */
1540  const Array<OneD, const NekDouble> &inarray,
1541  Array<OneD,NekDouble> &outarray)
1542  {
1543  v_FwdTrans_IterPerExp(inarray,outarray);
1544  }
1545 
1546  /**
1547  *
1548  */
1550  {
1551  v_SmoothField(field);
1552  }
1553 
1554  /**
1555  *
1556  */
1557  inline void ExpList::BwdTrans (
1558  const Array<OneD, const NekDouble> &inarray,
1559  Array<OneD, NekDouble> &outarray,
1560  CoeffState coeffstate)
1561  {
1562  v_BwdTrans(inarray,outarray,coeffstate);
1563  }
1564 
1565  /**
1566  *
1567  */
1569  const Array<OneD, const NekDouble> &inarray,
1570  Array<OneD, NekDouble> &outarray)
1571  {
1572  v_BwdTrans_IterPerExp(inarray,outarray);
1573  }
1574 
1575 
1576  /**
1577  *
1578  */
1580  const Array<OneD,const NekDouble> &inarray,
1581  Array<OneD, NekDouble> &outarray,
1582  CoeffState coeffstate)
1583  {
1584  v_MultiplyByInvMassMatrix(inarray,outarray,coeffstate);
1585  }
1586 
1587  /**
1588  *
1589  */
1590  inline void ExpList::HelmSolve(
1591  const Array<OneD, const NekDouble> &inarray,
1592  Array<OneD, NekDouble> &outarray,
1593  const FlagList &flags,
1594  const StdRegions::ConstFactorMap &factors,
1595  const StdRegions::VarCoeffMap &varcoeff,
1596  const Array<OneD, const NekDouble> &dirForcing)
1597  {
1598  v_HelmSolve(inarray, outarray, flags, factors, varcoeff, dirForcing);
1599  }
1600 
1601 
1602  /**
1603  *
1604  */
1606  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1607  const Array<OneD, const NekDouble> &inarray,
1608  Array<OneD, NekDouble> &outarray,
1609  const NekDouble lambda,
1610  CoeffState coeffstate,
1611  const Array<OneD, const NekDouble>& dirForcing)
1612  {
1613  v_LinearAdvectionDiffusionReactionSolve(velocity,inarray, outarray,
1614  lambda, coeffstate,dirForcing);
1615  }
1616 
1618  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1619  const Array<OneD, const NekDouble> &inarray,
1620  Array<OneD, NekDouble> &outarray,
1621  const NekDouble lambda,
1622  CoeffState coeffstate,
1623  const Array<OneD, const NekDouble>& dirForcing)
1624  {
1625  v_LinearAdvectionReactionSolve(velocity,inarray, outarray,
1626  lambda, coeffstate,dirForcing);
1627  }
1628 
1629  /**
1630  *
1631  */
1633  Array<OneD, NekDouble> &coord_1,
1634  Array<OneD, NekDouble> &coord_2)
1635 
1636  {
1637  v_GetCoords(coord_0,coord_1,coord_2);
1638  }
1639 
1640  /**
1641  *
1642  */
1644  Array<OneD, NekDouble> &out_d0,
1645  Array<OneD, NekDouble> &out_d1,
1646  Array<OneD, NekDouble> &out_d2)
1647  {
1648  v_PhysDeriv(inarray,out_d0,out_d1,out_d2);
1649  }
1650 
1651  /**
1652  *
1653  */
1654  inline void ExpList::PhysDeriv(
1655  const int dir,
1656  const Array<OneD, const NekDouble> &inarray,
1657  Array<OneD, NekDouble> &out_d)
1658  {
1659  v_PhysDeriv(dir,inarray,out_d);
1660  }
1661 
1662  inline void ExpList::PhysDeriv(
1663  Direction edir,
1664  const Array<OneD, const NekDouble> &inarray,
1665  Array<OneD, NekDouble> &out_d)
1666  {
1667  v_PhysDeriv(edir, inarray,out_d);
1668  }
1669 
1670  /**
1671  *
1672  */
1674  const Array<OneD, const NekDouble> &inarray,
1675  Array<OneD, NekDouble> &outarray,
1676  CoeffState coeffstate,
1677  bool Shuff,
1678  bool UnShuff)
1679  {
1680  v_HomogeneousFwdTrans(inarray,outarray,coeffstate,Shuff,UnShuff);
1681  }
1682 
1683  /**
1684  *
1685  */
1687  const Array<OneD, const NekDouble> &inarray,
1688  Array<OneD, NekDouble> &outarray,
1689  CoeffState coeffstate,
1690  bool Shuff,
1691  bool UnShuff)
1692  {
1693  v_HomogeneousBwdTrans(inarray,outarray,coeffstate,Shuff,UnShuff);
1694  }
1695 
1696  /**
1697  *
1698  */
1700  const Array<OneD, NekDouble> &inarray1,
1701  const Array<OneD, NekDouble> &inarray2,
1702  Array<OneD, NekDouble> &outarray,
1703  CoeffState coeffstate)
1704  {
1705  v_DealiasedProd(inarray1,inarray2,outarray,coeffstate);
1706  }
1707 
1708  /**
1709  *
1710  */
1712  Array<OneD, NekDouble> &BndVals,
1713  const Array<OneD, NekDouble> &TotField,
1714  int BndID)
1715  {
1716  v_GetBCValues(BndVals,TotField,BndID);
1717  }
1718 
1719  /**
1720  *
1721  */
1725  Array<OneD, NekDouble> &outarray,
1726  int BndID)
1727  {
1728  v_NormVectorIProductWRTBase(V1,V2,outarray,BndID);
1729  }
1730 
1731  /**
1732  * @param eid The index of the element to be checked.
1733  * @return The dimension of the coordinates of the specific element.
1734  */
1735  inline int ExpList::GetCoordim(int eid)
1736  {
1737  ASSERTL2(eid <= (*m_exp).size(),
1738  "eid is larger than number of elements");
1739  return (*m_exp)[eid]->GetCoordim();
1740  }
1741 
1742  /**
1743  * @param i The index of m_coeffs to be set
1744  * @param val The value which m_coeffs[i] is to be set to.
1745  */
1746  inline void ExpList::SetCoeff(int i, NekDouble val)
1747  {
1748  m_coeffs[i] = val;
1749  }
1750 
1751 
1752  /**
1753  * @param i The index of #m_coeffs to be set.
1754  * @param val The value which #m_coeffs[i] is to be set to.
1755  */
1756  inline void ExpList::SetCoeffs(int i, NekDouble val)
1757  {
1758  m_coeffs[i] = val;
1759  }
1760 
1761 
1763  {
1764  m_coeffs = inarray;
1765  }
1766 
1767  /**
1768  * As the function returns a constant reference to a
1769  * <em>const Array</em>, it is not possible to modify the
1770  * underlying data of the array #m_coeffs. In order to do
1771  * so, use the function #UpdateCoeffs instead.
1772  *
1773  * @return (A constant reference to) the array #m_coeffs.
1774  */
1776  {
1777  return m_coeffs;
1778  }
1779 
1781  Array<OneD,NekDouble>& outarray)
1782  {
1783  v_ImposeDirichletConditions(outarray);
1784  }
1785 
1787  {
1789  }
1790 
1791  inline void ExpList::LocalToGlobal(void)
1792  {
1793  v_LocalToGlobal();
1794  }
1795 
1796  inline void ExpList::GlobalToLocal(void)
1797  {
1798  v_GlobalToLocal();
1799  }
1800 
1801 
1802  /**
1803  * @param i The index of #m_coeffs to be returned
1804  * @return The NekDouble held in #m_coeffs[i].
1805  */
1807  {
1808  return m_coeffs[i];
1809  }
1810 
1811  /**
1812  * @param i The index of #m_coeffs to be returned
1813  * @return The NekDouble held in #m_coeffs[i].
1814  */
1816  {
1817  return m_coeffs[i];
1818  }
1819 
1820  /**
1821  * As the function returns a constant reference to a
1822  * <em>const Array</em> it is not possible to modify the
1823  * underlying data of the array #m_phys. In order to do
1824  * so, use the function #UpdatePhys instead.
1825  *
1826  * @return (A constant reference to) the array #m_phys.
1827  */
1829  {
1830  return m_phys;
1831  }
1832 
1833  /**
1834  * @return \f$N_{\mathrm{el}}\f$, the number of elements in the
1835  * expansion.
1836  */
1837  inline int ExpList::GetExpSize(void)
1838  {
1839  return (*m_exp).size();
1840  }
1841 
1842 
1843  /**
1844  * @param n The index of the element concerned.
1845  *
1846  * @return (A shared pointer to) the local expansion of the
1847  * \f$n^{\mathrm{th}}\f$ element.
1848  */
1850  {
1851  return (*m_exp)[n];
1852  }
1853 
1854  /**
1855  * @return (A const shared pointer to) the local expansion vector #m_exp
1856  */
1857  inline const boost::shared_ptr<LocalRegions::ExpansionVector>
1858  ExpList::GetExp(void) const
1859  {
1860  return m_exp;
1861  }
1862 
1863 
1864  /**
1865  *
1866  */
1867  inline int ExpList::GetCoeff_Offset(int n) const
1868  {
1869  return m_coeff_offset[n];
1870  }
1871 
1872  /**
1873  *
1874  */
1875  inline int ExpList::GetPhys_Offset(int n) const
1876  {
1877  return m_phys_offset[n];
1878  }
1879 
1880  /**
1881  *
1882  */
1883  inline int ExpList::GetOffset_Elmt_Id(int n) const
1884  {
1885  return m_offset_elmt_id[n];
1886  }
1887 
1888  /**
1889  * If one wants to get hold of the underlying data without modifying
1890  * them, rather use the function #GetCoeffs instead.
1891  *
1892  * @return (A reference to) the array #m_coeffs.
1893  */
1895  {
1896  return m_coeffs;
1897  }
1898 
1899  /**
1900  * If one wants to get hold of the underlying data without modifying
1901  * them, rather use the function #GetPhys instead.
1902  *
1903  * @return (A reference to) the array #m_phys.
1904  */
1906  {
1907  m_physState = true;
1908  return m_phys;
1909  }
1910 
1911 
1912  // functions associated with DisContField
1915  {
1916  return v_GetBndCondExpansions();
1917  }
1918 
1919  inline boost::shared_ptr<ExpList> &ExpList::UpdateBndCondExpansion(int i)
1920  {
1921  return v_UpdateBndCondExpansion(i);
1922  }
1923 
1924  inline void ExpList::Upwind(
1925  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
1926  const Array<OneD, const NekDouble> &Fwd,
1927  const Array<OneD, const NekDouble> &Bwd,
1928  Array<OneD, NekDouble> &Upwind)
1929  {
1930  v_Upwind(Vec, Fwd, Bwd, Upwind);
1931  }
1932 
1933  inline void ExpList::Upwind(
1934  const Array<OneD, const NekDouble> &Vn,
1935  const Array<OneD, const NekDouble> &Fwd,
1936  const Array<OneD, const NekDouble> &Bwd,
1937  Array<OneD, NekDouble> &Upwind)
1938  {
1939  v_Upwind(Vn, Fwd, Bwd, Upwind);
1940  }
1941 
1942  inline boost::shared_ptr<ExpList> &ExpList::GetTrace()
1943  {
1944  return v_GetTrace();
1945  }
1946 
1947  inline boost::shared_ptr<AssemblyMapDG> &ExpList::GetTraceMap()
1948  {
1949  return v_GetTraceMap();
1950  }
1951 
1953  {
1954  return v_GetTraceBndMap();
1955  }
1956 
1957  inline void ExpList::GetNormals(
1958  Array<OneD, Array<OneD, NekDouble> > &normals)
1959  {
1960  v_GetNormals(normals);
1961  }
1962 
1964  const Array<OneD, const NekDouble> &Fx,
1965  const Array<OneD, const NekDouble> &Fy,
1966  Array<OneD, NekDouble> &outarray)
1967  {
1968  v_AddTraceIntegral(Fx,Fy,outarray);
1969  }
1970 
1972  const Array<OneD, const NekDouble> &Fn,
1973  Array<OneD, NekDouble> &outarray)
1974  {
1975  v_AddTraceIntegral(Fn,outarray);
1976  }
1977 
1979  const Array<OneD, const NekDouble> &Fwd,
1980  const Array<OneD, const NekDouble> &Bwd,
1981  Array<OneD, NekDouble> &outarray)
1982  {
1983  v_AddFwdBwdTraceIntegral(Fwd,Bwd,outarray);
1984  }
1985 
1987  Array<OneD,NekDouble> &Fwd,
1988  Array<OneD,NekDouble> &Bwd)
1989  {
1990  v_GetFwdBwdTracePhys(Fwd,Bwd);
1991  }
1992 
1994  const Array<OneD,const NekDouble> &field,
1995  Array<OneD,NekDouble> &Fwd,
1996  Array<OneD,NekDouble> &Bwd)
1997  {
1998  v_GetFwdBwdTracePhys(field,Fwd,Bwd);
1999  }
2000 
2002  {
2003  v_ExtractTracePhys(outarray);
2004  }
2005 
2006 
2008  const Array<OneD, const NekDouble> &inarray,
2009  Array<OneD,NekDouble> &outarray)
2010  {
2011  v_ExtractTracePhys(inarray,outarray);
2012  }
2013 
2016  {
2017  return v_GetBndConditions();
2018  }
2019 
2020 
2023  {
2024  return v_UpdateBndConditions();
2025  }
2026 
2028  const NekDouble time,
2029  const std::string varName,
2030  const NekDouble x2_in,
2031  const NekDouble x3_in)
2032  {
2033  v_EvaluateBoundaryConditions(time, varName, x2_in, x3_in);
2034  }
2035 
2036  // Routines for continous matrix solution
2037  /**
2038  * This operation is equivalent to the evaluation of
2039  * \f$\underline{\boldsymbol{M}}^e\boldsymbol{\hat{u}}_l\f$, that is,
2040  * \f[ \left[
2041  * \begin{array}{cccc}
2042  * \boldsymbol{M}^1 & 0 & \hspace{3mm}0 \hspace{3mm}& 0 \\
2043  * 0 & \boldsymbol{M}^2 & 0 & 0 \\
2044  * 0 & 0 & \ddots & 0 \\
2045  * 0 & 0 & 0 & \boldsymbol{M}^{N_{\mathrm{el}}} \end{array} \right]
2046  *\left [ \begin{array}{c}
2047  * \boldsymbol{\hat{u}}^{1} \\
2048  * \boldsymbol{\hat{u}}^{2} \\
2049  * \vdots \\
2050  * \boldsymbol{\hat{u}}^{{{N_{\mathrm{el}}}}} \end{array} \right ]\f]
2051  * where \f$\boldsymbol{M}^e\f$ are the local matrices of type
2052  * specified by the key \a mkey. The decoupling of the local matrices
2053  * allows for a local evaluation of the operation. However, rather than
2054  * a local matrix-vector multiplication, the local operations are
2055  * evaluated as implemented in the function
2056  * StdRegions#StdExpansion#GeneralMatrixOp.
2057  *
2058  * @param mkey This key uniquely defines the type matrix
2059  * required for the operation.
2060  * @param inarray The vector \f$\boldsymbol{\hat{u}}_l\f$ of
2061  * size \f$N_{\mathrm{eof}}\f$.
2062  * @param outarray The resulting vector of size
2063  * \f$N_{\mathrm{eof}}\f$.
2064  */
2066  const GlobalMatrixKey &gkey,
2067  const Array<OneD,const NekDouble> &inarray,
2068  Array<OneD, NekDouble> &outarray,
2069  CoeffState coeffstate)
2070  {
2071  v_GeneralMatrixOp(gkey,inarray,outarray,coeffstate);
2072  }
2073 
2074 
2076  {
2078  }
2079 
2081  Array<OneD,int> &EdgeID)
2082  {
2083  v_GetBoundaryToElmtMap(ElmtID,EdgeID);
2084  }
2085 
2087 
2088  } //end of namespace
2089 } //end of namespace
2090 
2091 #endif // EXPLIST_H
2092 
virtual NekDouble v_L2(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
Definition: ExpList.cpp:1766
virtual void v_GlobalToLocal(void)
Definition: ExpList.cpp:2444
virtual void v_AddFwdBwdTraceIntegral(const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2298
virtual void v_SmoothField(Array< OneD, NekDouble > &field)
Definition: ExpList.cpp:594
virtual void v_PhysDeriv(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d0, Array< OneD, NekDouble > &out_d1, Array< OneD, NekDouble > &out_d2)
Definition: ExpList.cpp:448
void AddTraceIntegral(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:1963
virtual void v_GeneralMatrixOp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2480
const DNekScalBlkMatSharedPtr & GetBlockMatrix(const GlobalMatrixKey &gkey)
Definition: ExpList.cpp:772
void SetModifiedBasis(const bool modbasis)
Set Modified Basis for the stability analysis.
virtual void v_ExtractDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract data from raw field data into expansion list.
Definition: ExpList.cpp:2111
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:1775
void SmoothField(Array< OneD, NekDouble > &field)
Smooth a field across elements.
Definition: ExpList.h:1549
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
int GetExpIndex(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0, bool returnNearestElmt=false)
Definition: ExpList.cpp:1222
virtual void v_Reset()
Reset geometry information, metrics, matrix managers and geometry information.
Definition: ExpList.cpp:1394
boost::shared_ptr< GlobalMatrix > GenGlobalMatrix(const GlobalMatrixKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
Generates a global matrix from the given key and map.
Definition: ExpList.cpp:855
virtual void v_SetHomo1DSpecVanVisc(Array< OneD, NekDouble > visc)
Definition: ExpList.h:1328
void GeneralMatrixOp_IterPerExp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:787
virtual const Array< OneD, const boost::shared_ptr< ExpList > > & v_GetBndCondExpansions(void)
Definition: ExpList.cpp:2217
boost::shared_ptr< Transposition > TranspositionSharedPtr
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:1867
ExpList()
The default constructor.
Definition: ExpList.cpp:93
Local coefficients.
static ExpListSharedPtr NullExpListSharedPtr
Definition: ExpList.h:1343
static Array< OneD, NekDouble > NullNekDouble1DArray
void EvaluateBoundaryConditions(const NekDouble time=0.0, const std::string varName="", const NekDouble=NekConstants::kNekUnsetDouble, const NekDouble=NekConstants::kNekUnsetDouble)
Definition: ExpList.h:2027
std::map< int, vector< PeriodicEntity > > PeriodicMap
int Get1DScaledTotPoints(const NekDouble scale) const
Returns the total number of qudature points scaled by the factor scale on each 1D direction...
Definition: ExpList.h:1412
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:971
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:53
boost::shared_ptr< ExpList > & GetPlane(int n)
Definition: ExpList.h:859
boost::shared_ptr< BlockMatrixMap > BlockMatrixMapShPtr
A shared pointer to a BlockMatrixMap.
Definition: ExpList.h:97
LibUtilities::TranspositionSharedPtr GetTransposition(void)
This function returns the transposition class associaed with the homogeneous expansion.
Definition: ExpList.h:510
void MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1579
virtual boost::shared_ptr< AssemblyMapDG > & v_GetTraceMap()
Definition: ExpList.cpp:2261
int GetPhys_Offset(int n) const
Get the start offset position for a global list of m_phys correspoinding to element n...
Definition: ExpList.h:1875
Array< OneD, NekDouble > & UpdateCoeffs()
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:1894
virtual void v_FillBndCondFromField()
Definition: ExpList.cpp:2432
void ExtractFileBCs(const std::string &fileName, const std::string &varName, const boost::shared_ptr< ExpList > locExpList)
Definition: ExpList.cpp:1862
boost::shared_ptr< DNekMat > GenGlobalMatrixFull(const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
Definition: ExpList.cpp:992
NekDouble PhysIntegral(void)
This function integrates a function over the domain consisting of all the elements of the expansion...
Definition: ExpList.cpp:236
virtual void v_WriteTecplotZone(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1460
ExpansionType GetExpType(void)
Returns the type of the expansion.
Definition: ExpList.cpp:203
const Array< OneD, const int > & GetTraceBndMap(void)
Definition: ExpList.h:1952
boost::shared_ptr< FieldDefinitions > FieldDefinitionsSharedPtr
Definition: FieldIO.h:131
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2465
const Array< OneD, int > EvalBasisNumModesMaxPerExp(void) const
Returns the vector of the number of modes in the elemental basis order over all elements.
Definition: ExpList.h:1382
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:1858
std::vector< LibUtilities::FieldDefinitionsSharedPtr > GetFieldDefinitions()
Definition: ExpList.h:766
virtual void v_BwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1189
static BasisSharedPtr NullBasisSharedPtr
Definition: Basis.h:358
#define MULTI_REGIONS_EXPORT
boost::shared_ptr< ExpList > & UpdateBndCondExpansion(int i)
Definition: ExpList.h:1919
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1557
virtual void v_GetNormals(Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.cpp:2274
Array< OneD, const unsigned int > GetZIDs(void)
This function returns a vector containing the wave numbers in z-direction associated with the 3D homo...
Definition: ExpList.h:503
map< int, RobinBCInfoSharedPtr > GetRobinBCInfo()
Definition: ExpList.h:752
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:248
boost::shared_ptr< LibUtilities::SessionReader > GetSession()
Returns the session object.
Definition: ExpList.h:837
void ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Impose Dirichlet Boundary Conditions onto Array.
Definition: ExpList.h:1780
void GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Definition: ExpList.h:2080
void GeneralGetFieldDefinitions(std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef, int NumHomoDir=0, int NumHomoStrip=1, Array< OneD, LibUtilities::BasisSharedPtr > &HomoBasis=LibUtilities::NullBasisSharedPtr1DArray, std::vector< NekDouble > &HomoLen=LibUtilities::NullNekDoubleVector, std::vector< unsigned int > &HomoZIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoYIDs=LibUtilities::NullUnsignedIntVector)
Definition: ExpList.cpp:1937
BlockMatrixMapShPtr m_blockMat
Definition: ExpList.h:973
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:926
void LinearAdvectionDiffusionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Solve Advection Diffusion Reaction.
Definition: ExpList.h:1605
void LinearAdvectionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Solve Advection Diffusion Reaction.
Definition: ExpList.h:1617
void ApplyGeomInfo()
Apply geometry information to each expansion.
Definition: ExpList.cpp:1381
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:909
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
Definition: ExpList.cpp:1679
NekDouble Integral(const Array< OneD, const NekDouble > &inarray)
Definition: ExpList.h:480
virtual void v_NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
Definition: ExpList.cpp:2415
virtual LibUtilities::TranspositionSharedPtr v_GetTransposition(void)
Definition: ExpList.cpp:1818
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:1837
void PhysDeriv(Direction edir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
Definition: ExpList.h:1662
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:50
virtual boost::shared_ptr< ExpList > & v_UpdateBndCondExpansion(int i)
Definition: ExpList.cpp:2225
boost::shared_ptr< GlobalLinSys > GenGlobalBndLinSys(const GlobalLinSysKey &mkey, const AssemblyMapSharedPtr &locToGloMap)
Generate a GlobalLinSys from information provided by the key "mkey" and the mapping provided in LocTo...
Definition: ExpList.cpp:1152
void ExtractElmtDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract the data in fielddata into the coeffs using the basic ExpList Elemental expansions rather tha...
void PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function Galerkin projects the physical space points in inarray to outarray where inarray is ass...
Definition: ExpList.h:547
void MultiplyByElmtInvMass(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function elementally mulplies the coefficient space of Sin my the elemental inverse of the mass ...
Definition: ExpList.cpp:523
virtual void v_WriteTecplotHeader(std::ostream &outfile, std::string var="")
Definition: ExpList.cpp:1421
virtual NekDouble v_GetHomoLen(void)
Definition: ExpList.cpp:1826
void WriteVtkFooter(std::ostream &outfile)
Definition: ExpList.cpp:1673
boost::shared_ptr< GlobalOptParam > GlobalOptParamSharedPtr
Pointer to a GlobalOptParam object.
void SetCoeffs(int i, NekDouble val)
Set the i th coefficiient in m_coeffs to value val.
Definition: ExpList.h:1756
virtual void v_GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
Definition: ExpList.cpp:2500
const Array< OneD, const boost::shared_ptr< ExpList > > & GetBndCondExpansions()
Definition: ExpList.h:1914
void WriteTecplotField(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:379
map< GlobalMatrixKey, DNekScalBlkMatSharedPtr > BlockMatrixMap
A map between global matrix keys and their associated block matrices.
Definition: ExpList.h:95
boost::shared_ptr< ExpList > GetSharedThisPtr()
Returns a shared pointer to the current object.
Definition: ExpList.h:831
void Upwind(const Array< OneD, const Array< OneD, NekDouble > > &Vec, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
Definition: ExpList.h:1924
virtual Array< OneD, const NekDouble > v_HomogeneousEnergy(void)
Definition: ExpList.cpp:1810
virtual const Array< OneD, const int > & v_GetTraceBndMap()
Definition: ExpList.cpp:2269
virtual void v_WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var)
Definition: ExpList.cpp:1690
static std::vector< NekDouble > NullNekDoubleVector
Definition: FieldIO.h:61
void GeneralMatrixOp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
This function calculates the result of the multiplication of a matrix of type specified by mkey with ...
Definition: ExpList.h:2065
void NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
Definition: ExpList.h:1722
virtual void v_GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
Definition: ExpList.cpp:2407
void GetNormals(Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.h:1957
void Reset()
Reset geometry information and reset matrices.
Definition: ExpList.h:361
boost::shared_ptr< LibUtilities::Comm > GetComm()
Returns the comm object.
Definition: ExpList.h:843
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:958
LibUtilities::BasisSharedPtr GetHomogeneousBasis(void)
Definition: ExpList.h:854
const NekOptimize::GlobalOptParamSharedPtr & GetGlobalOptParam(void)
Definition: ExpList.h:747
int GetNumElmts(void)
This function returns the number of elements in the expansion which may be different for a homogeoeno...
Definition: ExpList.h:562
Base class for all multi-elemental spectral/hp expansions.
Definition: ExpList.h:101
boost::shared_ptr< Comm > CommSharedPtr
Pointer to a Communicator object.
Definition: Comm.h:53
boost::shared_ptr< ExpList > & GetTrace()
Definition: ExpList.h:1942
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 . ...
Definition: ExpList.h:1632
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2424
virtual LibUtilities::BasisSharedPtr v_GetHomogeneousBasis(void)
Definition: ExpList.h:1320
virtual void v_Upwind(const Array< OneD, const Array< OneD, NekDouble > > &Vec, const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &Upwind)
Definition: ExpList.cpp:2233
NekDouble Linf(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
This function calculates the error of the global spectral/hp element approximation.
Definition: ExpList.cpp:1727
boost::shared_ptr< LocalRegions::ExpansionVector > m_exp
The list of local expansions.
Definition: ExpList.h:947
int GetTotPoints(void) const
Returns the total number of quadrature points m_npoints .
Definition: ExpList.h:1401
void BwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function elementally evaluates the backward transformation of the global spectral/hp element exp...
Definition: ExpList.h:1568
virtual void v_FwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:562
virtual ~ExpList()
The default destructor.
Definition: ExpList.cpp:216
void IProductWRTDerivBase(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function calculates the inner product of a function with respect to the derivative (in directio...
Definition: ExpList.cpp:341
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:225
boost::shared_ptr< GlobalLinSys > GenGlobalLinSys(const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
This operation constructs the global linear system of type mkey.
Definition: ExpList.cpp:1133
void WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var="v")
Definition: ExpList.h:404
void AddFwdBwdTraceIntegral(const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:1978
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:961
void FillBndCondFromField(void)
Fill Bnd Condition expansion from the values stored in expansion.
Definition: ExpList.h:1786
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
Definition: ExpList.h:1340
void IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function calculates the inner product of a function with respect to all {local} expansion modes...
Definition: ExpList.h:1518
NekDouble GetHomoLen(void)
This function returns the Width of homogeneous direction associaed with the homogeneous expansion...
Definition: ExpList.h:517
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Definition: ExpList.cpp:1796
virtual std::vector< LibUtilities::FieldDefinitionsSharedPtr > v_GetFieldDefinitions(void)
Definition: ExpList.cpp:2047
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:935
virtual void v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1554
static PeriodicMap NullPeriodicMap
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:887
void WriteTecplotConnectivity(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:385
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:969
boost::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:74
void ExtractTracePhys(Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:2001
static std::vector< unsigned int > NullUnsignedIntVector
Definition: FieldIO.h:63
virtual void v_HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.cpp:2381
virtual void v_LocalToGlobal(void)
Definition: ExpList.cpp:2438
void SetWaveSpace(const bool wavespace)
Sets the wave space to the one of the possible configuration true or false.
Definition: ExpList.h:1442
void SetHomo1DSpecVanVisc(Array< OneD, NekDouble > visc)
This function sets the Spectral Vanishing Viscosity in homogeneous1D expansion.
Definition: ExpList.h:494
SpatialDomains::MeshGraphSharedPtr m_graph
Mesh associated with this expansion list.
Definition: ExpList.h:883
virtual void v_EvaluateBoundaryConditions(const NekDouble time=0.0, const std::string varName="", const NekDouble x2_in=NekConstants::kNekUnsetDouble, const NekDouble x3_in=NekConstants::kNekUnsetDouble)
Definition: ExpList.cpp:2600
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:880
SpatialDomains::MeshGraphSharedPtr GetGraph()
Definition: ExpList.h:848
double NekDouble
void WriteTecplotHeader(std::ostream &outfile, std::string var="")
Definition: ExpList.h:366
const DNekScalBlkMatSharedPtr GenBlockMatrix(const GlobalMatrixKey &gkey)
This function assembles the block diagonal matrix of local matrices of the type mtype.
Definition: ExpList.cpp:636
Defines a list of flags.
std::map< int, BoundaryConditionMapShPtr > BoundaryConditionCollection
Definition: Conditions.h:213
void HelmSolve(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const FlagList &flags, const StdRegions::ConstFactorMap &factors, const StdRegions::VarCoeffMap &varcoeff=StdRegions::NullVarCoeffMap, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Solve helmholtz problem.
Definition: ExpList.h:1590
void FwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function elementally evaluates the forward transformation of a function onto the global spectra...
Definition: ExpList.h:1539
void GetPeriodicEntities(PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces=NullPeriodicMap)
Definition: ExpList.h:757
virtual Array< OneD, SpatialDomains::BoundaryConditionShPtr > & v_UpdateBndConditions()
Definition: ExpList.cpp:2590
static const NekDouble kNekUnsetDouble
Describe a linear system.
void FwdTrans_BndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:572
virtual boost::shared_ptr< ExpList > & v_GetTrace()
Definition: ExpList.cpp:2253
void HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.h:1673
virtual void v_DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.cpp:2401
virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
Definition: ExpList.cpp:2186
virtual boost::shared_ptr< ExpList > & v_GetPlane(int n)
Definition: ExpList.cpp:2645
Describes a matrix with ordering defined by a local to global map.
static ExpList NullExpList
An empty ExpList object.
Definition: ExpList.h:1342
void LocalToGlobal(void)
Put the coefficients into global ordering using m_coeffs.
Definition: ExpList.h:1791
const Array< OneD, const NekDouble > & GetPhys() const
This function returns (a reference to) the array (implemented as m_phys) containing the function ev...
Definition: ExpList.h:1828
int GetNpoints(void) const
Returns the total number of quadrature points m_npoints .
Definition: ExpList.h:1433
Array< OneD, SpatialDomains::BoundaryConditionShPtr > & UpdateBndConditions()
Definition: ExpList.h:2022
virtual void v_GetPeriodicEntities(PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces)
Definition: ExpList.cpp:2622
virtual void v_PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1851
void GlobalEigenSystem(const boost::shared_ptr< DNekMat > &Gmat, Array< OneD, NekDouble > &EigValsReal, Array< OneD, NekDouble > &EigValsImag, Array< OneD, NekDouble > &EigVecs=NullNekDouble1DArray)
Array< OneD, NekDouble > & UpdatePhys()
This function returns (a reference to) the array (implemented as m_phys) containing the function ev...
Definition: ExpList.h:1905
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: ExpList.h:877
virtual void v_AddTraceIntegral(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2281
void GetFieldDefinitions(std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef)
Definition: ExpList.h:772
void ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
Extract the data from fromField using fromExpList the coeffs using the basic ExpList Elemental expans...
Definition: ExpList.cpp:2098
Array< OneD, const NekDouble > HomogeneousEnergy(void)
This function calculates the energy associated with each one of the modesof a 3D homogeneous nD expan...
Definition: ExpList.h:487
virtual void v_LinearAdvectionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Definition: ExpList.cpp:2369
virtual Array< OneD, const unsigned int > v_GetZIDs(void)
Definition: ExpList.cpp:1834
void DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1699
bool GetPhysState(void) const
This function indicates whether the array of physical values (implemented as m_phys) is filled or no...
Definition: ExpList.h:1499
static const Array< OneD, ExpListSharedPtr > NullExpListSharedPtrArray
Definition: ExpList.h:2086
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1528
virtual map< int, RobinBCInfoSharedPtr > v_GetRobinBCInfo(void)
Definition: ExpList.cpp:2612
static Array< OneD, BasisSharedPtr > NullBasisSharedPtr1DArray
Definition: Basis.h:359
#define ASSERTL2(condition, msg)
Assert Level 2 – Debugging which is used FULLDEBUG compilation mode. This level assert is designed t...
Definition: ErrorUtil.hpp:213
NekDouble GetCoeff(int i)
Get the i th value (coefficient) of m_coeffs.
Definition: ExpList.h:1806
NekDouble H1(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
Calculates the error of the global spectral/hp element approximation.
Definition: ExpList.cpp:1918
void SetCoeff(int i, NekDouble val)
Set the i th coefficiient in m_coeffs to value val.
Definition: ExpList.h:1746
int GetOffset_Elmt_Id(int n) const
Get the element id associated with the n th consecutive block of data in m_phys and m_coeffs...
Definition: ExpList.h:1883
virtual void v_IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:314
Array< OneD, const unsigned int > GetYIDs(void)
This function returns a vector containing the wave numbers in y-direction associated with the 3D homo...
Definition: ExpList.h:526
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
void IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1507
void SetPhysState(const bool physState)
This function manually sets whether the array of physical values (implemented as m_phys) is filled o...
Definition: ExpList.h:1490
void WriteVtkPieceFooter(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1684
void AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
Append the element data listed in elements fielddef->m_ElementIDs onto fielddata. ...
Definition: ExpList.h:781
int GetNcoeffs(void) const
Returns the total number of local degrees of freedom .
Definition: ExpList.h:1351
virtual void v_LinearAdvectionDiffusionReactionSolve(const Array< OneD, Array< OneD, NekDouble > > &velocity, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const NekDouble lambda, CoeffState coeffstate=eLocal, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Definition: ExpList.cpp:2357
virtual void v_HelmSolve(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const FlagList &flags, const StdRegions::ConstFactorMap &factors, const StdRegions::VarCoeffMap &varcoeff, const Array< OneD, const NekDouble > &dirForcing)
Definition: ExpList.cpp:2346
std::vector< int > m_coll_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:955
virtual void v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
Definition: ExpList.cpp:2061
static SpatialDomains::BoundaryConditionShPtr GetBoundaryCondition(const SpatialDomains::BoundaryConditionCollection &collection, unsigned int index, const std::string &variable)
Definition: ExpList.cpp:2631
boost::shared_ptr< AssemblyMapDG > & GetTraceMap(void)
Definition: ExpList.h:1947
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:2723
const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & GetBndConditions()
Definition: ExpList.h:2015
void WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip=0)
Definition: ExpList.h:394
void MultiplyByBlockMatrix(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:284
void SetCoeffsArray(Array< OneD, NekDouble > &inarray)
Set the m_coeffs array to inarray.
Definition: ExpList.h:1762
virtual void v_MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2337
int EvalBasisNumModesMax(void) const
Evaulates the maximum number of modes in the elemental basis order over all elements.
Definition: ExpList.h:1365
bool GetWaveSpace(void) const
This function returns the third direction expansion condition, which can be in wave space (coefficien...
Definition: ExpList.h:1450
void PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
This function interpolates the physical space points in inarray to outarray using the same points def...
Definition: ExpList.h:535
virtual void v_GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Definition: ExpList.cpp:2561
virtual Array< OneD, const unsigned int > v_GetYIDs(void)
Definition: ExpList.cpp:1842
virtual void v_SetUpPhysNormals()
Definition: ExpList.cpp:2553
void GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
Definition: ExpList.h:1711
int GetCoordim(int eid)
This function returns the dimension of the coordinates of the element eid.
Definition: ExpList.h:1735
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:206
boost::shared_ptr< Basis > BasisSharedPtr
void AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, Array< OneD, NekDouble > &coeffs)
Append the data in coeffs listed in elements fielddef->m_ElementIDs onto fielddata.
Definition: ExpList.h:791
virtual const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & v_GetBndConditions()
Definition: ExpList.cpp:2579
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:211
void SetPhysArray(Array< OneD, NekDouble > &inarray)
Sets the array m_phys.
Definition: ExpList.h:1481
void WriteTecplotZone(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:372
std::vector< int > m_coll_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:952
void WriteVtkHeader(std::ostream &outfile)
Definition: ExpList.cpp:1665
virtual void v_ExtractTracePhys(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2323
virtual void v_PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1857
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:432
void GlobalToLocal(void)
Put the coefficients into local ordering and place in m_coeffs.
Definition: ExpList.h:1796
NekDouble L2(const Array< OneD, const NekDouble > &inarray, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
This function calculates the error with respect to soln of the global spectral/hp element approximat...
Definition: ExpList.h:467
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1038
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2451
void GetFwdBwdTracePhys(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Definition: ExpList.h:1986
virtual void v_GetFwdBwdTracePhys(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Definition: ExpList.cpp:2307
void ExtractDataToCoeffs(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata, std::string &field, Array< OneD, NekDouble > &coeffs)
Extract the data in fielddata into the coeffs.
Definition: ExpList.cpp:2089
void HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.h:1686
void SetPhys(int i, NekDouble val)
Set the i th value of m_phys to value val.
Definition: ExpList.h:1456
virtual void v_WriteTecplotField(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1631
virtual int v_GetNumElmts(void)
Definition: ExpList.h:1026
Collections::CollectionVector m_collections
Definition: ExpList.h:949
virtual void v_ReadGlobalOptimizationParameters()
Definition: ExpList.cpp:2570
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:226
virtual void v_HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.cpp:2391
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2458