Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
ContField1D.cpp
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File ContField1D.cpp
4 //
5 // For more information, please see: http://www.nektar.info
6 //
7 // The MIT License
8 //
9 // Copyright (c) 2006 Division of Applied Mathematics, Brown University (USA),
10 // Department of Aeronautics, Imperial College London (UK), and Scientific
11 // Computing and Imaging Institute, University of Utah (USA).
12 //
13 // License for the specific language governing rights and limitations under
14 // Permission is hereby granted, free of charge, to any person obtaining a
15 // copy of this software and associated documentation files (the "Software"),
16 // to deal in the Software without restriction, including without limitation
17 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
18 // and/or sell copies of the Software, and to permit persons to whom the
19 // Software is furnished to do so, subject to the following conditions:
20 //
21 // The above copyright notice and this permission notice shall be included
22 // in all copies or substantial portions of the Software.
23 //
24 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
25 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
26 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
27 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
28 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
29 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
30 // DEALINGS IN THE SOFTWARE.
31 //
32 // Description: Field definition for 1D domain with boundary conditions
33 //
34 ///////////////////////////////////////////////////////////////////////////////
35 
38 
39 namespace Nektar
40 {
41  namespace MultiRegions
42  {
43  /**
44  * @class ContField1D
45  * As opposed to the class #ContExpList1D, the class #ContField1D is
46  * able to incorporate the boundary conditions imposed to the problem
47  * to be solved. Therefore, the class is equipped with three additional
48  * data members:
49  * - #m_bndCondExpansions
50  * - #m_bndTypes
51  * - #m_bndCondEquations
52  *
53  * The first data structure, #m_bndCondExpansions,
54  * contains the point Expansion on the boundary, #m_bndTypes
55  * stores information about the type of boundary condition on the
56  * different parts of the boundary while #m_bndCondEquations holds the
57  * equation of the imposed boundary conditions.
58  *
59  * Furthermore, in case of Dirichlet boundary conditions,
60  * this class is capable of lifting a known solution satisfying these
61  * boundary conditions. If we denote the unknown solution by
62  * \f$u^{\mathcal{H}}(\boldsymbol{x})\f$ and the known Dirichlet
63  * boundary conditions by \f$u^{\mathcal{D}}(\boldsymbol{x})\f$, the
64  * expansion then can be decomposed as
65  * \f[ u^{\delta}(\boldsymbol{x}_i)=u^{\mathcal{D}}(\boldsymbol{x}_i)+
66  * u^{\mathcal{H}}(\boldsymbol{x}_i)=\sum_{n=0}^{N^{\mathcal{D}}-1}
67  * \hat{u}_n^{\mathcal{D}}\Phi_n(\boldsymbol{x}_i)+
68  * \sum_{n={N^{\mathcal{D}}}}^{N_{\mathrm{dof}}
69  * -1}\hat{u}_n^{\mathcal{H}}
70  * \Phi_n(\boldsymbol{x}_i).\f]
71  * This lifting is accomplished by ordering the known global degrees of
72  * freedom, prescribed by the Dirichlet boundary conditions, first in
73  * the global array \f$\boldsymbol{\hat{u}}\f$, that is,
74  * \f[\boldsymbol{\hat{u}}=\left[ \begin{array}{c}
75  * \boldsymbol{\hat{u}}^{\mathcal{D}}\\
76  * \boldsymbol{\hat{u}}^{\mathcal{H}}
77  * \end{array} \right].\f]
78  * Such kind of expansions are also referred to as continuoous fields.
79  * This class should be used when solving 2D problems using a standard
80  * Galerkin approach.
81  */
82 
83  /**
84  * Constructs an empty 1D continuous field.
85  */
88  m_locToGloMap(),
89  m_globalLinSysManager(
90  boost::bind(&ContField1D::GenGlobalLinSys, this, _1),
91  std::string("GlobalLinSys"))
92  {
93  }
94 
95 
96  /**
97  * Given a mesh \a graph1D, containing information about the domain and
98  * the spectral/hp element expansion, this constructor fills the list
99  * of local expansions #m_exp with the proper expansions, calculates
100  * the total number of quadrature points \f$\boldsymbol{x}_i\f$ and
101  * local expansion coefficients \f$\hat{u}^e_n\f$ and allocates
102  * memory for the arrays #m_coeffs and #m_phys. Furthermore, it
103  * constructs the mapping array (contained in #m_locToGloMap) for the
104  * transformation between local elemental level and global level, it
105  * calculates the total number global expansion coefficients
106  * \f$\hat{u}_n\f$.
107  * The constructor also discretises the boundary conditions, specified
108  * by the argument \a bcs, by expressing them in terms of the
109  * coefficient of the expansion on the boundary.
110  *
111  * @param graph1D A 1D mesh, containing information about the
112  * domain and the spectral/hp element expansion.
113  * @param bcs The boundary conditions.
114  * @param variable An optional parameter to indicate for which
115  * variable the field should be constructed.
116  */
118  const LibUtilities::SessionReaderSharedPtr &pSession,
119  const SpatialDomains::MeshGraphSharedPtr &graph1D,
120  const std::string &variable,
121  const Collections::ImplementationType ImpType):
122  DisContField1D(pSession,graph1D,variable,false,ImpType),
123  m_locToGloMap(),
124  m_globalLinSysManager(
125  boost::bind(&ContField1D::GenGlobalLinSys, this, _1),
126  std::string("GlobalLinSys"))
127  {
128  SpatialDomains::BoundaryConditions bcs(pSession, graph1D);
129 
134  false,
135  variable,
137  }
138 
139 
140  /**
141  * Constructs a 1D continuous field as a copy of an existing field
142  * including all boundary conditions.
143  * @param In Existing continuous field to duplicate.
144  */
146  DisContField1D(In),
147  m_locToGloMap(In.m_locToGloMap),
148  m_globalLinSysManager(
149  boost::bind(&ContField1D::GenGlobalLinSys, this, _1),
150  std::string("GlobalLinSys"))
151  {
152  }
153 
154  /**
155  * Constructs a 1D continuous field as a copy of an existing explist 1D field
156  * and adding all the boundary conditions.
157  * @param In Existing explist1D field .
158  */
160  DisContField1D(In),
161  m_locToGloMap(),
162  m_globalLinSysManager(
163  boost::bind(&ContField1D::GenGlobalLinSys, this, _1),
164  std::string("GlobalLinSys"))
165  {
167  ::AllocateSharedPtr(pSession, m_ncoeffs, In);
168 
169  }
170 
171  /**
172  *
173  */
175  {
176  }
177 
178 
179  /**
180  * Given a function \f$f(x)\f$ defined at the quadrature
181  * points, this function determines the unknown global coefficients
182  * \f$\boldsymbol{\hat{u}}^{\mathcal{H}}\f$ employing a discrete
183  * Galerkin projection from physical space to coefficient
184  * space. The operation is evaluated by the function #GlobalSolve using
185  * the global mass matrix.
186  *
187  * The values of the function \f$f(x)\f$ evaluated at the
188  * quadrature points \f$x_i\f$ should be contained in the
189  * variable #m_phys of the ExpList object \a Sin. The resulting global
190  * coefficients \f$\hat{u}_g\f$ are stored in the array #m_coeffs.
191  *
192  * @param inarray Discrete \f$f(x)\f$.
193  * @param outarray Storage for result.
194  * @param coeffstate
195  */
197  Array<OneD, NekDouble> &outarray,
198  CoeffState coeffstate)
199  {
200  // Inner product of forcing
202  IProductWRTBase(inarray,wsp,eGlobal);
203 
204  // Solve the system
206 
207  GlobalSolve(key,wsp,outarray);
208  if(coeffstate == eLocal)
209  {
210  GlobalToLocal(outarray,outarray);
211  }
212  }
213 
214 
215  /**
216  * Given the coefficients of an expansion, this function evaluates the
217  * spectral/hp expansion \f$u^{\delta}(x)\f$ at the quadrature
218  * points \f$x_i\f$. This operation is evaluated locally by the
219  * function ExpList#BwdTrans.
220  *
221  * The coefficients of the expansion should be contained in the
222  * variable #m_coeffs of the ExpList object \a In. The resulting
223  * physical values at the quadrature points \f$u^{\delta}(x_i)\f$ are
224  * stored in the array #m_phys.
225  *
226  * @param In An ExpList, containing the local
227  * coefficients \f$\hat{u}_n^e\f$ in its array
228  * #m_coeffs.
229  */
231  const Array<OneD, const NekDouble> &inarray,
232  Array<OneD, NekDouble> &outarray,
233  CoeffState coeffstate)
234  {
235  Array<OneD, NekDouble> tmpinarray;
236  if(coeffstate != eLocal)
237  {
238  tmpinarray = Array<OneD, NekDouble>(inarray);
239  GlobalToLocal(inarray,tmpinarray);
240  }
241  else
242  {
243  tmpinarray = inarray;
244  }
245 
246  BwdTrans_IterPerExp(tmpinarray,outarray);
247  }
248 
249 
250  /**
251  *
252  */
254  const Array<OneD, const NekDouble> &inarray,
255  Array<OneD, NekDouble> &outarray,
256  CoeffState coeffstate)
257  {
259  if(coeffstate == eGlobal)
260  {
261  if(inarray.data() == outarray.data())
262  {
263  Array<OneD, NekDouble> tmp(inarray);
264  GlobalSolve(key,tmp,outarray);
265  }
266  else
267  {
268  GlobalSolve(key,inarray,outarray);
269  }
270  }
271  else
272  {
273  Array<OneD, NekDouble> globaltmp(m_ncoeffs,0.0);
274 
275  if(inarray.data() == outarray.data())
276  {
277  Array<OneD,NekDouble> tmp(inarray);
278  Assemble(tmp,outarray);
279  }
280  else
281  {
282  Assemble(inarray,outarray);
283  }
284 
285  GlobalSolve(key,outarray,globaltmp);
286  GlobalToLocal(globaltmp,outarray);
287  }
288  }
289 
290 
291  /**
292  * Given a linear system specified by the key \a key,
293  * \f[\boldsymbol{M}\boldsymbol{\hat{u}}_g=\boldsymbol{\hat{f}},\f]
294  * this function solves this linear system taking into account the
295  * boundary conditions specified in the data member
296  * #m_bndCondExpansions.
297  * Therefore, it adds an array \f$\boldsymbol{\hat{g}}\f$ which
298  * represents the non-zero surface integral resulting from the weak
299  * boundary conditions (e.g. Neumann boundary conditions) to the right
300  * hand side, that is,
301  * \f[\boldsymbol{M}\boldsymbol{\hat{u}}_g=\boldsymbol{\hat{f}}+
302  * \boldsymbol{\hat{g}}.\f]
303  * Furthermore, it lifts the known degrees of freedom which are
304  * prescribed by the Dirichlet boundary conditions. As these known
305  * coefficients \f$\boldsymbol{\hat{u}}^{\mathcal{D}}\f$ are numbered
306  * first in the global coefficient array \f$\boldsymbol{\hat{u}}_g\f$,
307  * the linear system can be decomposed as,
308  * \f[\left[\begin{array}{cc}
309  * \boldsymbol{M}^{\mathcal{DD}}&\boldsymbol{M}^{\mathcal{DH}}\\
310  * \boldsymbol{M}^{\mathcal{HD}}&\boldsymbol{M}^{\mathcal{HH}}
311  * \end{array}\right]
312  * \left[\begin{array}{c}
313  * \boldsymbol{\hat{u}}^{\mathcal{D}}\\
314  * \boldsymbol{\hat{u}}^{\mathcal{H}}
315  * \end{array}\right]=
316  * \left[\begin{array}{c}
317  * \boldsymbol{\hat{f}}^{\mathcal{D}}\\
318  * \boldsymbol{\hat{f}}^{\mathcal{H}}
319  * \end{array}\right]+
320  * \left[\begin{array}{c}
321  * \boldsymbol{\hat{g}}^{\mathcal{D}}\\
322  * \boldsymbol{\hat{g}}^{\mathcal{H}}
323  * \end{array}\right]
324  * \f]
325  * which will then be solved for the unknown coefficients
326  * \f$\boldsymbol{\hat{u}}^{\mathcal{H}}\f$ as,
327  * \f[
328  * \boldsymbol{M}^{\mathcal{HH}}\boldsymbol{\hat{u}}^{\mathcal{H}}
329  * = \boldsymbol{\hat{f}}^{\mathcal{H}}
330  * +\boldsymbol{\hat{g}}^{\mathcal{H}}
331  * -\boldsymbol{M}^{\mathcal{HD}}\boldsymbol{\hat{u}}^{\mathcal{D}}\f]
332  *
333  * @param key Specifes the linear system to solve.
334  * @param rhs Forcing term \f$\boldsymbol{f}\f$.
335  * @param inout Solution vector \f$\boldsymbol{\hat{u}}\f$.
336  * @param dirForcing .
337  */
339  const Array<OneD, const NekDouble>& rhs,
340  Array<OneD, NekDouble>& inout,
341  const Array<OneD, const NekDouble>& dirForcing)
342  {
343  int NumDirBcs = m_locToGloMap->GetNumGlobalDirBndCoeffs();
344  int contNcoeffs = m_locToGloMap->GetNumGlobalCoeffs();
345 
346  // STEP 1: SET THE DIRICHLET DOFS TO THE RIGHT VALUE
347  // IN THE SOLUTION ARRAY
349 
350  // STEP 2: CALCULATE THE HOMOGENEOUS COEFFICIENTS
351  if(contNcoeffs - NumDirBcs > 0)
352  {
354  LinSys->Solve(rhs,inout,m_locToGloMap,dirForcing);
355  }
356  }
357 
358 
359  /**
360  * The function searches the map #m_globalLinSys to see if the global
361  * matrix has been created before. If not, it calls the function
362  * #GenglobalLinSys to generate the requested global system.
363  *
364  * @param mkey Key specifying the linear system.
365  * @returns Pointer to the required linear system.
366  */
368  const GlobalLinSysKey &mkey)
369  {
370  return m_globalLinSysManager[mkey];
371  }
372 
374  const GlobalLinSysKey &mkey)
375  {
377  "To use method must have a AssemblyMap "
378  "attached to key");
380  }
381 
382 
383  /**
384  * The operation is evaluated locally (i.e. with respect to all local
385  * expansion modes) by the function ExpList#IProductWRTBase. The inner
386  * product with respect to the global expansion modes is than obtained
387  * by a global assembly operation.
388  *
389  * The values of the function \f$f(x)\f$ evaluated at the quadrature
390  * points \f$x_i\f$ should be contained in the variable #m_phys of the
391  * ExpList object \a in. The result is stored in the array
392  * #m_coeffs.
393  *
394  * @param In An ExpList, containing the discrete evaluation
395  * of \f$f(x)\f$ at the quadrature points in its
396  * array #m_phys.
397  */
399  const Array<OneD, const NekDouble> &inarray,
400  Array<OneD, NekDouble> &outarray,
401  CoeffState coeffstate)
402  {
403  if(coeffstate == eGlobal)
404  {
406  IProductWRTBase_IterPerExp(inarray,wsp);
407  Assemble(wsp,outarray);
408  }
409  else
410  {
411  IProductWRTBase_IterPerExp(inarray,outarray);
412  }
413  }
414 
415 
417  const Array<OneD, const NekDouble> &inarray,
418  Array<OneD, NekDouble> &outarray,
419  CoeffState coeffstate)
420  {
421  FwdTrans(inarray,outarray,coeffstate);
422  }
423 
425  const Array<OneD, const NekDouble> &inarray,
426  Array<OneD, NekDouble> &outarray,
427  CoeffState coeffstate)
428  {
429  MultiplyByInvMassMatrix(inarray,outarray,coeffstate);
430  }
431 
433  {
434  for(int i = 0; i < m_bndCondExpansions.num_elements(); ++i)
435  {
436  if(m_bndConditions[i]->GetBoundaryConditionType() == SpatialDomains::eDirichlet)
437  {
438  outarray[m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap(i)]
439  = m_bndCondExpansions[i]->GetCoeff(0);
440  }
441  }
442  }
443 
444  /**
445  * This operation is evaluated as:
446  * \f{tabbing}
447  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
448  * > > Do \= $i=$ $0,N_m^e-1$ \\
449  * > > > $\boldsymbol{\hat{u}}_g[\mbox{map}[e][i]] =
450  * \mbox{sign}[e][i] \cdot \boldsymbol{\hat{u}}^{e}[i]$\\
451  * > > continue\\
452  * > continue
453  * \f}
454  * where \a map\f$[e][i]\f$ is the mapping array and \a
455  * sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
456  * correct modal connectivity between the different elements (both
457  * these arrays are contained in the data member #m_locToGloMap). This
458  * operation is equivalent to the gather operation
459  * \f$\boldsymbol{\hat{u}}_g=\mathcal{A}^{-1}\boldsymbol{\hat{u}}_l\f$,
460  * where \f$\mathcal{A}\f$ is the
461  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
462  *
463  */
465  const Array<OneD, const NekDouble> &inarray,
466  Array<OneD,NekDouble> &outarray, bool useComm)
467  {
468  m_locToGloMap->LocalToGlobal(inarray, outarray, useComm);
469  }
470 
471 
472  void ContField1D::v_LocalToGlobal(bool useComm)
473  {
474  m_locToGloMap->LocalToGlobal(m_coeffs,m_coeffs, useComm);
475  }
476 
477  /**
478  * This operation is evaluated as:
479  * \f{tabbing}
480  * \hspace{1cm} \= Do \= $e=$ $1, N_{\mathrm{el}}$ \\
481  * > > Do \= $i=$ $0,N_m^e-1$ \\
482  * > > > $\boldsymbol{\hat{u}}^{e}[i] = \mbox{sign}[e][i] \cdot
483  * \boldsymbol{\hat{u}}_g[\mbox{map}[e][i]]$ \\
484  * > > continue \\
485  * > continue
486  * \f}
487  * where \a map\f$[e][i]\f$ is the mapping array and
488  * \a sign\f$[e][i]\f$ is an array of similar dimensions ensuring the
489  * correct modal connectivity between the different elements (both
490  * these arrays are contained in the data member #m_locToGloMap). This
491  * operation is equivalent to the scatter operation
492  * \f$\boldsymbol{\hat{u}}_l=\mathcal{A}\boldsymbol{\hat{u}}_g\f$, where
493  * \f$\mathcal{A}\f$ is the
494  * \f$N_{\mathrm{eof}}\times N_{\mathrm{dof}}\f$ permutation matrix.
495  *
496  */
498  const Array<OneD, const NekDouble> &inarray,
499  Array<OneD,NekDouble> &outarray)
500  {
501  m_locToGloMap->GlobalToLocal(inarray, outarray);
502  }
503 
505  {
506  m_locToGloMap->GlobalToLocal(m_coeffs,m_coeffs);
507  }
508 
509  /**
510  * Consider the one dimensional Helmholtz equation,
511  * \f[\frac{d^2u}{dx^2}-\lambda u(x) = f(x),\f]
512  * supplemented with appropriate boundary conditions (which are
513  * contained in the data member #m_bndCondExpansions). Applying a
514  * \f$C^0\f$ continuous Galerkin discretisation, this equation leads to
515  * the following linear system:
516  * \f[\left( \boldsymbol{M}+\lambda\boldsymbol{L}\right)
517  * \boldsymbol{\hat{u}}_g=\boldsymbol{\hat{f}}\f]
518  * where \f$\boldsymbol{M}\f$ and \f$\boldsymbol{L}\f$ are the mass and
519  * Laplacian matrix respectively. This function solves the system above
520  * for the global coefficients \f$\boldsymbol{\hat{u}}\f$ by a call to
521  * the function #GlobalSolve.
522  *
523  * The values of the function \f$f(x)\f$ evaluated at the
524  * quadrature points \f$\boldsymbol{x}_i\f$ should be contained in the
525  * variable #m_phys of the ExpList object \a inarray. The resulting
526  * global coefficients \f$\boldsymbol{\hat{u}}_g\f$ are stored in the
527  * array #m_coeffs.
528  *
529  * @param inarray Input containing forcing function
530  * \f$\boldsymbol{f}\f$ at the quadrature points.
531  * @param outarray Output containing the coefficients
532  * \f$\boldsymbol{u}_g\f$
533  * @param lambda Parameter value.
534  * @param Sigma Coefficients of lambda.
535  * @param varcoeff Variable diffusivity coefficients.
536  * @param coeffstate
537  * @param dirForcing Dirichlet Forcing.
538  */
540  const Array<OneD, const NekDouble> &inarray,
541  Array<OneD, NekDouble> &outarray,
542  const FlagList &flags,
543  const StdRegions::ConstFactorMap &factors,
544  const StdRegions::VarCoeffMap &varcoeff,
545  const Array<OneD, const NekDouble> &dirForcing,
546  const bool PhysSpaceForcing)
547  {
548  // Inner product of forcing
549  int contNcoeffs = m_locToGloMap->GetNumGlobalCoeffs();
550  Array<OneD,NekDouble> wsp(contNcoeffs);
551  if(PhysSpaceForcing)
552  {
553  IProductWRTBase(inarray,wsp,eGlobal);
554  }
555  else
556  {
557  Assemble(inarray,wsp);
558  }
559  // Note -1.0 term necessary to invert forcing function to
560  // be consistent with matrix definition
561  Vmath::Neg(contNcoeffs, wsp, 1);
562 
563  // Forcing function with weak boundary conditions
564  int i;
565  for(i = 0; i < m_bndCondExpansions.num_elements(); ++i)
566  {
567  if(m_bndConditions[i]->GetBoundaryConditionType() != SpatialDomains::eDirichlet)
568  {
569  wsp[m_locToGloMap->GetBndCondCoeffsToGlobalCoeffsMap(i)]
570  += m_bndCondExpansions[i]->GetCoeff(0);
571  }
572  }
573 
574  // Solve the system
576  m_locToGloMap,factors,varcoeff);
577 
578  if(flags.isSet(eUseGlobal))
579  {
580  GlobalSolve(key,wsp,outarray,dirForcing);
581  }
582  else
583  {
584  Array<OneD,NekDouble> tmp(contNcoeffs,0.0);
585  GlobalSolve(key,wsp,tmp,dirForcing);
586  GlobalToLocal(tmp,outarray);
587  }
588  }
589 
592  {
593  return GetBndConditions();
594  }
595 
597  const Array<OneD, const NekDouble> &inarray,
598  Array<OneD, NekDouble> &outarray,
599  CoeffState coeffstate)
600  {
601  BwdTrans(inarray,outarray,coeffstate);
602  }
603 
605  const Array<OneD, const NekDouble> &inarray,
606  Array<OneD, NekDouble> &outarray,
607  CoeffState coeffstate)
608  {
609  IProductWRTBase(inarray,outarray,coeffstate);
610  }
611 
612  /**
613  * This is equivalent to the operation:
614  * \f[\boldsymbol{M\hat{u}}_g\f]
615  * where \f$\boldsymbol{M}\f$ is the global matrix of type specified by
616  * \a mkey. After scattering the global array \a inarray to local
617  * level, this operation is evaluated locally by the function
618  * ExpList#GeneralMatrixOp. The global result is then obtained by a
619  * global assembly procedure.
620  *
621  * @param mkey This key uniquely defines the type matrix
622  * required for the operation.
623  * @param inarray The vector \f$\boldsymbol{\hat{u}}_g\f$ of size
624  * \f$N_{\mathrm{dof}}\f$.
625  * @param outarray The resulting vector of size
626  * \f$N_{\mathrm{dof}}\f$.
627  */
629  const GlobalMatrixKey &gkey,
630  const Array<OneD,const NekDouble> &inarray,
631  Array<OneD, NekDouble> &outarray,
632  CoeffState coeffstate)
633  {
634  if(coeffstate == eGlobal)
635  {
637  Array<OneD,NekDouble> tmp2(tmp1+m_ncoeffs);
638  GlobalToLocal(inarray,tmp1);
639  GeneralMatrixOp_IterPerExp(gkey,tmp1,tmp2);
640  Assemble(tmp2,outarray);
641  }
642  else
643  {
644  GeneralMatrixOp_IterPerExp(gkey,inarray,outarray);
645  }
646  }
647 
648 
649 
650  /**
651  * Reset the GlobalLinSys Manager
652  */
654  {
655  m_globalLinSysManager.ClearManager("GlobalLinSys");
656  }
657 
658  } // end of namespace
659 } //end of namespace
void GeneralMatrixOp_IterPerExp(const GlobalMatrixKey &gkey, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: ExpList.cpp:929
Local coefficients.
virtual void v_MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
virtual void v_BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
Array< OneD, SpatialDomains::BoundaryConditionShPtr > m_bndConditions
An array which contains the information about the boundary condition on the different boundary region...
STL namespace.
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:252
virtual void v_FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
Perform a forward transform.
Array< OneD, NekDouble > m_coeffs
Concatenation of all local expansion coefficients.
Definition: ExpList.h:998
void MultiplyByInvMassMatrix(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
boost::shared_ptr< SessionReader > SessionReaderSharedPtr
Definition: MeshPartition.h:51
virtual void v_IProductWRTBase(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate)
GlobalLinSysSharedPtr GetGlobalLinSys(const GlobalLinSysKey &mkey)
Returns the linear system specified by mkey.
Global coefficients.
bool isSet(const FlagType &key) const
Abstraction of a global continuous one-dimensional spectral/hp element expansion which approximates t...
Definition: ContField1D.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
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 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
int m_ncoeffs
The total number of local degrees of freedom. m_ncoeffs .
Definition: ExpList.h:976
LibUtilities::NekManager< GlobalLinSysKey, GlobalLinSys > m_globalLinSysManager
A manager which collects all the global linear systems being assembled, such that they should be cons...
Definition: ContField1D.h:156
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 ...
This class is the abstraction of a global discontinuous two- dimensional spectral/hp element expansio...
LibUtilities::SessionReaderSharedPtr m_session
Session.
Definition: ExpList.h:969
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:396
Defines a list of flags.
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)
Describe a linear system.
virtual void v_ImposeDirichletConditions(Array< OneD, NekDouble > &outarray)
Impose the Dirichlet Boundary Conditions on outarray.
Describes a matrix with ordering defined by a local to global map.
const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & GetBndConditions()
Definition: ContField1D.h:249
void GlobalSolve(const GlobalLinSysKey &key, const Array< OneD, const NekDouble > &rhs, Array< OneD, NekDouble > &inout, const Array< OneD, const NekDouble > &dirForcing=NullNekDouble1DArray)
Solve the linear system specified by the key key.
virtual void v_ClearGlobalLinSysManager(void)
virtual void v_GlobalToLocal(void)
ContField1D()
Default constructor.
Definition: ContField1D.cpp:86
Array< OneD, MultiRegions::ExpListSharedPtr > m_bndCondExpansions
Discretised boundary conditions.
bool LocToGloMapIsDefined() const
Returns true if a local to global map is defined.
void Assemble()
Assembles the global coefficients from the local coefficients .
Definition: ContField1D.h:276
void BwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
This function performs the backward transformation of the spectral/hp element expansion.
This class is the abstraction of a one-dimensional multi-elemental expansions which is merely a colle...
Definition: ExpList1D.h:61
boost::shared_ptr< GlobalLinSys > GlobalLinSysSharedPtr
Pointer to a GlobalLinSys object.
Definition: GlobalLinSys.h:52
virtual ~ContField1D()
Destructor.
virtual void v_LocalToGlobal(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, bool useComm)
Gathers the global coefficients from the local coefficients .
GlobalLinSysSharedPtr GenGlobalLinSys(const GlobalLinSysKey &mkey)
void FwdTrans(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, CoeffState coeffstate=eLocal)
Perform global forward transformation of a function ,.
PeriodicMap m_periodicVerts
A map which identifies groups of periodic vertices.
#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
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...
AssemblyMapCGSharedPtr m_locToGloMap
(A shared pointer to) the object which contains all the required information for the transformation f...
Definition: ContField1D.h:141
virtual const Array< OneD, const SpatialDomains::BoundaryConditionShPtr > & v_GetBndConditions()