Nektar++
Expansion.h
Go to the documentation of this file.
1 ///////////////////////////////////////////////////////////////////////////////
2 //
3 // File Expansion.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 // Permission is hereby granted, free of charge, to any person obtaining a
14 // copy of this software and associated documentation files (the "Software"),
15 // to deal in the Software without restriction, including without limitation
16 // the rights to use, copy, modify, merge, publish, distribute, sublicense,
17 // and/or sell copies of the Software, and to permit persons to whom the
18 // Software is furnished to do so, subject to the following conditions:
19 //
20 // The above copyright notice and this permission notice shall be included
21 // in all copies or substantial portions of the Software.
22 //
23 // THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
24 // OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
25 // FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
26 // THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
27 // LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
28 // FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
29 // DEALINGS IN THE SOFTWARE.
30 //
31 // Description: Header file for Expansion routines
32 //
33 ///////////////////////////////////////////////////////////////////////////////
34 
35 #ifndef EXPANSION_H
36 #define EXPANSION_H
37 
43 #include <map>
44 #include <memory>
45 #include <vector>
46 
47 namespace Nektar
48 {
49 namespace LocalRegions
50 {
51 
52 class Expansion;
53 class MatrixKey;
54 
56 
58 {
66 };
67 
68 typedef std::shared_ptr<Expansion> ExpansionSharedPtr;
69 typedef std::weak_ptr<Expansion> ExpansionWeakPtr;
70 typedef std::vector<ExpansionSharedPtr> ExpansionVector;
71 typedef std::map<MetricType, Array<OneD, NekDouble>> MetricMap;
72 
73 class Expansion : virtual public StdRegions::StdExpansion
74 {
75 public:
77  SpatialDomains::GeometrySharedPtr pGeom); // default constructor.
78  LOCAL_REGIONS_EXPORT Expansion(const Expansion &pSrc); // copy constructor.
80 
81  LOCAL_REGIONS_EXPORT void SetTraceExp(const int traceid,
84 
87 
89  const StdRegions::MatrixType mtype,
90  const StdRegions::ConstFactorMap &factors =
93 
95 
97 
99  CreateIndexMap(const IndexMapKey &ikey);
100 
102  CreateStaticCondMatrix(const MatrixKey &mkey);
103 
105  const;
106 
109  const StdRegions::MatrixType matrixType);
110 
113 
115  const NekDouble *data, const std::vector<unsigned int> &nummodes,
116  const int nmodes_offset, NekDouble *coeffs,
117  std::vector<LibUtilities::BasisType> &fromType);
118 
120  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
123  Array<OneD, NekDouble> &outarray);
125  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
127  Array<OneD, NekDouble> &outarray);
129  const int face, const std::shared_ptr<Expansion> &FaceExp,
131  Array<OneD, NekDouble> &outarray);
133  const int dir, const Array<OneD, const NekDouble> &inarray,
136  Array<OneD, NekDouble> &outarray);
139 
141  Array<OneD, Array<OneD, NekDouble>> &factors,
142  Array<OneD, Array<OneD, NekDouble>> &d0factors,
143  Array<OneD, Array<OneD, NekDouble>> &d1factors);
144 
146  {
147  return m_indexMapManager[ikey];
148  }
149 
151  const int dir, const Array<OneD, const NekDouble> &inarray,
152  Array<OneD, Array<OneD, NekDouble>> &outarray)
153  {
154  v_AlignVectorToCollapsedDir(dir, inarray, outarray);
155  }
156 
158 
160 
161  inline int GetLeftAdjacentElementTrace() const;
162 
163  inline int GetRightAdjacentElementTrace() const;
164 
165  inline void SetAdjacentElementExp(int traceid, ExpansionSharedPtr &e);
166 
168  {
169  return v_GetTraceOrient(trace);
170  }
171 
174  Array<OneD, NekDouble> &outarray)
175  {
176  v_SetCoeffsToOrientation(dir, inarray, outarray);
177  }
178 
179  /// Divided by the metric jacobi and quadrature weights
181  const Array<OneD, const NekDouble> &inarray,
182  Array<OneD, NekDouble> &outarray)
183  {
184  v_DivideByQuadratureMetric(inarray, outarray);
185  }
186 
187  /**
188  * @brief Extract the metric factors to compute the contravariant
189  * fluxes along edge \a edge and stores them into \a outarray
190  * following the local edge orientation (i.e. anticlockwise
191  * convention).
192  */
193  inline void GetTraceQFactors(const int trace,
194  Array<OneD, NekDouble> &outarray)
195  {
196  v_GetTraceQFactors(trace, outarray);
197  }
198 
199  inline void GetTracePhysVals(
200  const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp,
201  const Array<OneD, const NekDouble> &inarray,
202  Array<OneD, NekDouble> &outarray,
204  {
205  v_GetTracePhysVals(trace, TraceExp, inarray, outarray, orient);
206  }
207 
208  inline void GetTracePhysMap(const int edge, Array<OneD, int> &outarray)
209  {
210  v_GetTracePhysMap(edge, outarray);
211  }
212 
214  Array<OneD, int> &idmap, const int nq0,
215  const int nq1)
216  {
217  v_ReOrientTracePhysMap(orient, idmap, nq0, nq1);
218  }
219 
220  LOCAL_REGIONS_EXPORT const NormalVector &GetTraceNormal(const int id);
221 
222  inline void ComputeTraceNormal(const int id)
223  {
225  }
226 
228  {
229  return v_GetPhysNormals();
230  }
231 
233  {
234  v_SetPhysNormals(normal);
235  }
236 
237  inline void SetUpPhysNormals(const int trace)
238  {
239  v_SetUpPhysNormals(trace);
240  }
241 
242  inline void AddRobinMassMatrix(
243  const int traceid, const Array<OneD, const NekDouble> &primCoeffs,
244  DNekMatSharedPtr &inoutmat)
245  {
246  v_AddRobinMassMatrix(traceid, primCoeffs, inoutmat);
247  }
248 
249  inline void TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
250  {
251  v_TraceNormLen(traceid, h, p);
252  }
253 
255  const int traceid, const Array<OneD, const NekDouble> &primCoeffs,
256  const Array<OneD, NekDouble> &incoeffs, Array<OneD, NekDouble> &coeffs)
257  {
258  v_AddRobinTraceContribution(traceid, primCoeffs, incoeffs, coeffs);
259  }
260 
262  &GetElmtBndNormDirElmtLen(const int nbnd) const;
263 
266 
267 protected:
270 
271  std::map<int, ExpansionWeakPtr> m_traceExp;
275  std::map<int, NormalVector> m_traceNormals;
280 
281  /// the element length in each element boundary(Vertex, edge
282  /// or face) normal direction calculated based on the local
283  /// m_metricinfo times the standard element length (which is
284  /// 2.0)
285  std::map<int, Array<OneD, NekDouble>> m_elmtBndNormDirElmtLen;
286 
287  void ComputeLaplacianMetric();
290  const Array<OneD, const NekDouble> &direction,
292 
293  virtual void v_MultiplyByQuadratureMetric(
294  const Array<OneD, const NekDouble> &inarray,
295  Array<OneD, NekDouble> &outarray);
296 
297  virtual void v_DivideByQuadratureMetric(
298  const Array<OneD, const NekDouble> &inarray,
299  Array<OneD, NekDouble> &outarray);
300 
302  {
303  }
304 
305  virtual void v_GetCoords(Array<OneD, NekDouble> &coords_1,
306  Array<OneD, NekDouble> &coords_2,
307  Array<OneD, NekDouble> &coords_3);
308 
309  Array<OneD, NekDouble> v_GetMF(const int dir, const int shapedim,
310  const StdRegions::VarCoeffMap &varcoeffs);
311 
312  Array<OneD, NekDouble> v_GetMFDiv(const int dir,
313  const StdRegions::VarCoeffMap &varcoeffs);
314 
315  Array<OneD, NekDouble> v_GetMFMag(const int dir,
316  const StdRegions::VarCoeffMap &varcoeffs);
317 
319  const LocalRegions::MatrixKey &mkey);
320 
322  const DNekScalMatSharedPtr &r_bnd,
323  const StdRegions::MatrixType matrixType);
324 
326  const DNekScalMatSharedPtr &r_bnd);
327 
328  virtual void v_ExtractDataToCoeffs(
329  const NekDouble *data, const std::vector<unsigned int> &nummodes,
330  const int nmodes_offset, NekDouble *coeffs,
331  std::vector<LibUtilities::BasisType> &fromType);
332 
333  virtual void v_AddEdgeNormBoundaryInt(
334  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
337  Array<OneD, NekDouble> &outarray);
338  virtual void v_AddEdgeNormBoundaryInt(
339  const int edge, const std::shared_ptr<Expansion> &EdgeExp,
341  Array<OneD, NekDouble> &outarray);
342  virtual void v_AddFaceNormBoundaryInt(
343  const int face, const std::shared_ptr<Expansion> &FaceExp,
345  Array<OneD, NekDouble> &outarray);
346  virtual void v_DGDeriv(const int dir,
347  const Array<OneD, const NekDouble> &inarray,
350  Array<OneD, NekDouble> &outarray);
351  virtual NekDouble v_VectorFlux(
352  const Array<OneD, Array<OneD, NekDouble>> &vec);
353 
354  virtual void v_NormalTraceDerivFactors(
355  Array<OneD, Array<OneD, NekDouble>> &factors,
356  Array<OneD, Array<OneD, NekDouble>> &d0factors,
357  Array<OneD, Array<OneD, NekDouble>> &d1factors);
358 
359  virtual void v_AlignVectorToCollapsedDir(
360  const int dir, const Array<OneD, const NekDouble> &inarray,
361  Array<OneD, Array<OneD, NekDouble>> &outarray);
362 
363  virtual StdRegions::Orientation v_GetTraceOrient(int trace);
364 
367  Array<OneD, NekDouble> &outarray);
368 
369  virtual void v_GetTraceQFactors(const int trace,
370  Array<OneD, NekDouble> &outarray);
371 
372  virtual void v_GetTracePhysVals(
373  const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp,
374  const Array<OneD, const NekDouble> &inarray,
376 
377  virtual void v_GetTracePhysMap(const int edge, Array<OneD, int> &outarray);
378 
379  virtual void v_ReOrientTracePhysMap(const StdRegions::Orientation orient,
380  Array<OneD, int> &idmap, const int nq0,
381  const int nq1 = -1);
382 
383  virtual void v_ComputeTraceNormal(const int id);
384 
385  virtual const Array<OneD, const NekDouble> &v_GetPhysNormals(void);
386 
387  virtual void v_SetPhysNormals(Array<OneD, const NekDouble> &normal);
388 
389  virtual void v_SetUpPhysNormals(const int id);
390 
391  virtual void v_AddRobinMassMatrix(
392  const int face, const Array<OneD, const NekDouble> &primCoeffs,
393  DNekMatSharedPtr &inoutmat);
394 
395  virtual void v_AddRobinTraceContribution(
396  const int traceid, const Array<OneD, const NekDouble> &primCoeffs,
397  const Array<OneD, NekDouble> &incoeffs, Array<OneD, NekDouble> &coeffs);
398 
399  virtual void v_TraceNormLen(const int traceid, NekDouble &h, NekDouble &p);
400 
401  virtual void v_GenTraceExp(const int traceid, ExpansionSharedPtr &exp);
402 
403 private:
404 };
405 
407 {
408  ASSERTL1(traceid < GetNtraces(), "Trace is out of range.");
409 
410  ExpansionSharedPtr returnval;
411 
412  if (m_traceExp.count(traceid))
413  {
414  // Use stored value
415  returnval = m_traceExp[traceid].lock();
416  }
417  else
418  {
419  // Generate trace exp
420  v_GenTraceExp(traceid, returnval);
421  }
422 
423  return returnval;
424 }
425 
426 inline void Expansion::SetTraceExp(const int traceid, ExpansionSharedPtr &exp)
427 {
428  ASSERTL1(traceid < GetNtraces(), "Trace out of range.");
429 
430  m_traceExp[traceid] = exp;
431 }
432 
434 {
435  ASSERTL1(m_elementLeft.lock().get(), "Left adjacent element not set.");
436  return m_elementLeft.lock();
437 }
438 
440 {
441  ASSERTL1(m_elementLeft.lock().get(), "Right adjacent element not set.");
442 
443  return m_elementRight.lock();
444 }
445 
447 {
448  return m_elementTraceLeft;
449 }
450 
452 {
453  return m_elementTraceRight;
454 }
455 
456 inline void Expansion::SetAdjacentElementExp(int traceid,
457  ExpansionSharedPtr &exp)
458 {
459  if (m_elementLeft.lock().get())
460  {
461  m_elementRight = exp;
462  m_elementTraceRight = traceid;
463  }
464  else
465  {
466  m_elementLeft = exp;
467  m_elementTraceLeft = traceid;
468  }
469 }
470 
471 } // namespace LocalRegions
472 } // namespace Nektar
473 
474 #endif
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249
#define LOCAL_REGIONS_EXPORT
std::map< int, NormalVector > m_traceNormals
Definition: Expansion.h:275
std::map< int, Array< OneD, NekDouble > > m_elmtBndNormDirElmtLen
the element length in each element boundary(Vertex, edge or face) normal direction calculated based o...
Definition: Expansion.h:285
virtual void v_SetPhysNormals(Array< OneD, const NekDouble > &normal)
Definition: Expansion.cpp:873
SpatialDomains::GeometrySharedPtr GetGeom() const
Definition: Expansion.cpp:166
void SetTraceExp(const int traceid, ExpansionSharedPtr &f)
Definition: Expansion.h:426
virtual void v_TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
Definition: Expansion.cpp:901
IndexMapValuesSharedPtr CreateIndexMap(const IndexMapKey &ikey)
Definition: Expansion.cpp:180
SpatialDomains::GeometrySharedPtr m_geom
Definition: Expansion.h:272
Array< OneD, NekDouble > v_GetMF(const int dir, const int shapedim, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:627
void GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.h:208
const SpatialDomains::GeomFactorsSharedPtr & GetMetricInfo() const
Definition: Expansion.cpp:245
ExpansionSharedPtr GetLeftAdjacentElementExp() const
Definition: Expansion.h:433
DNekMatSharedPtr BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:100
void ComputeGmatcdotMF(const Array< TwoD, const NekDouble > &df, const Array< OneD, const NekDouble > &direction, Array< OneD, Array< OneD, NekDouble >> &dfdir)
Definition: Expansion.cpp:597
virtual void v_DivideByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:428
virtual void v_GenTraceExp(const int traceid, ExpansionSharedPtr &exp)
Definition: Expansion.cpp:854
NekDouble VectorFlux(const Array< OneD, Array< OneD, NekDouble >> &vec)
Definition: Expansion.cpp:144
virtual void v_AddRobinTraceContribution(const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
Definition: Expansion.cpp:893
ExpansionWeakPtr m_elementRight
Definition: Expansion.h:277
virtual void v_AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Expansion.cpp:916
virtual void v_GetTraceQFactors(const int trace, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:820
virtual void v_ComputeTraceNormal(const int id)
Definition: Expansion.cpp:861
void GetTraceQFactors(const int trace, Array< OneD, NekDouble > &outarray)
Extract the metric factors to compute the contravariant fluxes along edge edge and stores them into o...
Definition: Expansion.h:193
DNekScalBlkMatSharedPtr CreateStaticCondMatrix(const MatrixKey &mkey)
Definition: Expansion.cpp:272
virtual void v_DGDeriv(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble >> &coeffs, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:775
void SetUpPhysNormals(const int trace)
Definition: Expansion.h:237
void SetAdjacentElementExp(int traceid, ExpansionSharedPtr &e)
Definition: Expansion.h:456
virtual DNekMatSharedPtr v_BuildVertexMatrix(const DNekScalMatSharedPtr &r_bnd)
Definition: Expansion.cpp:733
virtual void v_SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:812
virtual DNekScalMatSharedPtr v_GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:264
int GetLeftAdjacentElementTrace() const
Definition: Expansion.h:446
void DivideByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Divided by the metric jacobi and quadrature weights.
Definition: Expansion.h:180
void GetTracePhysVals(const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient=StdRegions::eNoOrientation)
Definition: Expansion.h:199
virtual void v_AddFaceNormBoundaryInt(const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:767
void NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble >> &factors, Array< OneD, Array< OneD, NekDouble >> &d0factors, Array< OneD, Array< OneD, NekDouble >> &d1factors)
Definition: Expansion.cpp:149
Array< OneD, NekDouble > v_GetMFMag(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:703
virtual DNekMatSharedPtr v_BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:725
const Array< OneD, const NekDouble > & GetElmtBndNormDirElmtLen(const int nbnd) const
Definition: Expansion.cpp:907
virtual void v_AddRobinMassMatrix(const int face, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
Definition: Expansion.cpp:885
std::map< int, ExpansionWeakPtr > m_traceExp
Definition: Expansion.h:271
Array< OneD, NekDouble > v_GetMFDiv(const int dir, const StdRegions::VarCoeffMap &varcoeffs)
Definition: Expansion.cpp:680
void AddEdgeNormBoundaryInt(const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:113
virtual NekDouble v_VectorFlux(const Array< OneD, Array< OneD, NekDouble >> &vec)
Definition: Expansion.cpp:785
ExpansionWeakPtr m_elementLeft
Definition: Expansion.h:276
void AddFaceNormBoundaryInt(const int face, const std::shared_ptr< Expansion > &FaceExp, const Array< OneD, const NekDouble > &Fn, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:128
void DGDeriv(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, ExpansionSharedPtr > &EdgeExp, Array< OneD, Array< OneD, NekDouble >> &coeffs, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:135
void AddRobinMassMatrix(const int traceid, const Array< OneD, const NekDouble > &primCoeffs, DNekMatSharedPtr &inoutmat)
Definition: Expansion.h:242
virtual const Array< OneD, const NekDouble > & v_GetPhysNormals(void)
Definition: Expansion.cpp:867
SpatialDomains::GeomFactorsSharedPtr m_metricinfo
Definition: Expansion.h:273
virtual StdRegions::Orientation v_GetTraceOrient(int trace)
Definition: Expansion.cpp:806
void StdDerivBaseOnTraceMat(Array< OneD, DNekMatSharedPtr > &DerivMat)
Definition: Expansion.cpp:468
void ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
Definition: Expansion.cpp:105
ExpansionSharedPtr GetTraceExp(const int traceid)
Definition: Expansion.h:406
void SetPhysNormals(Array< OneD, const NekDouble > &normal)
Definition: Expansion.h:232
virtual void v_ReOrientTracePhysMap(const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1=-1)
Definition: Expansion.cpp:845
void ReOrientTracePhysMap(const StdRegions::Orientation orient, Array< OneD, int > &idmap, const int nq0, const int nq1)
Definition: Expansion.h:213
virtual void v_MultiplyByQuadratureMetric(const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:413
virtual void v_ComputeLaplacianMetric()
Definition: Expansion.h:301
virtual void v_GetCoords(Array< OneD, NekDouble > &coords_1, Array< OneD, NekDouble > &coords_2, Array< OneD, NekDouble > &coords_3)
Definition: Expansion.cpp:524
virtual void v_ExtractDataToCoeffs(const NekDouble *data, const std::vector< unsigned int > &nummodes, const int nmodes_offset, NekDouble *coeffs, std::vector< LibUtilities::BasisType > &fromType)
Definition: Expansion.cpp:741
ExpansionSharedPtr GetRightAdjacentElementExp() const
Definition: Expansion.h:439
StdRegions::Orientation GetTraceOrient(int trace)
Definition: Expansion.h:167
virtual void v_AddEdgeNormBoundaryInt(const int edge, const std::shared_ptr< Expansion > &EdgeExp, const Array< OneD, const NekDouble > &Fx, const Array< OneD, const NekDouble > &Fy, Array< OneD, NekDouble > &outarray)
Definition: Expansion.cpp:750
IndexMapValuesSharedPtr GetIndexMap(const IndexMapKey &ikey)
Definition: Expansion.h:145
DNekScalMatSharedPtr GetLocMatrix(const LocalRegions::MatrixKey &mkey)
Definition: Expansion.cpp:88
virtual void v_SetUpPhysNormals(const int id)
Definition: Expansion.cpp:879
virtual void v_GetTracePhysMap(const int edge, Array< OneD, int > &outarray)
Definition: Expansion.cpp:838
virtual void AddRobinTraceContribution(const int traceid, const Array< OneD, const NekDouble > &primCoeffs, const Array< OneD, NekDouble > &incoeffs, Array< OneD, NekDouble > &coeffs)
Definition: Expansion.h:254
int GetRightAdjacentElementTrace() const
Definition: Expansion.h:451
void TraceNormLen(const int traceid, NekDouble &h, NekDouble &p)
Definition: Expansion.h:249
const NormalVector & GetTraceNormal(const int id)
Definition: Expansion.cpp:250
virtual void v_NormalTraceDerivFactors(Array< OneD, Array< OneD, NekDouble >> &factors, Array< OneD, Array< OneD, NekDouble >> &d0factors, Array< OneD, Array< OneD, NekDouble >> &d1factors)
Definition: Expansion.cpp:795
virtual void v_GetTracePhysVals(const int trace, const StdRegions::StdExpansionSharedPtr &TraceExp, const Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray, StdRegions::Orientation orient)
Definition: Expansion.cpp:828
void AlignVectorToCollapsedDir(const int dir, const Array< OneD, const NekDouble > &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
Definition: Expansion.h:150
Expansion(SpatialDomains::GeometrySharedPtr pGeom)
Definition: Expansion.cpp:47
LibUtilities::NekManager< IndexMapKey, IndexMapValues, IndexMapKey::opLess > m_indexMapManager
Definition: Expansion.h:269
void ComputeTraceNormal(const int id)
Definition: Expansion.h:222
void SetCoeffsToOrientation(StdRegions::Orientation dir, Array< OneD, const NekDouble > &inarray, Array< OneD, NekDouble > &outarray)
Definition: Expansion.h:172
DNekMatSharedPtr BuildTransformationMatrix(const DNekScalMatSharedPtr &r_bnd, const StdRegions::MatrixType matrixType)
Definition: Expansion.cpp:94
const Array< OneD, const NekDouble > & GetPhysNormals(void)
Definition: Expansion.h:227
The base class for all shapes.
Definition: StdExpansion.h:71
int GetNtraces() const
Returns the number of trace elements connected to this element.
Definition: StdExpansion.h:359
std::shared_ptr< Expansion > ExpansionSharedPtr
Definition: Expansion.h:68
std::shared_ptr< IndexMapValues > IndexMapValuesSharedPtr
Definition: IndexMapKey.h:128
std::weak_ptr< Expansion > ExpansionWeakPtr
Definition: Expansion.h:69
Array< OneD, Array< OneD, NekDouble > > NormalVector
Definition: Expansion.h:53
std::vector< ExpansionSharedPtr > ExpansionVector
Definition: Expansion.h:70
std::map< MetricType, Array< OneD, NekDouble > > MetricMap
Definition: Expansion.h:71
std::shared_ptr< GeomFactors > GeomFactorsSharedPtr
Pointer to a GeomFactors object.
Definition: GeomFactors.h:62
std::shared_ptr< Geometry > GeometrySharedPtr
Definition: Geometry.h:53
std::shared_ptr< StdExpansion > StdExpansionSharedPtr
std::map< StdRegions::VarCoeffType, Array< OneD, NekDouble > > VarCoeffMap
Definition: StdRegions.hpp:240
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:282
static ConstFactorMap NullConstFactorMap
Definition: StdRegions.hpp:283
static VarCoeffMap NullVarCoeffMap
Definition: StdRegions.hpp:241
The above copyright notice and this permission notice shall be included.
Definition: CoupledSolver.h:1
std::shared_ptr< DNekScalMat > DNekScalMatSharedPtr
std::shared_ptr< DNekScalBlkMat > DNekScalBlkMatSharedPtr
Definition: NekTypeDefs.hpp:79
std::shared_ptr< DNekMat > DNekMatSharedPtr
Definition: NekTypeDefs.hpp:75
double NekDouble