Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
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 std::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  /// Constructor copying only elements defined in eIds.
118  const ExpList &in,
119  const std::vector<unsigned int> &eIDs,
120  const bool DeclareCoeffPhysArrays = true);
121 
122  /// The copy constructor.
124  const ExpList &in,
125  const bool DeclareCoeffPhysArrays = true);
126 
127  /// The default destructor.
128  MULTI_REGIONS_EXPORT virtual ~ExpList();
129 
130  /// Returns the total number of local degrees of freedom
131  /// \f$N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_m\f$.
132  inline int GetNcoeffs(void) const;
133 
134  /// Returns the total number of local degrees of freedom
135  /// for element eid
136  MULTI_REGIONS_EXPORT int GetNcoeffs(const int eid) const;
137 
138  /// Returns the type of the expansion
140 
141  /// Returns the type of the expansion
143 
144  /// Evaulates the maximum number of modes in the elemental basis
145  /// order over all elements
146  inline int EvalBasisNumModesMax(void) const;
147 
148  /// Returns the vector of the number of modes in the elemental
149  /// basis order over all elements.
151  EvalBasisNumModesMaxPerExp(void) const;
152 
153  /// Returns the total number of quadrature points #m_npoints
154  /// \f$=Q_{\mathrm{tot}}\f$.
155  inline int GetTotPoints(void) const;
156 
157  /// Returns the total number of quadrature points for eid's element
158  /// \f$=Q_{\mathrm{tot}}\f$.
159  inline int GetTotPoints(const int eid) const;
160 
161  /// Returns the total number of quadrature points #m_npoints
162  /// \f$=Q_{\mathrm{tot}}\f$.
163  inline int GetNpoints(void) const;
164 
165 
166  /// Returns the total number of qudature points scaled by
167  /// the factor scale on each 1D direction
168  inline int Get1DScaledTotPoints(const NekDouble scale) const;
169 
170  /// Sets the wave space to the one of the possible configuration
171  /// true or false
172  inline void SetWaveSpace(const bool wavespace);
173 
174 
175  ///Set Modified Basis for the stability analysis
176  inline void SetModifiedBasis(const bool modbasis);
177 
178  /// Set the \a i th value of \a m_phys to value \a val
179  inline void SetPhys(int i, NekDouble val);
180 
181  /// This function returns the third direction expansion condition,
182  /// which can be in wave space (coefficient) or not
183  /// It is stored in the variable m_WaveSpace.
184  inline bool GetWaveSpace(void) const;
185 
186  /// Fills the array #m_phys
187  inline void SetPhys(const Array<OneD, const NekDouble> &inarray);
188 
189  /// Sets the array #m_phys
190  inline void SetPhysArray(Array<OneD, NekDouble> &inarray);
191 
192  /// This function manually sets whether the array of physical
193  /// values \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) is
194  /// filled or not.
195  inline void SetPhysState(const bool physState);
196 
197  /// This function indicates whether the array of physical values
198  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) is filled or
199  /// not.
200  inline bool GetPhysState(void) const;
201 
202  /// This function integrates a function \f$f(\boldsymbol{x})\f$
203  /// over the domain consisting of all the elements of the expansion.
205 
206  /// This function integrates a function \f$f(\boldsymbol{x})\f$
207  /// over the domain consisting of all the elements of the expansion.
209  const Array<OneD,
210  const NekDouble> &inarray);
211 
212  /// This function calculates the inner product of a function
213  /// \f$f(\boldsymbol{x})\f$ with respect to all \emph{local}
214  /// expansion modes \f$\phi_n^e(\boldsymbol{x})\f$.
215  inline void IProductWRTBase_IterPerExp(
216  const Array<OneD, const NekDouble> &inarray,
217  Array<OneD, NekDouble> &outarray);
218 
219  ///
220  inline void IProductWRTBase(
221  const Array<OneD, const NekDouble> &inarray,
222  Array<OneD, NekDouble> &outarray,
223  CoeffState coeffstate = eLocal);
224 
225  /// This function calculates the inner product of a function
226  /// \f$f(\boldsymbol{x})\f$ with respect to the derivative (in
227  /// direction \param dir) of all \emph{local} expansion modes
228  /// \f$\phi_n^e(\boldsymbol{x})\f$.
230  (const int dir,
231  const Array<OneD, const NekDouble> &inarray,
232  Array<OneD, NekDouble> &outarray);
233 
234  /// This function calculates the inner product of a function
235  /// \f$f(\boldsymbol{x})\f$ with respect to the derivative (in
236  /// direction \param dir) of all \emph{local} expansion modes
237  /// \f$\phi_n^e(\boldsymbol{x})\f$.
239  (const Array<OneD, const Array<OneD, NekDouble> > &inarray,
240  Array<OneD, NekDouble> &outarray);
241 
242  /// This function elementally evaluates the forward transformation
243  /// of a function \f$u(\boldsymbol{x})\f$ onto the global
244  /// spectral/hp expansion.
245  inline void FwdTrans_IterPerExp (
246  const Array<OneD,
247  const NekDouble> &inarray,
248  Array<OneD,NekDouble> &outarray);
249 
250  ///
251  inline void FwdTrans(
252  const Array<OneD,
253  const NekDouble> &inarray,
254  Array<OneD, NekDouble> &outarray,
255  CoeffState coeffstate = eLocal);
256 
257  /// This function elementally mulplies the coefficient space of
258  /// Sin my the elemental inverse of the mass matrix.
260  const Array<OneD,
261  const NekDouble> &inarray,
262  Array<OneD, NekDouble> &outarray);
263 
264  ///
265  inline void MultiplyByInvMassMatrix(
266  const Array<OneD,const NekDouble> &inarray,
267  Array<OneD, NekDouble> &outarray,
268  CoeffState coeffstate = eLocal);
269 
270  /// Smooth a field across elements
271  inline void SmoothField(Array<OneD,NekDouble> &field);
272 
273  /// Solve helmholtz problem
274  inline void HelmSolve(
275  const Array<OneD, const NekDouble> &inarray,
276  Array<OneD, NekDouble> &outarray,
277  const FlagList &flags,
278  const StdRegions::ConstFactorMap &factors,
279  const StdRegions::VarCoeffMap &varcoeff =
281  const Array<OneD, const NekDouble> &dirForcing =
283  const bool PhysSpaceForcing = true);
284 
285  /// Solve Advection Diffusion Reaction
287  const Array<OneD, Array<OneD, NekDouble> > &velocity,
288  const Array<OneD, const NekDouble> &inarray,
289  Array<OneD, NekDouble> &outarray,
290  const NekDouble lambda,
291  CoeffState coeffstate = eLocal,
293  dirForcing = NullNekDouble1DArray);
294 
295 
296  /// Solve Advection Diffusion Reaction
297  inline void LinearAdvectionReactionSolve(
298  const Array<OneD, Array<OneD, NekDouble> > &velocity,
299  const Array<OneD, const NekDouble> &inarray,
300  Array<OneD, NekDouble> &outarray,
301  const NekDouble lambda,
302  CoeffState coeffstate = eLocal,
304  dirForcing = NullNekDouble1DArray);
305 
306  ///
308  const Array<OneD, const NekDouble> &inarray,
309  Array<OneD, NekDouble> &outarray);
310 
311 
312  /// This function elementally evaluates the backward transformation
313  /// of the global spectral/hp element expansion.
314  inline void BwdTrans_IterPerExp (
315  const Array<OneD, const NekDouble> &inarray,
316  Array<OneD,NekDouble> &outarray);
317 
318  ///
319  inline void BwdTrans (
320  const Array<OneD,
321  const NekDouble> &inarray,
322  Array<OneD,NekDouble> &outarray,
323  CoeffState coeffstate = eLocal);
324 
325  /// This function calculates the coordinates of all the elemental
326  /// quadrature points \f$\boldsymbol{x}_i\f$.
327  inline void GetCoords(
328  Array<OneD, NekDouble> &coord_0,
331 
332  // Homogeneous transforms
333  inline void HomogeneousFwdTrans(
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 HomogeneousBwdTrans(
341  const Array<OneD, const NekDouble> &inarray,
342  Array<OneD, NekDouble> &outarray,
343  CoeffState coeffstate = eLocal,
344  bool Shuff = true,
345  bool UnShuff = true);
346 
347  inline void DealiasedProd(
348  const Array<OneD, NekDouble> &inarray1,
349  const Array<OneD, NekDouble> &inarray2,
350  Array<OneD, NekDouble> &outarray,
351  CoeffState coeffstate = eLocal);
352 
353  inline void DealiasedDotProd(
354  const Array<OneD, Array<OneD, NekDouble> > &inarray1,
355  const Array<OneD, Array<OneD, NekDouble> > &inarray2,
356  Array<OneD, Array<OneD, NekDouble> > &outarray,
357  CoeffState coeffstate = eLocal);
358 
359  inline void GetBCValues(
360  Array<OneD, NekDouble> &BndVals,
361  const Array<OneD, NekDouble> &TotField,
362  int BndID);
363 
364  inline void NormVectorIProductWRTBase(
367  Array<OneD, NekDouble> &outarray,
368  int BndID);
369 
370  inline void NormVectorIProductWRTBase(
372  Array<OneD, NekDouble> &outarray);
373 
374  /// Apply geometry information to each expansion.
376 
377  /// Reset geometry information and reset matrices
379  {
380  v_Reset();
381  }
382 
383  void WriteTecplotHeader(std::ostream &outfile,
384  std::string var = "")
385  {
386  v_WriteTecplotHeader(outfile, var);
387  }
388 
390  std::ostream &outfile,
391  int expansion = -1)
392  {
393  v_WriteTecplotZone(outfile, expansion);
394  }
395 
396  void WriteTecplotField(std::ostream &outfile,
397  int expansion = -1)
398  {
399  v_WriteTecplotField(outfile, expansion);
400  }
401 
402  void WriteTecplotConnectivity(std::ostream &outfile,
403  int expansion = -1)
404  {
405  v_WriteTecplotConnectivity(outfile, expansion);
406  }
407 
408  MULTI_REGIONS_EXPORT void WriteVtkHeader(std::ostream &outfile);
409  MULTI_REGIONS_EXPORT void WriteVtkFooter(std::ostream &outfile);
410 
411  void WriteVtkPieceHeader(std::ostream &outfile, int expansion,
412  int istrip = 0)
413  {
414  v_WriteVtkPieceHeader(outfile, expansion, istrip);
415  }
416 
418  std::ostream &outfile,
419  int expansion);
420 
422  std::ostream &outfile,
423  int expansion,
424  std::string var = "v")
425  {
426  v_WriteVtkPieceData(outfile, expansion, var);
427  }
428 
429  /// This function returns the dimension of the coordinates of the
430  /// element \a eid.
431  // inline
432  MULTI_REGIONS_EXPORT int GetCoordim(int eid);
433 
434  /// Set the \a i th coefficiient in \a m_coeffs to value \a val
435  inline void SetCoeff(int i, NekDouble val);
436 
437  /// Set the \a i th coefficiient in #m_coeffs to value \a val
438  inline void SetCoeffs(int i, NekDouble val);
439 
440  /// Set the #m_coeffs array to inarray
441  inline void SetCoeffsArray(Array<OneD, NekDouble> &inarray);
442 
443  /// This function returns (a reference to) the array
444  /// \f$\boldsymbol{\hat{u}}_l\f$ (implemented as #m_coeffs)
445  /// containing all local expansion coefficients.
446  inline const Array<OneD, const NekDouble> &GetCoeffs() const;
447 
448  /// Impose Dirichlet Boundary Conditions onto Array
449  inline void ImposeDirichletConditions(
450  Array<OneD,NekDouble>& outarray);
451 
452 
453  /// Fill Bnd Condition expansion from the values stored in expansion
454  inline void FillBndCondFromField(void);
455 
456  /// Fill Bnd Condition expansion in nreg from the values stored in expansion
457  inline void FillBndCondFromField(const int nreg);
458 
459  /// Gathers the global coefficients \f$\boldsymbol{\hat{u}}_g\f$
460  /// from the local coefficients \f$\boldsymbol{\hat{u}}_l\f$.
461  // inline
462  MULTI_REGIONS_EXPORT inline void LocalToGlobal(bool useComm = true);
463 
465  const Array<OneD, const NekDouble> &inarray,
466  Array<OneD,NekDouble> &outarray,
467  bool useComm = true);
468 
469  /// Scatters from the global coefficients
470  /// \f$\boldsymbol{\hat{u}}_g\f$ to the local coefficients
471  /// \f$\boldsymbol{\hat{u}}_l\f$.
472  // inline
473  MULTI_REGIONS_EXPORT inline void GlobalToLocal(void);
474 
475  /**
476  * This operation is evaluated as:
477  * \f{tabbing}
478  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
479  * > > Do \= $i=$ $0,N_m^e-1$ \\
480  * > > > $\boldsymbol{\hat{u}}^{e}[i] = \mbox{sign}[e][i] \cdot
481  * \boldsymbol{\hat{u}}_g[\mbox{map}[e][i]]$ \\
482  * > > continue \\
483  * > continue
484  * \f}
485  * where \a map\f$[e][i]\f$ is the mapping array and \a
486  * sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
487  * correct modal connectivity between the different elements (both
488  * these arrays are contained in the data member #m_locToGloMap). This
489  * operation is equivalent to the scatter operation
490  * \f$\boldsymbol{\hat{u}}_l=\mathcal{A}\boldsymbol{\hat{u}}_g\f$,
491  * where \f$\mathcal{A}\f$ is the
492  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
493  *
494  * @param inarray An array of size \f$N_\mathrm{dof}\f$
495  * containing the global degrees of freedom
496  * \f$\boldsymbol{x}_g\f$.
497  * @param outarray The resulting local degrees of freedom
498  * \f$\boldsymbol{x}_l\f$ will be stored in this
499  * array of size \f$N_\mathrm{eof}\f$.
500  */
502  const Array<OneD, const NekDouble> &inarray,
503  Array<OneD,NekDouble> &outarray);
504 
505  /// Get the \a i th value (coefficient) of #m_coeffs
506  inline NekDouble GetCoeff(int i);
507 
508  /// Get the \a i th value (coefficient) of #m_coeffs
509  inline NekDouble GetCoeffs(int i);
510 
511  /// This function returns (a reference to) the array
512  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) containing the
513  /// function \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the
514  /// quadrature points.
515  // inline
517  &GetPhys() const;
518 
519  /// This function calculates the \f$L_\infty\f$ error of the global
520  /// spectral/hp element approximation.
522  const Array<OneD, const NekDouble> &inarray,
524 
525  /// This function calculates the \f$L_2\f$ error with
526  /// respect to soln of the global
527  /// spectral/hp element approximation.
529  const Array<OneD, const NekDouble> &inarray,
531  {
532  return v_L2(inarray, soln);
533  }
534 
535  /// Calculates the \f$H^1\f$ error of the global spectral/hp
536  /// element approximation.
538  const Array<OneD, const NekDouble> &inarray,
540 
542  {
543  return v_Integral(inarray);
544  }
545 
546  /// This function calculates the energy associated with
547  /// each one of the modesof a 3D homogeneous nD expansion
549  {
550  return v_HomogeneousEnergy();
551  }
552 
553  /// This function sets the Spectral Vanishing Viscosity
554  /// in homogeneous1D expansion.
556  {
558  }
559 
560  /// This function returns a vector containing the wave
561  /// numbers in z-direction associated
562  /// with the 3D homogenous expansion. Required if a
563  /// parellelisation is applied in the Fourier direction
565  {
566  return v_GetZIDs();
567  }
568 
569  /// This function returns the transposition class
570  /// associaed with the homogeneous expansion.
572  {
573  return v_GetTransposition();
574  }
575 
576  /// This function returns the Width of homogeneous direction
577  /// associaed with the homogeneous expansion.
579  {
580  return v_GetHomoLen();
581  }
582 
583  /// This function returns a vector containing the wave
584  /// numbers in y-direction associated
585  /// with the 3D homogenous expansion. Required if a
586  /// parellelisation is applied in the Fourier direction
588  {
589  return v_GetYIDs();
590  }
591 
592  /// This function interpolates the physical space points in
593  /// \a inarray to \a outarray using the same points defined in the
594  /// expansion but where the number of points are rescaled
595  /// by \a 1DScale
597  const NekDouble scale,
598  const Array<OneD, NekDouble> &inarray,
599  Array<OneD, NekDouble> &outarray)
600  {
601  v_PhysInterp1DScaled(scale, inarray,outarray);
602  }
603 
604  /// This function Galerkin projects the physical space points in
605  /// \a inarray to \a outarray where inarray is assumed to
606  /// be defined in the expansion but where the number of
607  /// points are rescaled by \a 1DScale
609  const NekDouble scale,
610  const Array<OneD, NekDouble> &inarray,
611  Array<OneD, NekDouble> &outarray)
612  {
613  v_PhysGalerkinProjection1DScaled(scale, inarray, outarray);
614  }
615 
616  /// This function returns the number of elements in the expansion.
617  inline int GetExpSize(void);
618 
619 
620  /// This function returns the number of elements in the
621  /// expansion which may be different for a homogeoenous extended
622  /// expansionp.
623  inline int GetNumElmts(void)
624  {
625  return v_GetNumElmts();
626  }
627 
628  /// This function returns the vector of elements in the expansion.
629  inline const boost::shared_ptr<LocalRegions::ExpansionVector>
630  GetExp() const;
631 
632  /// This function returns (a shared pointer to) the local elemental
633  /// expansion of the \f$n^{\mathrm{th}}\f$ element.
634  inline LocalRegions::ExpansionSharedPtr& GetExp(int n) const;
635 
636  /// This function returns (a shared pointer to) the local elemental
637  /// expansion containing the arbitrary point given by \a gloCoord.
639  const Array<OneD, const NekDouble> &gloCoord);
640 
641  /** This function returns the index of the local elemental
642  * expansion containing the arbitrary point given by \a gloCoord.
643  **/
645  const Array<OneD, const NekDouble> &gloCoord,
646  NekDouble tol = 0.0,
647  bool returnNearestElmt = false);
648 
649  /** This function returns the index and the Local
650  * Cartesian Coordinates \a locCoords of the local
651  * elemental expansion containing the arbitrary point
652  * given by \a gloCoords.
653  **/
655  const Array<OneD, const NekDouble> &gloCoords,
656  Array<OneD, NekDouble> &locCoords,
657  NekDouble tol = 0.0,
658  bool returnNearestElmt = false);
659 
660  /// Get the start offset position for a global list of #m_coeffs
661  /// correspoinding to element n.
662  inline int GetCoeff_Offset(int n) const;
663 
664  /// Get the start offset position for a global list of m_phys
665  /// correspoinding to element n.
666  inline int GetPhys_Offset(int n) const;
667 
668  /// This function returns (a reference to) the array
669  /// \f$\boldsymbol{\hat{u}}_l\f$ (implemented as #m_coeffs)
670  /// containing all local expansion coefficients.
672 
673  /// This function returns (a reference to) the array
674  /// \f$\boldsymbol{u}_l\f$ (implemented as #m_phys) containing the
675  /// function \f$u^{\delta}(\boldsymbol{x})\f$ evaluated at the
676  /// quadrature points.
678 
679  inline void PhysDeriv(
680  Direction edir,
681  const Array<OneD, const NekDouble> &inarray,
682  Array<OneD, NekDouble> &out_d);
683 
684  /// This function discretely evaluates the derivative of a function
685  /// \f$f(\boldsymbol{x})\f$ on the domain consisting of all
686  /// elements of the expansion.
687  inline void PhysDeriv(
688  const Array<OneD, const NekDouble> &inarray,
689  Array<OneD, NekDouble> &out_d0,
692 
693  inline void PhysDeriv(
694  const int dir,
695  const Array<OneD, const NekDouble> &inarray,
696  Array<OneD, NekDouble> &out_d);
697 
698  inline void CurlCurl(
701 
702  // functions associated with DisContField
705 
706  inline boost::shared_ptr<ExpList> &UpdateBndCondExpansion(int i);
707 
708  inline void Upwind(
709  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
710  const Array<OneD, const NekDouble> &Fwd,
711  const Array<OneD, const NekDouble> &Bwd,
713 
714  inline void Upwind(
716  const Array<OneD, const NekDouble> &Fwd,
717  const Array<OneD, const NekDouble> &Bwd,
719 
720  /**
721  * Return a reference to the trace space associated with this
722  * expansion list.
723  */
724  inline boost::shared_ptr<ExpList> &GetTrace();
725 
726  inline boost::shared_ptr<AssemblyMapDG> &GetTraceMap(void);
727 
728  inline const Array<OneD, const int> &GetTraceBndMap(void);
729 
730  inline void GetNormals(Array<OneD, Array<OneD, NekDouble> > &normals);
731 
732  inline void AddTraceIntegral(
735  Array<OneD, NekDouble> &outarray);
736 
737  inline void AddTraceIntegral(
739  Array<OneD, NekDouble> &outarray);
740 
741  inline void AddFwdBwdTraceIntegral(
742  const Array<OneD, const NekDouble> &Fwd,
743  const Array<OneD, const NekDouble> &Bwd,
744  Array<OneD, NekDouble> &outarray);
745 
746  inline void GetFwdBwdTracePhys(
748  Array<OneD,NekDouble> &Bwd);
749 
750  inline void GetFwdBwdTracePhys(
751  const Array<OneD,const NekDouble> &field,
753  Array<OneD,NekDouble> &Bwd);
754 
755  inline const std::vector<bool> &GetLeftAdjacentFaces(void) const;
756 
757  inline void ExtractTracePhys(Array<OneD,NekDouble> &outarray);
758 
759  inline void ExtractTracePhys(
760  const Array<OneD, const NekDouble> &inarray,
761  Array<OneD,NekDouble> &outarray);
762 
763  inline const Array<OneD, const SpatialDomains::
765 
766  inline Array<OneD, SpatialDomains::
768 
769  inline void EvaluateBoundaryConditions(
770  const NekDouble time = 0.0,
771  const std::string varName = "",
774 
775  // Routines for continous matrix solution
776  /// This function calculates the result of the multiplication of a
777  /// matrix of type specified by \a mkey with a vector given by \a
778  /// inarray.
779  inline void GeneralMatrixOp(
780  const GlobalMatrixKey &gkey,
781  const Array<OneD,const NekDouble> &inarray,
782  Array<OneD, NekDouble> &outarray,
783  CoeffState coeffstate = eLocal);
784 
786  const GlobalMatrixKey &gkey,
787  const Array<OneD,const NekDouble> &inarray,
788  Array<OneD, NekDouble> &outarray);
789 
790  inline void SetUpPhysNormals();
791 
792  inline void GetBoundaryToElmtMap(Array<OneD, int> &ElmtID,
793  Array<OneD,int> &EdgeID);
794 
795  inline void GetBndElmtExpansion(int i,
796  boost::shared_ptr<ExpList> &result,
797  const bool DeclareCoeffPhysArrays = true);
798 
799  inline void ExtractElmtToBndPhys(int i,
801  Array<OneD, NekDouble> &boundary);
802 
803  inline void ExtractPhysToBndElmt(int i,
804  const Array<OneD, const NekDouble> &phys,
805  Array<OneD, NekDouble> &bndElmt);
806 
807  inline void ExtractPhysToBnd(int i,
808  const Array<OneD, const NekDouble> &phys,
810 
811  inline void GetBoundaryNormals(int i,
812  Array<OneD, Array<OneD, NekDouble> > &normals);
813 
815  std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef,
816  int NumHomoDir = 0,
819  std::vector<NekDouble> &HomoLen =
821  bool homoStrips = false,
822  std::vector<unsigned int> &HomoSIDs =
824  std::vector<unsigned int> &HomoZIDs =
826  std::vector<unsigned int> &HomoYIDs =
828 
829 
831  {
832  return m_globalOptParam;
833  }
834 
835  std::map<int, RobinBCInfoSharedPtr> GetRobinBCInfo()
836  {
837  return v_GetRobinBCInfo();
838  }
839 
841  PeriodicMap &periodicVerts,
842  PeriodicMap &periodicEdges,
843  PeriodicMap &periodicFaces = NullPeriodicMap)
844  {
845  v_GetPeriodicEntities(periodicVerts, periodicEdges, periodicFaces);
846  }
847 
848  std::vector<LibUtilities::FieldDefinitionsSharedPtr>
850  {
851  return v_GetFieldDefinitions();
852  }
853 
854 
855  void GetFieldDefinitions(std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef)
856  {
857  v_GetFieldDefinitions(fielddef);
858  }
859 
860 
861 
862  /// Append the element data listed in elements
863  /// fielddef->m_ElementIDs onto fielddata
866  std::vector<NekDouble> &fielddata)
867  {
868  v_AppendFieldData(fielddef,fielddata);
869  }
870 
871 
872  /// Append the data in coeffs listed in elements
873  /// fielddef->m_ElementIDs onto fielddata
876  std::vector<NekDouble> &fielddata,
877  Array<OneD, NekDouble> &coeffs)
878  {
879  v_AppendFieldData(fielddef,fielddata,coeffs);
880  }
881 
882 
883  /** \brief Extract the data in fielddata into the coeffs
884  * using the basic ExpList Elemental expansions rather
885  * than planes in homogeneous case
886  */
889  std::vector<NekDouble> &fielddata,
890  std::string &field,
891  Array<OneD, NekDouble> &coeffs);
892 
893 
894  /** \brief Extract the data from fromField using
895  * fromExpList the coeffs using the basic ExpList
896  * Elemental expansions rather than planes in homogeneous
897  * case
898  */
900  const boost::shared_ptr<ExpList> &fromExpList,
901  const Array<OneD, const NekDouble> &fromCoeffs,
902  Array<OneD, NekDouble> &toCoeffs);
903 
904 
905  //Extract data in fielddata into the m_coeffs_list for the 3D stability analysis (base flow is 2D)
908  std::vector<NekDouble> &fielddata,
909  std::string &field,
910  Array<OneD, NekDouble> &coeffs);
911 
912 
913  /// Returns a shared pointer to the current object.
914  boost::shared_ptr<ExpList> GetSharedThisPtr()
915  {
916  return shared_from_this();
917  }
918 
919  /// Returns the session object
920  boost::shared_ptr<LibUtilities::SessionReader> GetSession() const
921  {
922  return m_session;
923  }
924 
925  /// Returns the comm object
926  boost::shared_ptr<LibUtilities::Comm> GetComm()
927  {
928  return m_comm;
929  }
930 
932  {
933  return m_graph;
934  }
935 
936  // Wrapper functions for Homogeneous Expansions
938  {
939  return v_GetHomogeneousBasis();
940  }
941 
942  boost::shared_ptr<ExpList> &GetPlane(int n)
943  {
944  return v_GetPlane(n);
945  }
946 
947  //expansion type
949 
953 
955 
956  protected:
957  /// Definition of the total number of degrees of freedom and
958  /// quadrature points and offsets to access data
959  void SetCoeffPhysOffsets();
960 
961  boost::shared_ptr<DNekMat> GenGlobalMatrixFull(
962  const GlobalLinSysKey &mkey,
963  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
964 
965  /// Communicator
967 
968  /// Session
970 
971  /// Mesh associated with this expansion list.
973 
974  /// The total number of local degrees of freedom. #m_ncoeffs
975  /// \f$=N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_l\f$
977 
978  /** The total number of quadrature points. #m_npoints
979  *\f$=Q_{\mathrm{tot}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_Q\f$
980  **/
982 
983  /**
984  * \brief Concatenation of all local expansion coefficients.
985  *
986  * The array of length #m_ncoeffs\f$=N_{\mathrm{eof}}\f$ which is
987  * the concatenation of the local expansion coefficients
988  * \f$\hat{u}_n^e\f$ over all \f$N_{\mathrm{el}}\f$ elements
989  * \f[\mathrm{\texttt{m\_coeffs}}=\boldsymbol{\hat{u}}_{l} =
990  * \underline{\boldsymbol{\hat{u}}}^e = \left [ \begin{array}{c}
991  * \boldsymbol{\hat{u}}^{1} \ \
992  * \boldsymbol{\hat{u}}^{2} \ \
993  * \vdots \ \
994  * \boldsymbol{\hat{u}}^{{{N_{\mathrm{el}}}}} \end{array} \right ],
995  * \quad
996  * \mathrm{where}\quad \boldsymbol{\hat{u}}^{e}[n]=\hat{u}_n^{e}\f]
997  */
999 
1000  /**
1001  * \brief The global expansion evaluated at the quadrature points
1002  *
1003  * The array of length #m_npoints\f$=Q_{\mathrm{tot}}\f$ containing
1004  * the evaluation of \f$u^{\delta}(\boldsymbol{x})\f$ at the
1005  * quadrature points \f$\boldsymbol{x}_i\f$.
1006  * \f[\mathrm{\texttt{m\_phys}}=\boldsymbol{u}_{l} =
1007  * \underline{\boldsymbol{u}}^e = \left [ \begin{array}{c}
1008  * \boldsymbol{u}^{1} \ \
1009  * \boldsymbol{u}^{2} \ \
1010  * \vdots \ \
1011  * \boldsymbol{u}^{{{N_{\mathrm{el}}}}} \end{array} \right ],\quad
1012  * \mathrm{where}\quad
1013  * \boldsymbol{u}^{e}[i]=u^{\delta}(\boldsymbol{x}_i)\f]
1014  */
1016 
1017  /**
1018  * \brief The state of the array #m_phys.
1019  *
1020  * Indicates whether the array #m_phys, created to contain the
1021  * evaluation of \f$u^{\delta}(\boldsymbol{x})\f$ at the quadrature
1022  * points \f$\boldsymbol{x}_i\f$, is filled with these values.
1023  */
1025 
1026  /**
1027  * \brief The list of local expansions.
1028  *
1029  * The (shared pointer to the) vector containing (shared pointers
1030  * to) all local expansions. The fact that the local expansions are
1031  * all stored as a (pointer to a) #StdExpansion, the abstract base
1032  * class for all local expansions, allows a general implementation
1033  * where most of the routines for the derived classes are defined
1034  * in the #ExpList base class.
1035  */
1036  boost::shared_ptr<LocalRegions::ExpansionVector> m_exp;
1037 
1039 
1040  /// Offset of elemental data into the array #m_coeffs
1041  std::vector<int> m_coll_coeff_offset;
1042 
1043  /// Offset of elemental data into the array #m_phys
1044  std::vector<int> m_coll_phys_offset;
1045 
1046  /// Offset of elemental data into the array #m_coeffs
1048 
1049  /// Offset of elemental data into the array #m_phys
1051 
1053 
1054  BlockMatrixMapShPtr m_blockMat;
1055 
1056  //@todo should this be in ExpList or ExpListHomogeneous1D.cpp
1057  // it's a bool which determine if the expansion is in the wave space (coefficient space)
1058  // or not
1060 
1061  /// Mapping from geometry ID of element to index inside #m_exp
1062  boost::unordered_map<int, int> m_elmtToExpId;
1063 
1064  /// This function assembles the block diagonal matrix of local
1065  /// matrices of the type \a mtype.
1067  const GlobalMatrixKey &gkey);
1068 
1070  const GlobalMatrixKey &gkey);
1071 
1072  void MultiplyByBlockMatrix(
1073  const GlobalMatrixKey &gkey,
1074  const Array<OneD,const NekDouble> &inarray,
1075  Array<OneD, NekDouble> &outarray);
1076 
1077  /// Generates a global matrix from the given key and map.
1078  boost::shared_ptr<GlobalMatrix> GenGlobalMatrix(
1079  const GlobalMatrixKey &mkey,
1080  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
1081 
1082 
1083  void GlobalEigenSystem(
1084  const boost::shared_ptr<DNekMat> &Gmat,
1085  Array<OneD, NekDouble> &EigValsReal,
1086  Array<OneD, NekDouble> &EigValsImag,
1087  Array<OneD, NekDouble> &EigVecs
1089 
1090 
1091  /// This operation constructs the global linear system of type \a
1092  /// mkey.
1093  boost::shared_ptr<GlobalLinSys> GenGlobalLinSys(
1094  const GlobalLinSysKey &mkey,
1095  const boost::shared_ptr<AssemblyMapCG> &locToGloMap);
1096 
1097  /// Generate a GlobalLinSys from information provided by the key
1098  /// "mkey" and the mapping provided in LocToGloBaseMap.
1099  boost::shared_ptr<GlobalLinSys> GenGlobalBndLinSys(
1100  const GlobalLinSysKey &mkey,
1101  const AssemblyMapSharedPtr &locToGloMap);
1102 
1104  {
1106  }
1107 
1108  // Virtual prototypes
1109 
1110  virtual int v_GetNumElmts(void)
1111  {
1112  return (*m_exp).size();
1113  }
1114 
1116  &v_GetBndCondExpansions(void);
1117 
1118  virtual boost::shared_ptr<ExpList> &v_UpdateBndCondExpansion(int i);
1119 
1120  virtual void v_Upwind(
1121  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
1122  const Array<OneD, const NekDouble> &Fwd,
1123  const Array<OneD, const NekDouble> &Bwd,
1125 
1126  virtual void v_Upwind(
1127  const Array<OneD, const NekDouble> &Vn,
1128  const Array<OneD, const NekDouble> &Fwd,
1129  const Array<OneD, const NekDouble> &Bwd,
1131 
1132  virtual boost::shared_ptr<ExpList> &v_GetTrace();
1133 
1134  virtual boost::shared_ptr<AssemblyMapDG> &v_GetTraceMap();
1135 
1136  virtual const Array<OneD, const int> &v_GetTraceBndMap();
1137 
1138  virtual void v_GetNormals(
1139  Array<OneD, Array<OneD, NekDouble> > &normals);
1140 
1141  virtual void v_AddTraceIntegral(
1142  const Array<OneD, const NekDouble> &Fx,
1143  const Array<OneD, const NekDouble> &Fy,
1144  Array<OneD, NekDouble> &outarray);
1145 
1146  virtual void v_AddTraceIntegral(
1147  const Array<OneD, const NekDouble> &Fn,
1148  Array<OneD, NekDouble> &outarray);
1149 
1150  virtual void v_AddFwdBwdTraceIntegral(
1151  const Array<OneD, const NekDouble> &Fwd,
1152  const Array<OneD, const NekDouble> &Bwd,
1153  Array<OneD, NekDouble> &outarray);
1154 
1155  virtual void v_GetFwdBwdTracePhys(
1156  Array<OneD,NekDouble> &Fwd,
1157  Array<OneD,NekDouble> &Bwd);
1158 
1159  virtual void v_GetFwdBwdTracePhys(
1160  const Array<OneD,const NekDouble> &field,
1161  Array<OneD,NekDouble> &Fwd,
1162  Array<OneD,NekDouble> &Bwd);
1163 
1164  virtual const std::vector<bool> &v_GetLeftAdjacentFaces(void) const;
1165 
1166  virtual void v_ExtractTracePhys(
1167  Array<OneD,NekDouble> &outarray);
1168 
1169  virtual void v_ExtractTracePhys(
1170  const Array<OneD, const NekDouble> &inarray,
1171  Array<OneD,NekDouble> &outarray);
1172 
1173  virtual void v_MultiplyByInvMassMatrix(
1174  const Array<OneD,const NekDouble> &inarray,
1175  Array<OneD, NekDouble> &outarray,
1176  CoeffState coeffstate);
1177 
1178  virtual void v_HelmSolve(
1179  const Array<OneD, const NekDouble> &inarray,
1180  Array<OneD, NekDouble> &outarray,
1181  const FlagList &flags,
1182  const StdRegions::ConstFactorMap &factors,
1183  const StdRegions::VarCoeffMap &varcoeff,
1184  const Array<OneD, const NekDouble> &dirForcing,
1185  const bool PhysSpaceForcing);
1186 
1188  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1189  const Array<OneD, const NekDouble> &inarray,
1190  Array<OneD, NekDouble> &outarray,
1191  const NekDouble lambda,
1192  CoeffState coeffstate = eLocal,
1194  dirForcing = NullNekDouble1DArray);
1195 
1196  virtual void v_LinearAdvectionReactionSolve(
1197  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1198  const Array<OneD, const NekDouble> &inarray,
1199  Array<OneD, NekDouble> &outarray,
1200  const NekDouble lambda,
1201  CoeffState coeffstate = eLocal,
1203  dirForcing = NullNekDouble1DArray);
1204 
1205  // wrapper functions about virtual functions
1206  virtual void v_ImposeDirichletConditions(Array<OneD,NekDouble>& outarray);
1207 
1208  virtual void v_FillBndCondFromField();
1209 
1210  virtual void v_FillBndCondFromField(const int nreg);
1211 
1212  virtual void v_Reset();
1213 
1214  virtual void v_LocalToGlobal(bool UseComm);
1215 
1216  virtual void v_LocalToGlobal(
1217  const Array<OneD, const NekDouble> &inarray,
1218  Array<OneD,NekDouble> &outarray,
1219  bool UseComm);
1220 
1221  virtual void v_GlobalToLocal(void);
1222 
1223  virtual void v_GlobalToLocal(
1224  const Array<OneD, const NekDouble> &inarray,
1225  Array<OneD,NekDouble> &outarray);
1226 
1227  virtual void v_BwdTrans(
1228  const Array<OneD,const NekDouble> &inarray,
1229  Array<OneD, NekDouble> &outarray,
1230  CoeffState coeffstate);
1231 
1232  virtual void v_BwdTrans_IterPerExp(
1233  const Array<OneD,const NekDouble> &inarray,
1234  Array<OneD,NekDouble> &outarray);
1235 
1236  virtual void v_FwdTrans(
1237  const Array<OneD,const NekDouble> &inarray,
1238  Array<OneD, NekDouble> &outarray,
1239  CoeffState coeffstate);
1240 
1241  virtual void v_FwdTrans_IterPerExp(
1242  const Array<OneD,const NekDouble> &inarray,
1243  Array<OneD,NekDouble> &outarray);
1244 
1245  virtual void v_SmoothField(Array<OneD,NekDouble> &field);
1246 
1247  virtual void v_IProductWRTBase(
1248  const Array<OneD, const NekDouble> &inarray,
1249  Array<OneD, NekDouble> &outarray,
1250  CoeffState coeffstate);
1251 
1252  virtual void v_IProductWRTBase_IterPerExp(
1253  const Array<OneD,const NekDouble> &inarray,
1254  Array<OneD, NekDouble> &outarray);
1255 
1256  virtual void v_GeneralMatrixOp(
1257  const GlobalMatrixKey &gkey,
1258  const Array<OneD,const NekDouble> &inarray,
1259  Array<OneD, NekDouble> &outarray,
1260  CoeffState coeffstate);
1261 
1262  virtual void v_GetCoords(
1263  Array<OneD, NekDouble> &coord_0,
1264  Array<OneD, NekDouble> &coord_1,
1266 
1267  virtual void v_PhysDeriv(
1268  const Array<OneD, const NekDouble> &inarray,
1269  Array<OneD, NekDouble> &out_d0,
1270  Array<OneD, NekDouble> &out_d1,
1271  Array<OneD, NekDouble> &out_d2);
1272 
1273  virtual void v_PhysDeriv(
1274  const int dir,
1275  const Array<OneD, const NekDouble> &inarray,
1276  Array<OneD, NekDouble> &out_d);
1277 
1278  virtual void v_PhysDeriv(
1279  Direction edir,
1280  const Array<OneD, const NekDouble> &inarray,
1281  Array<OneD, NekDouble> &out_d);
1282 
1283  virtual void v_CurlCurl(
1284  Array<OneD, Array<OneD, NekDouble> > &Vel,
1285  Array<OneD, Array<OneD, NekDouble> > &Q);
1286 
1287  virtual void v_HomogeneousFwdTrans(
1288  const Array<OneD, const NekDouble> &inarray,
1289  Array<OneD, NekDouble> &outarray,
1290  CoeffState coeffstate = eLocal,
1291  bool Shuff = true,
1292  bool UnShuff = true);
1293 
1294  virtual void v_HomogeneousBwdTrans(
1295  const Array<OneD, const NekDouble> &inarray,
1296  Array<OneD, NekDouble> &outarray,
1297  CoeffState coeffstate = eLocal,
1298  bool Shuff = true,
1299  bool UnShuff = true);
1300 
1301  virtual void v_DealiasedProd(
1302  const Array<OneD, NekDouble> &inarray1,
1303  const Array<OneD, NekDouble> &inarray2,
1304  Array<OneD, NekDouble> &outarray,
1305  CoeffState coeffstate = eLocal);
1306 
1307  virtual void v_DealiasedDotProd(
1308  const Array<OneD, Array<OneD, NekDouble> > &inarray1,
1309  const Array<OneD, Array<OneD, NekDouble> > &inarray2,
1310  Array<OneD, Array<OneD, NekDouble> > &outarray,
1311  CoeffState coeffstate = eLocal);
1312 
1313  virtual void v_GetBCValues(
1314  Array<OneD, NekDouble> &BndVals,
1315  const Array<OneD, NekDouble> &TotField,
1316  int BndID);
1317 
1318  virtual void v_NormVectorIProductWRTBase(
1321  Array<OneD, NekDouble> &outarray,
1322  int BndID);
1323 
1324  virtual void v_NormVectorIProductWRTBase(
1325  Array<OneD, Array<OneD, NekDouble> > &V,
1326  Array<OneD, NekDouble> &outarray);
1327 
1328  virtual void v_SetUpPhysNormals();
1329 
1330  virtual void v_GetBoundaryToElmtMap(Array<OneD, int> &ElmtID,
1331  Array<OneD,int> &EdgeID);
1332 
1333  virtual void v_GetBndElmtExpansion(int i,
1334  boost::shared_ptr<ExpList> &result,
1335  const bool DeclareCoeffPhysArrays);
1336 
1337  virtual void v_ExtractElmtToBndPhys(int i,
1338  Array<OneD, NekDouble> &elmt,
1339  Array<OneD, NekDouble> &boundary);
1340 
1341  virtual void v_ExtractPhysToBndElmt(int i,
1342  const Array<OneD, const NekDouble> &phys,
1343  Array<OneD, NekDouble> &bndElmt);
1344 
1345  virtual void v_ExtractPhysToBnd(int i,
1346  const Array<OneD, const NekDouble> &phys,
1347  Array<OneD, NekDouble> &bnd);
1348 
1349  virtual void v_GetBoundaryNormals(int i,
1350  Array<OneD, Array<OneD, NekDouble> > &normals);
1351 
1352  virtual void v_ReadGlobalOptimizationParameters();
1353 
1354  virtual std::vector<LibUtilities::FieldDefinitionsSharedPtr>
1355  v_GetFieldDefinitions(void);
1356 
1357  virtual void v_GetFieldDefinitions(
1358  std::vector<LibUtilities::FieldDefinitionsSharedPtr> &fielddef);
1359 
1360  virtual void v_AppendFieldData(
1362  std::vector<NekDouble> &fielddata);
1363 
1364  virtual void v_AppendFieldData(
1366  std::vector<NekDouble> &fielddata,
1367  Array<OneD, NekDouble> &coeffs);
1368 
1369  virtual void v_ExtractDataToCoeffs(
1371  std::vector<NekDouble> &fielddata, std::string &field,
1372  Array<OneD, NekDouble> &coeffs);
1373 
1374  virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr<ExpList> &fromExpList, const Array<OneD, const NekDouble> &fromCoeffs, Array<OneD, NekDouble> &toCoeffs);
1375 
1376  virtual void v_WriteTecplotHeader(std::ostream &outfile,
1377  std::string var = "");
1378  virtual void v_WriteTecplotZone(std::ostream &outfile,
1379  int expansion);
1380  virtual void v_WriteTecplotField(std::ostream &outfile,
1381  int expansion);
1382  virtual void v_WriteTecplotConnectivity(std::ostream &outfile,
1383  int expansion);
1384  virtual void v_WriteVtkPieceHeader(
1385  std::ostream &outfile,
1386  int expansion,
1387  int istrip);
1388 
1389  virtual void v_WriteVtkPieceData(
1390  std::ostream &outfile,
1391  int expansion,
1392  std::string var);
1393 
1394  virtual NekDouble v_L2(
1395  const Array<OneD, const NekDouble> &phys,
1397 
1398  virtual NekDouble v_Integral (
1399  const Array<OneD, const NekDouble> &inarray);
1400 
1403  virtual NekDouble v_GetHomoLen(void);
1406 
1407  // 1D Scaling and projection
1408  virtual void v_PhysInterp1DScaled(
1409  const NekDouble scale, const Array<OneD, NekDouble> &inarray,
1410  Array<OneD, NekDouble> &outarray);
1411 
1412  virtual void v_PhysGalerkinProjection1DScaled(
1413  const NekDouble scale,
1414  const Array<OneD, NekDouble> &inarray,
1415  Array<OneD, NekDouble> &outarray);
1416 
1417  virtual void v_ClearGlobalLinSysManager(void);
1418 
1419  void ExtractFileBCs(const std::string &fileName,
1421  const std::string &varName,
1422  const boost::shared_ptr<ExpList> locExpList);
1423 
1424  // Utility function for a common case of retrieving a
1425  // BoundaryCondition from a boundary condition collection.
1428  GetBoundaryCondition(const SpatialDomains::
1429  BoundaryConditionCollection& collection,
1430  unsigned int index, const std::string& variable);
1431 
1432 
1433  private:
1435 
1438 
1439  virtual void v_EvaluateBoundaryConditions(
1440  const NekDouble time = 0.0,
1441  const std::string varName = "",
1444 
1445  virtual std::map<int, RobinBCInfoSharedPtr> v_GetRobinBCInfo(void);
1446 
1447 
1448  virtual void v_GetPeriodicEntities(
1449  PeriodicMap &periodicVerts,
1450  PeriodicMap &periodicEdges,
1451  PeriodicMap &periodicFaces);
1452 
1453  // Homogeneous direction wrapper functions.
1455  {
1456  ASSERTL0(false,
1457  "This method is not defined or valid for this class type");
1459  }
1460 
1461  // wrapper function to set viscosity for Homo1D expansion
1463  {
1464  ASSERTL0(false,
1465  "This method is not defined or valid for this class type");
1466  }
1467 
1468 
1469  virtual boost::shared_ptr<ExpList> &v_GetPlane(int n);
1470  };
1471 
1472 
1473  /// Shared pointer to an ExpList object.
1474  typedef boost::shared_ptr<ExpList> ExpListSharedPtr;
1475  /// An empty ExpList object.
1477  static ExpListSharedPtr NullExpListSharedPtr;
1478 
1479  // Inline routines follow.
1480 
1481  /**
1482  * Returns the total number of local degrees of freedom
1483  * \f$N_{\mathrm{eof}}=\sum_{e=1}^{{N_{\mathrm{el}}}}N^{e}_m\f$.
1484  */
1485  inline int ExpList::GetNcoeffs() const
1486  {
1487  return m_ncoeffs;
1488  }
1489 
1490  inline int ExpList::GetNcoeffs(const int eid) const
1491  {
1492  return (*m_exp)[eid]->GetNcoeffs();
1493  }
1494 
1495  /**
1496  * Evaulates the maximum number of modes in the elemental basis
1497  * order over all elements
1498  */
1500  {
1501  unsigned int i;
1502  int returnval = 0;
1503 
1504  for(i= 0; i < (*m_exp).size(); ++i)
1505  {
1506  returnval = (std::max)(returnval,
1507  (*m_exp)[i]->EvalBasisNumModesMax());
1508  }
1509 
1510  return returnval;
1511  }
1512 
1513  /**
1514  *
1515  */
1517  const
1518  {
1519  unsigned int i;
1520  Array<OneD,int> returnval((*m_exp).size(),0);
1521 
1522  for(i= 0; i < (*m_exp).size(); ++i)
1523  {
1524  returnval[i]
1525  = (std::max)(returnval[i],(*m_exp)[i]->EvalBasisNumModesMax());
1526  }
1527 
1528  return returnval;
1529  }
1530 
1531 
1532  /**
1533  *
1534  */
1535  inline int ExpList::GetTotPoints() const
1536  {
1537  return m_npoints;
1538  }
1539 
1540  inline int ExpList::GetTotPoints(const int eid) const
1541  {
1542  return (*m_exp)[eid]->GetTotPoints();
1543  }
1544 
1545 
1546  inline int ExpList::Get1DScaledTotPoints(const NekDouble scale) const
1547  {
1548  int returnval = 0;
1549  int cnt;
1550  int nbase = (*m_exp)[0]->GetNumBases();
1551 
1552  for(int i = 0; i < (*m_exp).size(); ++i)
1553  {
1554  cnt = 1;
1555  for(int j = 0; j < nbase; ++j)
1556  {
1557  cnt *= (int)(scale*((*m_exp)[i]->GetNumPoints(j)));
1558  }
1559  returnval += cnt;
1560  }
1561  return returnval;
1562  }
1563 
1564  /**
1565  *
1566  */
1567  inline int ExpList::GetNpoints() const
1568  {
1569  return m_npoints;
1570  }
1571 
1572 
1573  /**
1574  *
1575  */
1576  inline void ExpList::SetWaveSpace(const bool wavespace)
1577  {
1578  m_WaveSpace = wavespace;
1579  }
1580 
1581  /**
1582  *
1583  */
1584  inline bool ExpList::GetWaveSpace() const
1585  {
1586  return m_WaveSpace;
1587  }
1588 
1589  /// Set the \a i th value of\a m_phys to value \a val
1590  inline void ExpList::SetPhys(int i, NekDouble val)
1591  {
1592  m_phys[i] = val;
1593  }
1594 
1595 
1596  /**
1597  * This function fills the array \f$\boldsymbol{u}_l\f$, the evaluation
1598  * of the expansion at the quadrature points (implemented as #m_phys),
1599  * with the values of the array \a inarray.
1600  *
1601  * @param inarray The array containing the values where
1602  * #m_phys should be filled with.
1603  */
1604  inline void ExpList::SetPhys(
1605  const Array<OneD, const NekDouble> &inarray)
1606  {
1607  ASSERTL0(inarray.num_elements() == m_npoints,
1608  "Input array does not have correct number of elements.");
1609 
1610  Vmath::Vcopy(m_npoints,&inarray[0],1,&m_phys[0],1);
1611  m_physState = true;
1612  }
1613 
1614 
1616  {
1617  m_phys = inarray;
1618  }
1619 
1620 
1621  /**
1622  * @param physState \a true (=filled) or \a false (=not filled).
1623  */
1624  inline void ExpList::SetPhysState(const bool physState)
1625  {
1626  m_physState = physState;
1627  }
1628 
1629 
1630  /**
1631  * @return physState \a true (=filled) or \a false (=not filled).
1632  */
1633  inline bool ExpList::GetPhysState() const
1634  {
1635  return m_physState;
1636  }
1637 
1638  /**
1639  *
1640  */
1642  const Array<OneD, const NekDouble> &inarray,
1643  Array<OneD, NekDouble> &outarray,
1644  CoeffState coeffstate)
1645  {
1646  v_IProductWRTBase(inarray,outarray, coeffstate);
1647  }
1648 
1649  /**
1650  *
1651  */
1653  const Array<OneD, const NekDouble> &inarray,
1654  Array<OneD, NekDouble> &outarray)
1655  {
1656  v_IProductWRTBase_IterPerExp(inarray,outarray);
1657  }
1658 
1659  /**
1660  *
1661  */
1662  inline void ExpList::FwdTrans(
1663  const Array<OneD, const NekDouble> &inarray,
1664  Array<OneD, NekDouble> &outarray,
1665  CoeffState coeffstate)
1666  {
1667  v_FwdTrans(inarray,outarray,coeffstate);
1668  }
1669 
1670  /**
1671  *
1672  */
1674  const Array<OneD, const NekDouble> &inarray,
1675  Array<OneD,NekDouble> &outarray)
1676  {
1677  v_FwdTrans_IterPerExp(inarray,outarray);
1678  }
1679 
1680  /**
1681  *
1682  */
1684  {
1685  v_SmoothField(field);
1686  }
1687 
1688  /**
1689  *
1690  */
1691  inline void ExpList::BwdTrans (
1692  const Array<OneD, const NekDouble> &inarray,
1693  Array<OneD, NekDouble> &outarray,
1694  CoeffState coeffstate)
1695  {
1696  v_BwdTrans(inarray,outarray,coeffstate);
1697  }
1698 
1699  /**
1700  *
1701  */
1703  const Array<OneD, const NekDouble> &inarray,
1704  Array<OneD, NekDouble> &outarray)
1705  {
1706  v_BwdTrans_IterPerExp(inarray,outarray);
1707  }
1708 
1709 
1710  /**
1711  *
1712  */
1714  const Array<OneD,const NekDouble> &inarray,
1715  Array<OneD, NekDouble> &outarray,
1716  CoeffState coeffstate)
1717  {
1718  v_MultiplyByInvMassMatrix(inarray,outarray,coeffstate);
1719  }
1720 
1721  /**
1722  *
1723  */
1724  inline void ExpList::HelmSolve(
1725  const Array<OneD, const NekDouble> &inarray,
1726  Array<OneD, NekDouble> &outarray,
1727  const FlagList &flags,
1728  const StdRegions::ConstFactorMap &factors,
1729  const StdRegions::VarCoeffMap &varcoeff,
1730  const Array<OneD, const NekDouble> &dirForcing,
1731  const bool PhysSpaceForcing)
1732 
1733  {
1734  v_HelmSolve(inarray, outarray, flags, factors, varcoeff,
1735  dirForcing, PhysSpaceForcing);
1736  }
1737 
1738 
1739  /**
1740  *
1741  */
1743  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1744  const Array<OneD, const NekDouble> &inarray,
1745  Array<OneD, NekDouble> &outarray,
1746  const NekDouble lambda,
1747  CoeffState coeffstate,
1748  const Array<OneD, const NekDouble>& dirForcing)
1749  {
1750  v_LinearAdvectionDiffusionReactionSolve(velocity,inarray, outarray,
1751  lambda, coeffstate,dirForcing);
1752  }
1753 
1755  const Array<OneD, Array<OneD, NekDouble> > &velocity,
1756  const Array<OneD, const NekDouble> &inarray,
1757  Array<OneD, NekDouble> &outarray,
1758  const NekDouble lambda,
1759  CoeffState coeffstate,
1760  const Array<OneD, const NekDouble>& dirForcing)
1761  {
1762  v_LinearAdvectionReactionSolve(velocity,inarray, outarray,
1763  lambda, coeffstate,dirForcing);
1764  }
1765 
1766  /**
1767  *
1768  */
1770  Array<OneD, NekDouble> &coord_1,
1771  Array<OneD, NekDouble> &coord_2)
1772 
1773  {
1774  v_GetCoords(coord_0,coord_1,coord_2);
1775  }
1776 
1777  /**
1778  *
1779  */
1781  Array<OneD, NekDouble> &out_d0,
1782  Array<OneD, NekDouble> &out_d1,
1783  Array<OneD, NekDouble> &out_d2)
1784  {
1785  v_PhysDeriv(inarray,out_d0,out_d1,out_d2);
1786  }
1787 
1788  /**
1789  *
1790  */
1791  inline void ExpList::PhysDeriv(
1792  const int dir,
1793  const Array<OneD, const NekDouble> &inarray,
1794  Array<OneD, NekDouble> &out_d)
1795  {
1796  v_PhysDeriv(dir,inarray,out_d);
1797  }
1798 
1799  inline void ExpList::PhysDeriv(
1800  Direction edir,
1801  const Array<OneD, const NekDouble> &inarray,
1802  Array<OneD, NekDouble> &out_d)
1803  {
1804  v_PhysDeriv(edir, inarray,out_d);
1805  }
1806 
1807  inline void ExpList::CurlCurl(
1810  {
1811  v_CurlCurl(Vel, Q);
1812  }
1813 
1814  /**
1815  *
1816  */
1818  const Array<OneD, const NekDouble> &inarray,
1819  Array<OneD, NekDouble> &outarray,
1820  CoeffState coeffstate,
1821  bool Shuff,
1822  bool UnShuff)
1823  {
1824  v_HomogeneousFwdTrans(inarray,outarray,coeffstate,Shuff,UnShuff);
1825  }
1826 
1827  /**
1828  *
1829  */
1831  const Array<OneD, const NekDouble> &inarray,
1832  Array<OneD, NekDouble> &outarray,
1833  CoeffState coeffstate,
1834  bool Shuff,
1835  bool UnShuff)
1836  {
1837  v_HomogeneousBwdTrans(inarray,outarray,coeffstate,Shuff,UnShuff);
1838  }
1839 
1840  /**
1841  *
1842  */
1844  const Array<OneD, NekDouble> &inarray1,
1845  const Array<OneD, NekDouble> &inarray2,
1846  Array<OneD, NekDouble> &outarray,
1847  CoeffState coeffstate)
1848  {
1849  v_DealiasedProd(inarray1,inarray2,outarray,coeffstate);
1850  }
1851 
1852  /**
1853  *
1854  */
1856  const Array<OneD, Array<OneD, NekDouble> > &inarray1,
1857  const Array<OneD, Array<OneD, NekDouble> > &inarray2,
1858  Array<OneD, Array<OneD, NekDouble> > &outarray,
1859  CoeffState coeffstate)
1860  {
1861  v_DealiasedDotProd(inarray1,inarray2,outarray,coeffstate);
1862  }
1863 
1864  /**
1865  *
1866  */
1868  Array<OneD, NekDouble> &BndVals,
1869  const Array<OneD, NekDouble> &TotField,
1870  int BndID)
1871  {
1872  v_GetBCValues(BndVals,TotField,BndID);
1873  }
1874 
1875  /**
1876  *
1877  */
1881  Array<OneD, NekDouble> &outarray,
1882  int BndID)
1883  {
1884  v_NormVectorIProductWRTBase(V1,V2,outarray,BndID);
1885  }
1886 
1889  Array<OneD, NekDouble> &outarray)
1890  {
1891  v_NormVectorIProductWRTBase(V, outarray);
1892  }
1893 
1894  /**
1895  * @param eid The index of the element to be checked.
1896  * @return The dimension of the coordinates of the specific element.
1897  */
1898  inline int ExpList::GetCoordim(int eid)
1899  {
1900  ASSERTL2(eid <= (*m_exp).size(),
1901  "eid is larger than number of elements");
1902  return (*m_exp)[eid]->GetCoordim();
1903  }
1904 
1905  /**
1906  * @param i The index of m_coeffs to be set
1907  * @param val The value which m_coeffs[i] is to be set to.
1908  */
1909  inline void ExpList::SetCoeff(int i, NekDouble val)
1910  {
1911  m_coeffs[i] = val;
1912  }
1913 
1914 
1915  /**
1916  * @param i The index of #m_coeffs to be set.
1917  * @param val The value which #m_coeffs[i] is to be set to.
1918  */
1919  inline void ExpList::SetCoeffs(int i, NekDouble val)
1920  {
1921  m_coeffs[i] = val;
1922  }
1923 
1924 
1926  {
1927  m_coeffs = inarray;
1928  }
1929 
1930  /**
1931  * As the function returns a constant reference to a
1932  * <em>const Array</em>, it is not possible to modify the
1933  * underlying data of the array #m_coeffs. In order to do
1934  * so, use the function #UpdateCoeffs instead.
1935  *
1936  * @return (A constant reference to) the array #m_coeffs.
1937  */
1939  {
1940  return m_coeffs;
1941  }
1942 
1944  Array<OneD,NekDouble>& outarray)
1945  {
1946  v_ImposeDirichletConditions(outarray);
1947  }
1948 
1950  {
1952  }
1953 
1954  inline void ExpList::FillBndCondFromField(const int nreg)
1955  {
1956  v_FillBndCondFromField(nreg);
1957  }
1958 
1959  inline void ExpList::LocalToGlobal(bool useComm)
1960  {
1961  v_LocalToGlobal(useComm);
1962  }
1963 
1965  const Array<OneD, const NekDouble> &inarray,
1966  Array<OneD,NekDouble> &outarray,
1967  bool useComm)
1968  {
1969  v_LocalToGlobal(inarray, outarray,useComm);
1970  }
1971 
1972  inline void ExpList::GlobalToLocal(void)
1973  {
1974  v_GlobalToLocal();
1975  }
1976 
1977  /**
1978  * This operation is evaluated as:
1979  * \f{tabbing}
1980  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
1981  * > > Do \= $i=$ $0,N_m^e-1$ \\
1982  * > > > $\boldsymbol{\hat{u}}^{e}[i] = \mbox{sign}[e][i] \cdot
1983  * \boldsymbol{\hat{u}}_g[\mbox{map}[e][i]]$ \\
1984  * > > continue \\
1985  * > continue
1986  * \f}
1987  * where \a map\f$[e][i]\f$ is the mapping array and \a
1988  * sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
1989  * correct modal connectivity between the different elements (both
1990  * these arrays are contained in the data member #m_locToGloMap). This
1991  * operation is equivalent to the scatter operation
1992  * \f$\boldsymbol{\hat{u}}_l=\mathcal{A}\boldsymbol{\hat{u}}_g\f$, where
1993  * \f$\mathcal{A}\f$ is the
1994  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
1995  *
1996  * @param inarray An array of size \f$N_\mathrm{dof}\f$
1997  * containing the global degrees of freedom
1998  * \f$\boldsymbol{x}_g\f$.
1999  * @param outarray The resulting local degrees of freedom
2000  * \f$\boldsymbol{x}_l\f$ will be stored in this
2001  * array of size \f$N_\mathrm{eof}\f$.
2002  */
2004  const Array<OneD, const NekDouble> &inarray,
2005  Array<OneD,NekDouble> &outarray)
2006  {
2007  v_GlobalToLocal(inarray, outarray);
2008  }
2009 
2010 
2011  /**
2012  * @param i The index of #m_coeffs to be returned
2013  * @return The NekDouble held in #m_coeffs[i].
2014  */
2016  {
2017  return m_coeffs[i];
2018  }
2019 
2020  /**
2021  * @param i The index of #m_coeffs to be returned
2022  * @return The NekDouble held in #m_coeffs[i].
2023  */
2025  {
2026  return m_coeffs[i];
2027  }
2028 
2029  /**
2030  * As the function returns a constant reference to a
2031  * <em>const Array</em> it is not possible to modify the
2032  * underlying data of the array #m_phys. In order to do
2033  * so, use the function #UpdatePhys instead.
2034  *
2035  * @return (A constant reference to) the array #m_phys.
2036  */
2038  {
2039  return m_phys;
2040  }
2041 
2042  /**
2043  * @return \f$N_{\mathrm{el}}\f$, the number of elements in the
2044  * expansion.
2045  */
2046  inline int ExpList::GetExpSize(void)
2047  {
2048  return (*m_exp).size();
2049  }
2050 
2051 
2052  /**
2053  * @param n The index of the element concerned.
2054  *
2055  * @return (A shared pointer to) the local expansion of the
2056  * \f$n^{\mathrm{th}}\f$ element.
2057  */
2059  {
2060  return (*m_exp)[n];
2061  }
2062 
2063  /**
2064  * @return (A const shared pointer to) the local expansion vector #m_exp
2065  */
2066  inline const boost::shared_ptr<LocalRegions::ExpansionVector>
2067  ExpList::GetExp(void) const
2068  {
2069  return m_exp;
2070  }
2071 
2072 
2073  /**
2074  *
2075  */
2076  inline int ExpList::GetCoeff_Offset(int n) const
2077  {
2078  return m_coeff_offset[n];
2079  }
2080 
2081  /**
2082  *
2083  */
2084  inline int ExpList::GetPhys_Offset(int n) const
2085  {
2086  return m_phys_offset[n];
2087  }
2088 
2089  /**
2090  * If one wants to get hold of the underlying data without modifying
2091  * them, rather use the function #GetCoeffs instead.
2092  *
2093  * @return (A reference to) the array #m_coeffs.
2094  */
2096  {
2097  return m_coeffs;
2098  }
2099 
2100  /**
2101  * If one wants to get hold of the underlying data without modifying
2102  * them, rather use the function #GetPhys instead.
2103  *
2104  * @return (A reference to) the array #m_phys.
2105  */
2107  {
2108  m_physState = true;
2109  return m_phys;
2110  }
2111 
2112 
2113  // functions associated with DisContField
2116  {
2117  return v_GetBndCondExpansions();
2118  }
2119 
2120  inline boost::shared_ptr<ExpList> &ExpList::UpdateBndCondExpansion(int i)
2121  {
2122  return v_UpdateBndCondExpansion(i);
2123  }
2124 
2125  inline void ExpList::Upwind(
2126  const Array<OneD, const Array<OneD, NekDouble> > &Vec,
2127  const Array<OneD, const NekDouble> &Fwd,
2128  const Array<OneD, const NekDouble> &Bwd,
2129  Array<OneD, NekDouble> &Upwind)
2130  {
2131  v_Upwind(Vec, Fwd, Bwd, Upwind);
2132  }
2133 
2134  inline void ExpList::Upwind(
2135  const Array<OneD, const NekDouble> &Vn,
2136  const Array<OneD, const NekDouble> &Fwd,
2137  const Array<OneD, const NekDouble> &Bwd,
2138  Array<OneD, NekDouble> &Upwind)
2139  {
2140  v_Upwind(Vn, Fwd, Bwd, Upwind);
2141  }
2142 
2143  inline boost::shared_ptr<ExpList> &ExpList::GetTrace()
2144  {
2145  return v_GetTrace();
2146  }
2147 
2148  inline boost::shared_ptr<AssemblyMapDG> &ExpList::GetTraceMap()
2149  {
2150  return v_GetTraceMap();
2151  }
2152 
2154  {
2155  return v_GetTraceBndMap();
2156  }
2157 
2158  inline void ExpList::GetNormals(
2159  Array<OneD, Array<OneD, NekDouble> > &normals)
2160  {
2161  v_GetNormals(normals);
2162  }
2163 
2165  const Array<OneD, const NekDouble> &Fx,
2166  const Array<OneD, const NekDouble> &Fy,
2167  Array<OneD, NekDouble> &outarray)
2168  {
2169  v_AddTraceIntegral(Fx,Fy,outarray);
2170  }
2171 
2173  const Array<OneD, const NekDouble> &Fn,
2174  Array<OneD, NekDouble> &outarray)
2175  {
2176  v_AddTraceIntegral(Fn,outarray);
2177  }
2178 
2180  const Array<OneD, const NekDouble> &Fwd,
2181  const Array<OneD, const NekDouble> &Bwd,
2182  Array<OneD, NekDouble> &outarray)
2183  {
2184  v_AddFwdBwdTraceIntegral(Fwd,Bwd,outarray);
2185  }
2186 
2188  Array<OneD,NekDouble> &Fwd,
2189  Array<OneD,NekDouble> &Bwd)
2190  {
2191  v_GetFwdBwdTracePhys(Fwd,Bwd);
2192  }
2193 
2195  const Array<OneD,const NekDouble> &field,
2196  Array<OneD,NekDouble> &Fwd,
2197  Array<OneD,NekDouble> &Bwd)
2198  {
2199  v_GetFwdBwdTracePhys(field,Fwd,Bwd);
2200  }
2201 
2202  inline const std::vector<bool> &ExpList::GetLeftAdjacentFaces(void) const
2203  {
2204  return v_GetLeftAdjacentFaces();
2205  }
2206 
2208  {
2209  v_ExtractTracePhys(outarray);
2210  }
2211 
2212 
2214  const Array<OneD, const NekDouble> &inarray,
2215  Array<OneD,NekDouble> &outarray)
2216  {
2217  v_ExtractTracePhys(inarray,outarray);
2218  }
2219 
2222  {
2223  return v_GetBndConditions();
2224  }
2225 
2226 
2229  {
2230  return v_UpdateBndConditions();
2231  }
2232 
2234  const NekDouble time,
2235  const std::string varName,
2236  const NekDouble x2_in,
2237  const NekDouble x3_in)
2238  {
2239  v_EvaluateBoundaryConditions(time, varName, x2_in, x3_in);
2240  }
2241 
2242  // Routines for continous matrix solution
2243  /**
2244  * This operation is equivalent to the evaluation of
2245  * \f$\underline{\boldsymbol{M}}^e\boldsymbol{\hat{u}}_l\f$, that is,
2246  * \f[ \left[
2247  * \begin{array}{cccc}
2248  * \boldsymbol{M}^1 & 0 & \hspace{3mm}0 \hspace{3mm}& 0 \\
2249  * 0 & \boldsymbol{M}^2 & 0 & 0 \\
2250  * 0 & 0 & \ddots & 0 \\
2251  * 0 & 0 & 0 & \boldsymbol{M}^{N_{\mathrm{el}}} \end{array} \right]
2252  *\left [ \begin{array}{c}
2253  * \boldsymbol{\hat{u}}^{1} \\
2254  * \boldsymbol{\hat{u}}^{2} \\
2255  * \vdots \\
2256  * \boldsymbol{\hat{u}}^{{{N_{\mathrm{el}}}}} \end{array} \right ]\f]
2257  * where \f$\boldsymbol{M}^e\f$ are the local matrices of type
2258  * specified by the key \a mkey. The decoupling of the local matrices
2259  * allows for a local evaluation of the operation. However, rather than
2260  * a local matrix-vector multiplication, the local operations are
2261  * evaluated as implemented in the function
2262  * StdRegions#StdExpansion#GeneralMatrixOp.
2263  *
2264  * @param mkey This key uniquely defines the type matrix
2265  * required for the operation.
2266  * @param inarray The vector \f$\boldsymbol{\hat{u}}_l\f$ of
2267  * size \f$N_{\mathrm{eof}}\f$.
2268  * @param outarray The resulting vector of size
2269  * \f$N_{\mathrm{eof}}\f$.
2270  */
2272  const GlobalMatrixKey &gkey,
2273  const Array<OneD,const NekDouble> &inarray,
2274  Array<OneD, NekDouble> &outarray,
2275  CoeffState coeffstate)
2276  {
2277  v_GeneralMatrixOp(gkey,inarray,outarray,coeffstate);
2278  }
2279 
2280 
2282  {
2284  }
2285 
2287  Array<OneD,int> &EdgeID)
2288  {
2289  v_GetBoundaryToElmtMap(ElmtID,EdgeID);
2290  }
2291 
2292  inline void ExpList::GetBndElmtExpansion(int i,
2293  boost::shared_ptr<ExpList> &result,
2294  const bool DeclareCoeffPhysArrays)
2295  {
2296  v_GetBndElmtExpansion(i, result, DeclareCoeffPhysArrays);
2297  }
2298 
2300  Array<OneD, NekDouble> &elmt,
2301  Array<OneD, NekDouble> &boundary)
2302  {
2303  v_ExtractElmtToBndPhys(i, elmt, boundary);
2304  }
2305 
2307  const Array<OneD, const NekDouble> &phys,
2308  Array<OneD, NekDouble> &bndElmt)
2309  {
2310  v_ExtractPhysToBndElmt(i, phys, bndElmt);
2311  }
2312 
2313  inline void ExpList::ExtractPhysToBnd(int i,
2314  const Array<OneD, const NekDouble> &phys,
2316  {
2317  v_ExtractPhysToBnd(i, phys, bnd);
2318  }
2319 
2320  inline void ExpList::GetBoundaryNormals(int i,
2321  Array<OneD, Array<OneD, NekDouble> > &normals)
2322  {
2323  v_GetBoundaryNormals(i, normals);
2324  }
2325 
2327 
2328  } //end of namespace
2329 } //end of namespace
2330 
2331 #endif // EXPLIST_H
2332 
virtual NekDouble v_L2(const Array< OneD, const NekDouble > &phys, const Array< OneD, const NekDouble > &soln=NullNekDouble1DArray)
Definition: ExpList.cpp:1909
virtual void v_GlobalToLocal(void)
Definition: ExpList.cpp:2692
virtual void v_AddFwdBwdTraceIntegral(const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2459
virtual void v_SmoothField(Array< OneD, NekDouble > &field)
Definition: ExpList.cpp:736
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:514
void AddTraceIntegral(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:2164
virtual void v_GeneralMatrixOp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2736
const DNekScalBlkMatSharedPtr & GetBlockMatrix(const GlobalMatrixKey &gkey)
Definition: ExpList.cpp:914
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:2262
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
void SmoothField(Array< OneD, NekDouble > &field)
Smooth a field across elements.
Definition: ExpList.h:1683
void GetBndElmtExpansion(int i, boost::shared_ptr< ExpList > &result, const bool DeclareCoeffPhysArrays=true)
Definition: ExpList.h:2292
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
int GetExpIndex(const Array< OneD, const NekDouble > &gloCoord, NekDouble tol=0.0, bool returnNearestElmt=false)
Definition: ExpList.cpp:1364
virtual void v_Reset()
Reset geometry information, metrics, matrix managers and geometry information.
Definition: ExpList.cpp:1537
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:997
virtual void v_SetHomo1DSpecVanVisc(Array< OneD, NekDouble > visc)
Definition: ExpList.h:1462
virtual void v_GetBoundaryNormals(int i, Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.cpp:2944
void LocalToGlobal(bool useComm=true)
Gathers the global coefficients from the local coefficients .
Definition: ExpList.h:1959
void GeneralMatrixOp_IterPerExp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:929
virtual const Array< OneD, const boost::shared_ptr< ExpList > > & v_GetBndCondExpansions(void)
Definition: ExpList.cpp:2378
void ExtractFileBCs(const std::string &fileName, LibUtilities::CommSharedPtr comm, const std::string &varName, const boost::shared_ptr< ExpList > locExpList)
Definition: ExpList.cpp:2011
boost::shared_ptr< Transposition > TranspositionSharedPtr
void ExtractPhysToBnd(int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bnd)
Definition: ExpList.h:2313
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
ExpList()
The default constructor.
Definition: ExpList.cpp:95
Local coefficients.
static ExpListSharedPtr NullExpListSharedPtr
Definition: ExpList.h:1477
static Array< OneD, NekDouble > NullNekDouble1DArray
boost::unordered_map< int, int > m_elmtToExpId
Mapping from geometry ID of element to index inside m_exp.
Definition: ExpList.h:1062
void EvaluateBoundaryConditions(const NekDouble time=0.0, const std::string varName="", const NekDouble=NekConstants::kNekUnsetDouble, const NekDouble=NekConstants::kNekUnsetDouble)
Definition: ExpList.h:2233
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:1546
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1052
boost::shared_ptr< AssemblyMap > AssemblyMapSharedPtr
Definition: AssemblyMap.h:53
boost::shared_ptr< ExpList > & GetPlane(int n)
Definition: ExpList.h:942
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:571
void MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1713
virtual boost::shared_ptr< AssemblyMapDG > & v_GetTraceMap()
Definition: ExpList.cpp:2422
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:2084
virtual void v_ExtractPhysToBnd(int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bnd)
Definition: ExpList.cpp:2906
Array< OneD, NekDouble > & UpdateCoeffs()
This function returns (a reference to) the array (implemented as m_coeffs) containing all local expa...
Definition: ExpList.h:2095
virtual void v_FillBndCondFromField()
Definition: ExpList.cpp:2662
std::vector< Collection > CollectionVector
Definition: Collection.h:98
boost::shared_ptr< DNekMat > GenGlobalMatrixFull(const GlobalLinSysKey &mkey, const boost::shared_ptr< AssemblyMapCG > &locToGloMap)
Definition: ExpList.cpp:1134
NekDouble PhysIntegral(void)
This function integrates a function over the domain consisting of all the elements of the expansion...
Definition: ExpList.cpp:302
virtual void v_WriteTecplotZone(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1603
ExpansionType GetExpType(void)
Returns the type of the expansion.
Definition: ExpList.cpp:269
const Array< OneD, const int > & GetTraceBndMap(void)
Definition: ExpList.h:2153
boost::shared_ptr< FieldDefinitions > FieldDefinitionsSharedPtr
Definition: FieldIO.h:181
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2721
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:1516
const boost::shared_ptr< LocalRegions::ExpansionVector > GetExp() const
This function returns the vector of elements in the expansion.
Definition: ExpList.h:2067
std::vector< LibUtilities::FieldDefinitionsSharedPtr > GetFieldDefinitions()
Definition: ExpList.h:849
virtual void v_BwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1331
void GeneralGetFieldDefinitions(std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef, int NumHomoDir=0, Array< OneD, LibUtilities::BasisSharedPtr > &HomoBasis=LibUtilities::NullBasisSharedPtr1DArray, std::vector< NekDouble > &HomoLen=LibUtilities::NullNekDoubleVector, bool homoStrips=false, std::vector< unsigned int > &HomoSIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoZIDs=LibUtilities::NullUnsignedIntVector, std::vector< unsigned int > &HomoYIDs=LibUtilities::NullUnsignedIntVector)
Definition: ExpList.cpp:2090
static BasisSharedPtr NullBasisSharedPtr
Definition: Basis.h:358
#define MULTI_REGIONS_EXPORT
boost::shared_ptr< ExpList > & UpdateBndCondExpansion(int i)
Definition: ExpList.h:2120
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1691
virtual void v_GetNormals(Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.cpp:2435
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:564
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:252
void ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Impose Dirichlet Boundary Conditions onto Array.
Definition: ExpList.h:1943
void GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Definition: ExpList.h:2286
BlockMatrixMapShPtr m_blockMat
Definition: ExpList.h:1054
Array< OneD, NekDouble > m_phys
The global expansion evaluated at the quadrature points.
Definition: ExpList.h:1015
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:1742
virtual const std::vector< bool > & v_GetLeftAdjacentFaces(void) const
Definition: ExpList.cpp:2484
virtual void v_GetBndElmtExpansion(int i, boost::shared_ptr< ExpList > &result, const bool DeclareCoeffPhysArrays)
Definition: ExpList.cpp:2817
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:1754
void ApplyGeomInfo()
Apply geometry information to each expansion.
Definition: ExpList.cpp:1524
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:998
virtual void v_WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip)
Definition: ExpList.cpp:1822
NekDouble Integral(const Array< OneD, const NekDouble > &inarray)
Definition: ExpList.h:541
virtual void v_NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
Definition: ExpList.cpp:2599
virtual LibUtilities::TranspositionSharedPtr v_GetTransposition(void)
Definition: ExpList.cpp:1961
int GetExpSize(void)
This function returns the number of elements in the expansion.
Definition: ExpList.h:2046
void PhysDeriv(Direction edir, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &out_d)
Definition: ExpList.h:1799
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual boost::shared_ptr< ExpList > & v_UpdateBndCondExpansion(int i)
Definition: ExpList.cpp:2386
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:1294
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...
virtual void v_ExtractElmtToBndPhys(int i, Array< OneD, NekDouble > &elmt, Array< OneD, NekDouble > &boundary)
Definition: ExpList.cpp:2827
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:608
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:665
virtual void v_WriteTecplotHeader(std::ostream &outfile, std::string var="")
Definition: ExpList.cpp:1564
virtual NekDouble v_GetHomoLen(void)
Definition: ExpList.cpp:1969
void WriteVtkFooter(std::ostream &outfile)
Definition: ExpList.cpp:1816
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:1919
virtual void v_GetCoords(Array< OneD, NekDouble > &coord_0, Array< OneD, NekDouble > &coord_1, Array< OneD, NekDouble > &coord_2=NullNekDouble1DArray)
Definition: ExpList.cpp:2756
const Array< OneD, const boost::shared_ptr< ExpList > > & GetBndCondExpansions()
Definition: ExpList.h:2115
void WriteTecplotField(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:396
boost::shared_ptr< ExpList > GetSharedThisPtr()
Returns a shared pointer to the current object.
Definition: ExpList.h:914
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:2125
virtual Array< OneD, const NekDouble > v_HomogeneousEnergy(void)
Definition: ExpList.cpp:1953
virtual const Array< OneD, const int > & v_GetTraceBndMap()
Definition: ExpList.cpp:2430
virtual void v_WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var)
Definition: ExpList.cpp:1833
static std::vector< NekDouble > NullNekDoubleVector
void ExtractPhysToBndElmt(int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bndElmt)
Definition: ExpList.h:2306
virtual void v_ExtractPhysToBndElmt(int i, const Array< OneD, const NekDouble > &phys, Array< OneD, NekDouble > &bndElmt)
Definition: ExpList.cpp:2866
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:2271
void NormVectorIProductWRTBase(Array< OneD, const NekDouble > &V1, Array< OneD, const NekDouble > &V2, Array< OneD, NekDouble > &outarray, int BndID)
Definition: ExpList.h:1878
virtual void v_GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
Definition: ExpList.cpp:2591
void GetNormals(Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.h:2158
virtual void v_LocalToGlobal(bool UseComm)
Definition: ExpList.cpp:2676
void Reset()
Reset geometry information and reset matrices.
Definition: ExpList.h:378
boost::shared_ptr< LibUtilities::Comm > GetComm()
Returns the comm object.
Definition: ExpList.h:926
Array< OneD, int > m_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:1047
LibUtilities::BasisSharedPtr GetHomogeneousBasis(void)
Definition: ExpList.h:937
const NekOptimize::GlobalOptParamSharedPtr & GetGlobalOptParam(void)
Definition: ExpList.h:830
int GetNumElmts(void)
This function returns the number of elements in the expansion which may be different for a homogeoeno...
Definition: ExpList.h:623
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:55
boost::shared_ptr< ExpList > & GetTrace()
Definition: ExpList.h:2143
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:1769
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2654
virtual LibUtilities::BasisSharedPtr v_GetHomogeneousBasis(void)
Definition: ExpList.h:1454
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:2394
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:1870
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
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:1702
virtual void v_FwdTrans_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:704
virtual ~ExpList()
The default destructor.
Definition: ExpList.cpp:282
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:407
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:227
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:1275
void CurlCurl(Array< OneD, Array< OneD, NekDouble > > &Vel, Array< OneD, Array< OneD, NekDouble > > &Q)
Definition: ExpList.h:1807
void WriteVtkPieceData(std::ostream &outfile, int expansion, std::string var="v")
Definition: ExpList.h:421
void AddFwdBwdTraceIntegral(const Array< OneD, const NekDouble > &Fwd, const Array< OneD, const NekDouble > &Bwd, Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:2179
Array< OneD, int > m_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:1050
void FillBndCondFromField(void)
Fill Bnd Condition expansion from the values stored in expansion.
Definition: ExpList.h:1949
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:1652
NekDouble GetHomoLen(void)
This function returns the Width of homogeneous direction associaed with the homogeneous expansion...
Definition: ExpList.h:578
virtual NekDouble v_Integral(const Array< OneD, const NekDouble > &inarray)
Definition: ExpList.cpp:1939
virtual std::vector< LibUtilities::FieldDefinitionsSharedPtr > v_GetFieldDefinitions(void)
Definition: ExpList.cpp:2198
bool m_physState
The state of the array m_phys.
Definition: ExpList.h:1024
virtual void v_WriteTecplotConnectivity(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1697
static PeriodicMap NullPeriodicMap
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:976
void WriteTecplotConnectivity(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:402
boost::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:74
void ExtractTracePhys(Array< OneD, NekDouble > &outarray)
Definition: ExpList.h:2207
boost::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
static std::vector< unsigned int > NullUnsignedIntVector
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:2552
void SetWaveSpace(const bool wavespace)
Sets the wave space to the one of the possible configuration true or false.
Definition: ExpList.h:1576
void SetCoeffPhysOffsets()
Definition of the total number of degrees of freedom and quadrature points and offsets to access data...
Definition: ExpList.cpp:247
void SetHomo1DSpecVanVisc(Array< OneD, NekDouble > visc)
This function sets the Spectral Vanishing Viscosity in homogeneous1D expansion.
Definition: ExpList.h:555
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, const bool PhysSpaceForcing)
Definition: ExpList.cpp:2516
SpatialDomains::MeshGraphSharedPtr m_graph
Mesh associated with this expansion list.
Definition: ExpList.h:972
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:3028
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:969
boost::shared_ptr< LibUtilities::SessionReader > GetSession() const
Returns the session object.
Definition: ExpList.h:920
SpatialDomains::MeshGraphSharedPtr GetGraph()
Definition: ExpList.h:931
double NekDouble
void WriteTecplotHeader(std::ostream &outfile, std::string var="")
Definition: ExpList.h:383
const DNekScalBlkMatSharedPtr GenBlockMatrix(const GlobalMatrixKey &gkey)
This function assembles the block diagonal matrix of local matrices of the type mtype.
Definition: ExpList.cpp:778
boost::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
Defines a list of flags.
std::map< int, BoundaryConditionMapShPtr > BoundaryConditionCollection
Definition: Conditions.h:226
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:1673
void GetPeriodicEntities(PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces=NullPeriodicMap)
Definition: ExpList.h:840
virtual Array< OneD, SpatialDomains::BoundaryConditionShPtr > & v_UpdateBndConditions()
Definition: ExpList.cpp:3018
std::map< int, std::vector< PeriodicEntity > > PeriodicMap
virtual void v_ClearGlobalLinSysManager(void)
Definition: ExpList.cpp:2005
static const NekDouble kNekUnsetDouble
Describe a linear system.
void FwdTrans_BndConstrained(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:714
virtual boost::shared_ptr< ExpList > & v_GetTrace()
Definition: ExpList.cpp:2414
void HomogeneousFwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.h:1817
virtual void v_DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.cpp:2572
virtual void v_ExtractCoeffsToCoeffs(const boost::shared_ptr< ExpList > &fromExpList, const Array< OneD, const NekDouble > &fromCoeffs, Array< OneD, NekDouble > &toCoeffs)
Definition: ExpList.cpp:2353
virtual boost::shared_ptr< ExpList > & v_GetPlane(int n)
Definition: ExpList.cpp:3073
Describes a matrix with ordering defined by a local to global map.
const std::vector< bool > & GetLeftAdjacentFaces(void) const
Definition: ExpList.h:2202
static ExpList NullExpList
An empty ExpList object.
Definition: ExpList.h:1476
std::map< int, RobinBCInfoSharedPtr > GetRobinBCInfo()
Definition: ExpList.h:835
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:2037
int GetNpoints(void) const
Returns the total number of quadrature points m_npoints .
Definition: ExpList.h:1567
Array< OneD, SpatialDomains::BoundaryConditionShPtr > & UpdateBndConditions()
Definition: ExpList.h:2228
virtual void v_GetPeriodicEntities(PeriodicMap &periodicVerts, PeriodicMap &periodicEdges, PeriodicMap &periodicFaces)
Definition: ExpList.cpp:3050
virtual void v_PhysInterp1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:1994
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:2106
MultiRegions::Direction const DirCartesianMap[]
Definition: ExpList.h:86
LibUtilities::CommSharedPtr m_comm
Communicator.
Definition: ExpList.h:966
virtual void v_AddTraceIntegral(const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2442
void GetFieldDefinitions(std::vector< LibUtilities::FieldDefinitionsSharedPtr > &fielddef)
Definition: ExpList.h:855
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:2249
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:548
virtual void v_CurlCurl(Array< OneD, Array< OneD, NekDouble > > &Vel, Array< OneD, Array< OneD, NekDouble > > &Q)
Definition: ExpList.cpp:586
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:2540
virtual Array< OneD, const unsigned int > v_GetZIDs(void)
Definition: ExpList.cpp:1977
void ExtractElmtToBndPhys(int i, Array< OneD, NekDouble > &elmt, Array< OneD, NekDouble > &boundary)
Definition: ExpList.h:2299
std::map< GlobalMatrixKey, DNekScalBlkMatSharedPtr > BlockMatrixMap
A map between global matrix keys and their associated block matrices.
Definition: ExpList.h:95
void DealiasedProd(const Array< OneD, NekDouble > &inarray1, const Array< OneD, NekDouble > &inarray2, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1843
bool GetPhysState(void) const
This function indicates whether the array of physical values (implemented as m_phys) is filled or no...
Definition: ExpList.h:1633
static const Array< OneD, ExpListSharedPtr > NullExpListSharedPtrArray
Definition: ExpList.h:2326
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1662
virtual std::map< int, RobinBCInfoSharedPtr > v_GetRobinBCInfo(void)
Definition: ExpList.cpp:3040
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:250
virtual void v_DealiasedDotProd(const Array< OneD, Array< OneD, NekDouble > > &inarray1, const Array< OneD, Array< OneD, NekDouble > > &inarray2, Array< OneD, Array< OneD, NekDouble > > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.cpp:2581
NekDouble GetCoeff(int i)
Get the i th value (coefficient) of m_coeffs.
Definition: ExpList.h:2015
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:2071
void SetCoeff(int i, NekDouble val)
Set the i th coefficiient in m_coeffs to value val.
Definition: ExpList.h:1909
virtual void v_IProductWRTBase_IterPerExp(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:380
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:587
void IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1641
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:1624
void WriteVtkPieceFooter(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1827
void AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
Append the element data listed in elements fielddef->m_ElementIDs onto fielddata. ...
Definition: ExpList.h:864
int GetNcoeffs(void) const
Returns the total number of local degrees of freedom .
Definition: ExpList.h:1485
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:2528
std::vector< int > m_coll_phys_offset
Offset of elemental data into the array m_phys.
Definition: ExpList.h:1044
virtual void v_AppendFieldData(LibUtilities::FieldDefinitionsSharedPtr &fielddef, std::vector< NekDouble > &fielddata)
Definition: ExpList.cpp:2212
static SpatialDomains::BoundaryConditionShPtr GetBoundaryCondition(const SpatialDomains::BoundaryConditionCollection &collection, unsigned int index, const std::string &variable)
Definition: ExpList.cpp:3059
boost::shared_ptr< AssemblyMapDG > & GetTraceMap(void)
Definition: ExpList.h:2148
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:3151
const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & GetBndConditions()
Definition: ExpList.h:2221
void WriteVtkPieceHeader(std::ostream &outfile, int expansion, int istrip=0)
Definition: ExpList.h:411
void MultiplyByBlockMatrix(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:350
void SetCoeffsArray(Array< OneD, NekDouble > &inarray)
Set the m_coeffs array to inarray.
Definition: ExpList.h:1925
virtual void v_MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2507
void DealiasedDotProd(const Array< OneD, Array< OneD, NekDouble > > &inarray1, const Array< OneD, Array< OneD, NekDouble > > &inarray2, Array< OneD, Array< OneD, NekDouble > > &outarray, CoeffState coeffstate=eLocal)
Definition: ExpList.h:1855
int EvalBasisNumModesMax(void) const
Evaulates the maximum number of modes in the elemental basis order over all elements.
Definition: ExpList.h:1499
void GetBoundaryNormals(int i, Array< OneD, Array< OneD, NekDouble > > &normals)
Definition: ExpList.h:2320
bool GetWaveSpace(void) const
This function returns the third direction expansion condition, which can be in wave space (coefficien...
Definition: ExpList.h:1584
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:596
virtual void v_GetBoundaryToElmtMap(Array< OneD, int > &ElmtID, Array< OneD, int > &EdgeID)
Definition: ExpList.cpp:2989
virtual Array< OneD, const unsigned int > v_GetYIDs(void)
Definition: ExpList.cpp:1985
virtual void v_SetUpPhysNormals()
Definition: ExpList.cpp:2809
void GetBCValues(Array< OneD, NekDouble > &BndVals, const Array< OneD, NekDouble > &TotField, int BndID)
Definition: ExpList.h:1867
int GetCoordim(int eid)
This function returns the dimension of the coordinates of the element eid.
Definition: ExpList.h:1898
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:219
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:874
virtual const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & v_GetBndConditions()
Definition: ExpList.cpp:3007
void SetExpType(ExpansionType Type)
Returns the type of the expansion.
Definition: ExpList.cpp:277
void SetPhysArray(Array< OneD, NekDouble > &inarray)
Sets the array m_phys.
Definition: ExpList.h:1615
void WriteTecplotZone(std::ostream &outfile, int expansion=-1)
Definition: ExpList.h:389
std::vector< int > m_coll_coeff_offset
Offset of elemental data into the array m_coeffs.
Definition: ExpList.h:1041
void WriteVtkHeader(std::ostream &outfile)
Definition: ExpList.cpp:1808
virtual void v_ExtractTracePhys(Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2493
virtual void v_PhysGalerkinProjection1DScaled(const NekDouble scale, const Array< OneD, NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:2000
boost::shared_ptr< MeshGraph > MeshGraphSharedPtr
Definition: MeshGraph.h:442
void GlobalToLocal(void)
Scatters from the global coefficients to the local coefficients .
Definition: ExpList.h:1972
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:528
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2707
void GetFwdBwdTracePhys(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Definition: ExpList.h:2187
virtual void v_GetFwdBwdTracePhys(Array< OneD, NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd)
Definition: ExpList.cpp:2468
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:2240
void HomogeneousBwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal, bool Shuff=true, bool UnShuff=true)
Definition: ExpList.h:1830
void SetPhys(int i, NekDouble val)
Set the i th value of m_phys to value val.
Definition: ExpList.h:1590
virtual void v_WriteTecplotField(std::ostream &outfile, int expansion)
Definition: ExpList.cpp:1774
virtual int v_GetNumElmts(void)
Definition: ExpList.h:1110
Collections::CollectionVector m_collections
Definition: ExpList.h:1038
virtual void v_ReadGlobalOptimizationParameters()
Definition: ExpList.cpp:2998
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, const bool PhysSpaceForcing=true)
Solve helmholtz problem.
Definition: ExpList.h:1724
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:228
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:2562
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Definition: ExpList.cpp:2714