Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ContField2D.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ContField2D.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: Field definition in tow-dimensions
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
36 #ifndef NEKTAR_LIBS_MULTIREGIONS_CONTFIELD2D_H
37 #define NEKTAR_LIBS_MULTIREGIONS_CONTFIELD2D_H
38 
47 
48 
49 namespace Nektar
50 {
51  namespace MultiRegions
52  {
53  /// This class is the abstraction of a global continuous two-
54  /// dimensional spectral/hp element expansion which approximates the
55  /// solution of a set of partial differential equations.
57  {
58  public:
59  /// The default constructor.
61 
62  /// This constructor sets up global continuous field based on an
63  /// input mesh and boundary conditions.
67  const std::string &variable = "DefaultVar",
68  const bool DeclareCoeffPhysArrays = true,
69  const bool CheckIfSingularSystem = false,
72 
73  /// Construct a global continuous field with solution type based on
74  /// another field but using a separate input mesh and boundary
75  /// conditions.
76  MULTI_REGIONS_EXPORT ContField2D(const ContField2D &In,
78  const std::string &variable,
79  const bool DeclareCoeffPhysArrays = true,
80  const bool CheckIfSingularSystem = false);
81 
82  /// The copy constructor.
83  MULTI_REGIONS_EXPORT ContField2D(const ContField2D &In, bool DeclareCoeffPhysArrays = true);
84 
85  /// The default destructor.
87 
88  /// Assembles the global coefficients \f$\boldsymbol{\hat{u}}_g\f$
89  /// from the local coefficients \f$\boldsymbol{\hat{u}}_l\f$.
90  inline void Assemble();
91 
92  /// Assembles the global coefficients \f$\boldsymbol{\hat{u}}_g\f$
93  /// from the local coefficients \f$\boldsymbol{\hat{u}}_l\f$.
94  inline void Assemble(
95  const Array<OneD, const NekDouble> &inarray,
96  Array<OneD,NekDouble> &outarray) const;
97 
98  /// Returns the map from local to global level.
100  const;
101 
102 
103  /// Calculates the inner product of a function
104  /// \f$f(\boldsymbol{x})\f$ with respect to all <em>global</em>
105  /// expansion modes \f$\phi_n^e(\boldsymbol{x})\f$.
106  inline void IProductWRTBase(
107  const Array<OneD, const NekDouble> &inarray,
108  Array<OneD, NekDouble> &outarray,
109  CoeffState coeffstate = eLocal);
110 
111  /// Performs the global forward transformation of a function
112  /// \f$f(\boldsymbol{x})\f$, subject to the boundary conditions
113  /// specified.
115  Array<OneD, NekDouble> &outarray,
116  CoeffState coeffstate = eLocal);
117 
118  /// Performs the backward transformation of the spectral/hp
119  /// element expansion.
120  inline void BwdTrans(
121  const Array<OneD, const NekDouble> &inarray,
122  Array<OneD, NekDouble> &outarray,
123  CoeffState coeffstate = eLocal);
124 
125  /// Multiply a solution by the inverse mass matrix.
127  const Array<OneD, const NekDouble> &inarray,
128  Array<OneD, NekDouble> &outarray,
129  CoeffState coeffstate = eLocal);
130 
131  /// Solves the two-dimensional Laplace equation, subject to the
132  /// boundary conditions specified.
134  Array<OneD, NekDouble> &outarray,
135  const Array<OneD, const NekDouble> &dirForcing
138  variablecoeffs = NullNekDoubleArrayofArray,
139  NekDouble time = 0.0,
140  CoeffState coeffstate = eLocal);
141 
142 
143  /// Compute the eigenvalues of the linear advection operator.
145  const NekDouble ay,
150 
151  /// Returns the boundary conditions expansion.
154 
155  /// Returns the boundary conditions.
158 
159  inline int GetGlobalMatrixNnz(const GlobalMatrixKey &gkey);
160 
161  protected:
162 
163  private:
164  /// (A shared pointer to) the object which contains all the
165  /// required information for the transformation from local to
166  /// global degrees of freedom.
168 
169  /// (A shared pointer to) a list which collects all the global
170  /// matrices being assembled, such that they should be constructed
171  /// only once.
173 
174  /// A manager which collects all the global
175  /// linear systems being assembled, such that they should be
176  /// constructed only once.
178 
179  /// Solves the linear system specified by the key \a key.
181  const Array<OneD, const NekDouble> &rhs,
182  Array<OneD, NekDouble> &inout,
183  const Array<OneD, const NekDouble> &dirForcing
185 
186  /// Returns the global matrix specified by \a mkey.
188 
189  /// Returns the linear system specified by the key \a mkey.
191 
193 
194  /// Impose the Dirichlet Boundary Conditions on outarray
196 
198 
199  MULTI_REGIONS_EXPORT virtual void v_FillBndCondFromField(const int nreg);
200 
201  /// Gathers the global coefficients \f$\boldsymbol{\hat{u}}_g\f$
202  /// from the local coefficients \f$\boldsymbol{\hat{u}}_l\f$.
204  const Array<OneD, const NekDouble> &inarray,
205  Array<OneD,NekDouble> &outarray,
206  bool useComm);
207 
208  MULTI_REGIONS_EXPORT virtual void v_LocalToGlobal(bool useComm);
209 
210  /// Scatters from the global coefficients
211  /// \f$\boldsymbol{\hat{u}}_g\f$ to the local coefficients
212  /// \f$\boldsymbol{\hat{u}}_l\f$.
214  const Array<OneD, const NekDouble> &inarray,
215  Array<OneD,NekDouble> &outarray);
216 
217  MULTI_REGIONS_EXPORT virtual void v_GlobalToLocal(void);
218 
219  /// Template method virtual forwarder for FwdTrans().
220  MULTI_REGIONS_EXPORT virtual void v_BwdTrans(
221  const Array<OneD, const NekDouble> &inarray,
222  Array<OneD, NekDouble> &outarray,
223  CoeffState coeffstate);
224 
225 
226  /// Template method virtual forwarder for FwdTrans().
227  MULTI_REGIONS_EXPORT virtual void v_FwdTrans(
228  const Array<OneD, const NekDouble> &inarray,
229  Array<OneD, NekDouble> &outarray,
230  CoeffState coeffstate);
231 
232  /// Template method virtual forwarded for SmoothField().
233  MULTI_REGIONS_EXPORT virtual void v_SmoothField(
234  Array<OneD,NekDouble> &field);
235 
236  /// Template method virtual forwarder for MultiplyByInvMassMatrix().
238  const Array<OneD, const NekDouble> &inarray,
239  Array<OneD, NekDouble> &outarray,
240  CoeffState coeffstate);
241 
242  /// Solves the two-dimensional Helmholtz equation, subject to the
243  /// boundary conditions specified.
244  MULTI_REGIONS_EXPORT virtual void v_HelmSolve(
245  const Array<OneD, const NekDouble> &inarray,
246  Array<OneD, NekDouble> &outarray,
247  const FlagList &flags,
248  const StdRegions::ConstFactorMap &factors,
249  const StdRegions::VarCoeffMap &varcoeff,
250  const Array<OneD, const NekDouble> &dirForcing,
251  const bool PhysSpaceForcing);
252 
253  /// Calculates the result of the multiplication of a global
254  /// matrix of type specified by \a mkey with a vector given by \a
255  /// inarray.
256  virtual void v_GeneralMatrixOp(
257  const GlobalMatrixKey &gkey,
258  const Array<OneD,const NekDouble> &inarray,
259  Array<OneD, NekDouble> &outarray,
260  CoeffState coeffstate);
261 
262  // Solve the linear advection problem assuming that m_coeffs
263  // vector contains an intial estimate for solution
265  const Array<OneD, const NekDouble> &inarray,
266  Array<OneD, NekDouble> &outarray,
267  const NekDouble lambda,
268  CoeffState coeffstate = eLocal,
270  dirForcing = NullNekDouble1DArray);
271 
272  // Solve the linear advection problem assuming that m_coeff
273  // vector contains an intial estimate for solution
275  const Array<OneD, const NekDouble> &inarray,
276  Array<OneD, NekDouble> &outarray,
277  const NekDouble lambda,
278  CoeffState coeffstate = eLocal,
280  dirForcing = NullNekDouble1DArray);
281 
282  /// Template method virtual forwarder for GetBndConditions().
283  MULTI_REGIONS_EXPORT virtual const Array<OneD,const SpatialDomains
286 
287  };
288 
289  typedef boost::shared_ptr<ContField2D> ContField2DSharedPtr;
290 
291  /**
292  * This operation is evaluated as:
293  * \f{tabbing}
294  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
295  * > > Do \= $i=$ $0,N_m^e-1$ \\
296  * > > > $\boldsymbol{\hat{u}}_g[\mbox{map}[e][i]] =
297  * \boldsymbol{\hat{u}}_g[\mbox{map}[e][i]]+\mbox{sign}[e][i] \cdot
298  * \boldsymbol{\hat{u}}^{e}[i]$\\
299  * > > continue\\
300  * > continue
301  * \f}
302  * where \a map\f$[e][i]\f$ is the mapping array and \a
303  * sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
304  * correct modal connectivity between the different elements (both
305  * these arrays are contained in the data member #m_locToGloMap). This
306  * operation is equivalent to the gather operation
307  * \f$\boldsymbol{\hat{u}}_g=\mathcal{A}^{T}\boldsymbol{\hat{u}}_l\f$,
308  * where \f$\mathcal{A}\f$ is the
309  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
310  *
311  * @note The array #m_coeffs should be filled with the local
312  * coefficients \f$\boldsymbol{\hat{u}}_l\f$ and that the
313  * resulting global coefficients \f$\boldsymbol{\hat{u}}_g\f$
314  * will be stored in #m_coeffs.
315  */
316  inline void ContField2D::Assemble()
317  {
318  m_locToGloMap->Assemble(m_coeffs,m_coeffs);
319  }
320 
321  /**
322  * This operation is evaluated as:
323  * \f{tabbing}
324  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
325  * > > Do \= $i=$ $0,N_m^e-1$ \\
326  * > > > $\boldsymbol{\hat{u}}_g[\mbox{map}[e][i]] =
327  * \boldsymbol{\hat{u}}_g[\mbox{map}[e][i]]+\mbox{sign}[e][i] \cdot
328  * \boldsymbol{\hat{u}}^{e}[i]$\\
329  * > > continue\\
330  * > continue
331  * \f}
332  * where \a map\f$[e][i]\f$ is the mapping array and \a
333  * sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
334  * correct modal connectivity between the different elements (both
335  * these arrays are contained in the data member #m_locToGloMap). This
336  * operation is equivalent to the gather operation
337  * \f$\boldsymbol{\hat{u}}_g=\mathcal{A}^{T}\boldsymbol{\hat{u}}_l\f$,
338  * where \f$\mathcal{A}\f$ is the
339  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
340  *
341  * @param inarray An array of size \f$N_\mathrm{eof}\f$
342  * containing the local degrees of freedom
343  * \f$\boldsymbol{x}_l\f$.
344  * @param outarray The resulting global degrees of freedom
345  * \f$\boldsymbol{x}_g\f$ will be stored in this
346  * array of size \f$N_\mathrm{dof}\f$.
347  */
349  const Array<OneD, const NekDouble> &inarray,
350  Array<OneD,NekDouble> &outarray) const
351  {
352  m_locToGloMap->Assemble(inarray,outarray);
353  }
354 
355 
356  inline const AssemblyMapCGSharedPtr&
358  {
359  return m_locToGloMap;
360  }
361 
362 
363  /**
364  * The operation is evaluated locally (i.e. with respect to all local
365  * expansion modes) by the function ExpList#IProductWRTBase. The inner
366  * product with respect to the global expansion modes is than obtained
367  * by a global assembly operation.
368  *
369  * The values of the function \f$f(\boldsymbol{x})\f$ evaluated at the
370  * quadrature points \f$\boldsymbol{x}_i\f$ should be contained in the
371  * variable #m_phys of the ExpList object \a in. The result is stored
372  * in the array #m_coeffs.
373  *
374  * @param In An ExpList, containing the discrete evaluation
375  * of \f$f(\boldsymbol{x})\f$ at the quadrature
376  * points in its array #m_phys.
377  */
379  const Array<OneD, const NekDouble> &inarray,
380  Array<OneD, NekDouble> &outarray,
381  CoeffState coeffstate)
382 
383  {
384  if(coeffstate == eGlobal)
385  {
386  bool doGlobalOp = m_globalOptParam->DoGlobalMatOp(
388 
389  if(doGlobalOp)
390  {
392  m_locToGloMap);
394  mat->Multiply(inarray,outarray);
395  m_locToGloMap->UniversalAssemble(outarray);
396  }
397  else
398  {
400  IProductWRTBase_IterPerExp(inarray,wsp);
401  Assemble(wsp,outarray);
402  }
403  }
404  else
405  {
406  IProductWRTBase_IterPerExp(inarray,outarray);
407  }
408  }
409 
410  /**
411  * Given the coefficients of an expansion, this function evaluates the
412  * spectral/hp expansion \f$u^{\delta}(\boldsymbol{x})\f$ at the
413  * quadrature points \f$\boldsymbol{x}_i\f$. This operation is
414  * evaluated locally by the function ExpList#BwdTrans.
415  *
416  * The coefficients of the expansion should be contained in the variable
417  * #m_coeffs of the ExpList object \a In. The resulting physical values
418  * at the quadrature points \f$u^{\delta}(\boldsymbol{x}_i)\f$ are
419  * stored in the array #m_phys.
420  *
421  * @param In An ExpList, containing the local coefficients
422  * \f$\hat{u}_n^e\f$ in its array #m_coeffs.
423  */
425  const Array<OneD, const NekDouble> &inarray,
426  Array<OneD, NekDouble> &outarray,
427  CoeffState coeffstate)
428  {
429  if(coeffstate == eGlobal)
430  {
431  bool doGlobalOp = m_globalOptParam->DoGlobalMatOp(
433 
434  if(doGlobalOp)
435  {
438  mat->Multiply(inarray,outarray);
439  }
440  else
441  {
443  GlobalToLocal(inarray,wsp);
444  BwdTrans_IterPerExp(wsp,outarray);
445  }
446  }
447  else
448  {
449  BwdTrans_IterPerExp(inarray,outarray);
450  }
451  }
452 
455  {
456  return m_bndCondExpansions;
457  }
458 
461  {
462  return m_bndConditions;
463  }
464 
466  {
468  "To use method must have a AssemblyMap "
469  "attached to key");
470 
471  GlobalMatrixMap::iterator matrixIter = m_globalMat->find(gkey);
472 
473  if(matrixIter == m_globalMat->end())
474  {
475  return 0;
476  }
477  else
478  {
479  return matrixIter->second->GetNumNonZeroEntries();
480  }
481 
482  return 0;
483  }
484 
485  } //end of namespace
486 } //end of namespace
487 
488 #endif // MULTIERGIONS_CONTFIELD2D_H
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Template method virtual forwarder for FwdTrans().
void LaplaceSolve(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray, const Array< OneD, Array< OneD, NekDouble > > &variablecoeffs=NullNekDoubleArrayofArray, NekDouble time=0.0, CoeffState coeffstate=eLocal)
Solves the two-dimensional Laplace equation, subject to the boundary conditions specified.
AssemblyMapCGSharedPtr m_locToGloMap
(A shared pointer to) the object which contains all the required information for the transformation f...
Definition: ContField2D.h:167
Local coefficients.
static Array< OneD, NekDouble > NullNekDouble1DArray
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Performs the global forward transformation of a function , subject to the boundary conditions specifi...
NekOptimize::GlobalOptParamSharedPtr m_globalOptParam
Definition: ExpList.h:1052
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)
const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & GetBndConditions()
Returns the boundary conditions.
Definition: ContField2D.h:460
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Impose the Dirichlet Boundary Conditions on outarray.
GlobalLinSysSharedPtr GenGlobalLinSys(const GlobalLinSysKey &mkey)
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Performs the backward transformation of the spectral/hp element expansion.
Definition: ContField2D.h:424
#define MULTI_REGIONS_EXPORT
boost::shared_ptr< ContField2D > ContField2DSharedPtr
Definition: ContField2D.h:289
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:252
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:998
boost::shared_ptr< GlobalMatrix > GlobalMatrixSharedPtr
Shared pointer to a GlobalMatrix object.
Definition: GlobalMatrix.h:89
LibUtilities::NekManager< GlobalLinSysKey, GlobalLinSys > m_globalLinSysManager
A manager which collects all the global linear systems being assembled, such that they should be cons...
Definition: ContField2D.h:177
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
Global coefficients.
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool useComm)
Gathers the global coefficients from the local coefficients .
Array< OneD, MultiRegions::ExpListSharedPtr > m_bndCondExpansions
An object which contains the discretised boundary conditions.
This class is the abstraction of a global continuous two- dimensional spectral/hp element expansion w...
Definition: ContField2D.h:56
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
boost::shared_ptr< GlobalMatrixMap > GlobalMatrixMapShPtr
Shared pointer to a global matrix map.
Definition: GlobalMatrix.h:93
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:227
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
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Template method virtual forwarder for FwdTrans().
virtual void v_SmoothField(Array< OneD, NekDouble > &field)
Template method virtual forwarded for SmoothField().
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:976
GlobalMatrixMapShPtr m_globalMat
(A shared pointer to) a list which collects all the global matrices being assembled, such that they should be constructed only once.
Definition: ContField2D.h:172
void GlobalSolve(const GlobalLinSysKey &key, const Array< OneD, const NekDouble > &rhs, Array< OneD, NekDouble > &inout, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Solves the linear system specified by the key key.
GlobalMatrixSharedPtr GetGlobalMatrix(const GlobalMatrixKey &mkey)
Returns the global matrix specified by mkey.
virtual ~ContField2D()
The default destructor.
virtual void v_GeneralMatrixOp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Calculates the result of the multiplication of a global matrix of type specified by mkey with a vecto...
ContField2D()
The default constructor.
Definition: ContField2D.cpp:88
virtual void v_GlobalToLocal(void)
void Assemble()
Assembles the global coefficients from the local coefficients .
Definition: ContField2D.h:316
double NekDouble
Defines a list of flags.
Array< OneD, SpatialDomains::BoundaryConditionShPtr > m_bndConditions
An array which contains the information about the boundary condition on the different boundary region...
virtual const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & v_GetBndConditions()
Template method virtual forwarder for GetBndConditions().
Describe a linear system.
void IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Calculates the inner product of a function with respect to all global expansion modes ...
Definition: ContField2D.h:378
Describes a matrix with ordering defined by a local to global map.
const AssemblyMapCGSharedPtr & GetLocalToGlobalMap() const
Returns the map from local to global level.
Definition: ContField2D.h:357
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
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)
Solves the two-dimensional Helmholtz equation, subject to the boundary conditions specified...
bool LocToGloMapIsDefined() const
Returns true if a local to global map is defined.
virtual void v_MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Template method virtual forwarder for MultiplyByInvMassMatrix().
void MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Multiply a solution by the inverse mass matrix.
boost::shared_ptr< GlobalLinSys > GlobalLinSysSharedPtr
Pointer to a GlobalLinSys object.
Definition: GlobalLinSys.h:52
GlobalLinSysSharedPtr GetGlobalLinSys(const GlobalLinSysKey &mkey)
Returns the linear system specified by the key mkey.
void LinearAdvectionEigs(const NekDouble ax, const NekDouble ay, Array< OneD, NekDouble > &Real, Array< OneD, NekDouble > &Imag, Array< OneD, NekDouble > &Evecs=NullNekDouble1DArray)
Compute the eigenvalues of the linear advection operator.
virtual void v_ClearGlobalLinSysManager(void)
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayofArray
boost::shared_ptr< BoundaryConditionBase > BoundaryConditionShPtr
Definition: Conditions.h:219
int GetGlobalMatrixNnz(const GlobalMatrixKey &gkey)
Definition: ContField2D.h:465
boost::shared_ptr< AssemblyMapCG > AssemblyMapCGSharedPtr
Definition: AssemblyMapCG.h:52
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode...
Definition: ErrorUtil.hpp:228
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
const Array< OneD, const MultiRegions::ExpListSharedPtr > & GetBndCondExpansions()
Returns the boundary conditions expansion.
Definition: ContField2D.h:454
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)