Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Private Member Functions | Friends | List of all members
Nektar::SmoothShockCapture Class Reference

Smooth artificial diffusion for shock capture for compressible flow problems. More...

#include <SmoothShockCapture.h>

Inheritance diagram for Nektar::SmoothShockCapture:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SmoothShockCapture:
Collaboration graph
[legend]

Static Public Member Functions

static ArtificialDiffusionSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
 Creates an instance of this class. More...
 

Static Public Attributes

static std::string className
 Name of the class. More...
 

Protected Member Functions

virtual void v_DoArtificialDiffusion (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
virtual void v_GetArtificialViscosity (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &mu)
 
- Protected Member Functions inherited from Nektar::ArtificialDiffusion
 ArtificialDiffusion (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
 Constructor. More...
 

Private Member Functions

 SmoothShockCapture (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
 
virtual ~SmoothShockCapture (void)
 
void GetForcingTerm (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > outarrayForcing)
 

Friends

class MemoryManager< SmoothShockCapture >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::ArtificialDiffusion
virtual ~ArtificialDiffusion ()
 
void DoArtificialDiffusion (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 Apply the artificial diffusion. More...
 
void GetArtificialViscosity (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &mu)
 Calculate the artificial viscosity. More...
 
- Protected Attributes inherited from Nektar::ArtificialDiffusion
LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_fields
 Array of fields. More...
 
VariableConverterSharedPtr m_varConv
 Auxiliary object to convert variables. More...
 
SolverUtils::DiffusionSharedPtr m_diffusion
 LDG Diffusion operator. More...
 
NekDouble m_FacL
 Parameters. More...
 
NekDouble m_FacH
 
NekDouble m_hFactor
 
NekDouble m_C1
 
NekDouble m_C2
 
NekDouble m_mu0
 
NekDouble m_Skappa
 
NekDouble m_Kappa
 

Detailed Description

Smooth artificial diffusion for shock capture for compressible flow problems.

Definition at line 49 of file SmoothShockCapture.h.

Constructor & Destructor Documentation

Nektar::SmoothShockCapture::SmoothShockCapture ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const int  spacedim 
)
private

Definition at line 48 of file SmoothShockCapture.cpp.

References ASSERTL0, and Nektar::ArtificialDiffusion::m_fields.

51  : ArtificialDiffusion(pSession, pFields, spacedim)
52 {
53  ASSERTL0(m_fields.num_elements() == spacedim + 3,
54  "Not enough variables for smooth shock capturing; "
55  "make sure you have added eps to variable list.");
56 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
ArtificialDiffusion(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
Constructor.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
virtual Nektar::SmoothShockCapture::~SmoothShockCapture ( void  )
inlineprivatevirtual

Definition at line 84 of file SmoothShockCapture.h.

84 {};

Member Function Documentation

static ArtificialDiffusionSharedPtr Nektar::SmoothShockCapture::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const int  spacedim 
)
inlinestatic

Creates an instance of this class.

Definition at line 56 of file SmoothShockCapture.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr(), and CellMLToNektar.cellml_metadata::p.

60  {
62  AllocateSharedPtr(pSession, pFields, spacedim);
63  return p;
64  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< ArtificialDiffusion > ArtificialDiffusionSharedPtr
A shared pointer to a artificial diffusion object.
void Nektar::SmoothShockCapture::GetForcingTerm ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > >  outarrayForcing 
)
private

Definition at line 162 of file SmoothShockCapture.cpp.

References Nektar::ArtificialDiffusion::m_C1, Nektar::ArtificialDiffusion::m_fields, Nektar::ArtificialDiffusion::m_hFactor, Nektar::ArtificialDiffusion::m_varConv, Vmath::Vadd(), and Vmath::Vmax().

Referenced by v_DoArtificialDiffusion().

165 {
166  const int nPts = m_fields[0]->GetTotPoints();
167  const int nvariables = m_fields.num_elements();
168  const int nElements = m_fields[0]->GetExpSize();
169 
170  Array<OneD, NekDouble> Sensor(nPts, 0.0);
171  Array<OneD, NekDouble> SensorKappa(nPts, 0.0);
172  Array <OneD, NekDouble > Lambda(nPts, 0.0);
173  Array <OneD, NekDouble > Tau(nPts, 1.0);
174  Array <OneD, NekDouble > soundspeed(nPts, 0.0);
175  Array <OneD, NekDouble > pressure(nPts, 0.0);
176  Array <OneD, NekDouble > absVelocity(nPts, 0.0);
177 
178  Array<OneD,int> pOrderElmt = m_fields[0]->EvalBasisNumModesMaxPerExp();
179  Array<OneD, NekDouble> pOrder (nPts, 0.0);
180 
181  // Thermodynamic related quantities
182  m_varConv->GetPressure(inarray, pressure);
183  m_varConv->GetSoundSpeed(inarray, pressure, soundspeed);
184  m_varConv->GetAbsoluteVelocity(inarray, absVelocity);
185  m_varConv->GetSensor(m_fields[0], inarray, Sensor, SensorKappa);
186 
187  // Determine the maximum wavespeed
188  Vmath::Vadd(nPts, absVelocity, 1, soundspeed, 1, Lambda, 1);
189 
190  NekDouble LambdaMax = Vmath::Vmax(nPts, Lambda, 1);
191 
192  int PointCount = 0;
193 
194  for (int e = 0; e < nElements; e++)
195  {
196  int nQuadPointsElement = m_fields[0]->GetExp(e)->GetTotPoints();
197 
198  for (int n = 0; n < nQuadPointsElement; n++)
199  {
200  pOrder[n + PointCount] = pOrderElmt[e];
201 
202  // order 1.0e-06
203  Tau[n + PointCount] =
204  1.0 / (m_C1*pOrder[n + PointCount]*LambdaMax);
205 
206  outarrayForcing[nvariables-1][n + PointCount] =
207  1 / Tau[n + PointCount] * (m_hFactor * LambdaMax /
208  pOrder[n + PointCount] *
209  SensorKappa[n + PointCount] -
210  inarray[nvariables-1][n + PointCount]);
211  }
212  PointCount += nQuadPointsElement;
213  }
214 }
Definition: Tau.hpp:39
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:779
double NekDouble
VariableConverterSharedPtr m_varConv
Auxiliary object to convert variables.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:299
void Nektar::SmoothShockCapture::v_DoArtificialDiffusion ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protectedvirtual

Reimplemented from Nektar::ArtificialDiffusion.

Definition at line 58 of file SmoothShockCapture.cpp.

References GetForcingTerm(), Nektar::ArtificialDiffusion::m_C2, Nektar::ArtificialDiffusion::m_diffusion, Nektar::ArtificialDiffusion::m_fields, Nektar::ArtificialDiffusion::m_varConv, Vmath::Smul(), Vmath::Vadd(), and Vmath::Vmax().

61 {
62  int i;
63  int nvariables = inarray.num_elements();
64  int npoints = m_fields[0]->GetNpoints();
65 
66  Array<OneD, Array<OneD, NekDouble> > outarrayDiff(nvariables);
67 
68  for (i = 0; i < nvariables; ++i)
69  {
70  outarrayDiff[i] = Array<OneD, NekDouble>(npoints, 0.0);
71  }
72 
73  m_diffusion->Diffuse(nvariables, m_fields, inarray, outarrayDiff);
74 
75  const Array<OneD, int> ExpOrder = m_fields[0]->EvalBasisNumModesMaxPerExp();
76 
77  NekDouble pOrder = Vmath::Vmax(ExpOrder.num_elements(), ExpOrder, 1);
78 
79  Array <OneD, NekDouble > a_vel (npoints, 0.0);
80  Array <OneD, NekDouble > u_abs (npoints, 0.0);
81  Array <OneD, NekDouble > tmp (npoints, 0.0);
82 
83  m_varConv->GetPressure(inarray, tmp);
84  m_varConv->GetSoundSpeed(inarray, tmp, a_vel);
85  m_varConv->GetAbsoluteVelocity(inarray, u_abs);
86 
87  Vmath::Vadd(npoints, a_vel, 1, u_abs, 1, tmp, 1);
88 
89  NekDouble max_wave_sp = Vmath::Vmax(npoints, tmp, 1);
90 
91  NekDouble fac = m_C2*max_wave_sp*pOrder;
92 
93  Vmath::Smul(npoints,
94  fac,
95  outarrayDiff[nvariables-1], 1,
96  outarrayDiff[nvariables-1], 1);
97 
98  for (i = 0; i < nvariables; ++i)
99  {
100  Vmath::Vadd(npoints,
101  outarray[i], 1,
102  outarrayDiff[i], 1,
103  outarray[i], 1);
104  }
105 
106  Array<OneD, Array<OneD, NekDouble> > outarrayForcing(nvariables);
107 
108  for (i = 0; i < nvariables; ++i)
109  {
110  outarrayForcing[i] = Array<OneD, NekDouble>(npoints, 0.0);
111  }
112 
113  GetForcingTerm(inarray, outarrayForcing);
114 
115  for (i = 0; i < nvariables; ++i)
116  {
117  // Add Forcing Term
118  Vmath::Vadd(npoints,
119  outarray[i], 1,
120  outarrayForcing[i], 1,
121  outarray[i], 1);
122  }
123 }
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.cpp:779
SolverUtils::DiffusionSharedPtr m_diffusion
LDG Diffusion operator.
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:213
double NekDouble
void GetForcingTerm(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > outarrayForcing)
VariableConverterSharedPtr m_varConv
Auxiliary object to convert variables.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:299
void Nektar::SmoothShockCapture::v_GetArtificialViscosity ( const Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, NekDouble > &  mu 
)
protectedvirtual

Implements Nektar::ArtificialDiffusion.

Definition at line 125 of file SmoothShockCapture.cpp.

References Nektar::ArtificialDiffusion::m_FacH, Nektar::ArtificialDiffusion::m_FacL, Nektar::ArtificialDiffusion::m_fields, Nektar::ArtificialDiffusion::m_mu0, and Nektar::ArtificialDiffusion::m_varConv.

128 {
129  int nvariables = physfield.num_elements();
130  int nPts = m_fields[0]->GetTotPoints();
131 
132  Array<OneD, NekDouble > sensor (nPts, 0.0);
133  Array<OneD, NekDouble > SensorKappa(nPts, 0.0);
134 
135  // Calculate sensor
136  m_varConv->GetSensor(m_fields[0], physfield, sensor, SensorKappa);
137 
138  NekDouble ThetaH = m_FacH;
139  NekDouble ThetaL = m_FacL;
140 
141  NekDouble Phi0 = (ThetaH+ThetaL)/2;
142  NekDouble DeltaPhi = ThetaH-Phi0;
143 
144  for (int e = 0; e < mu.num_elements(); e++)
145  {
146  if (physfield[nvariables-1][e] <= (Phi0 - DeltaPhi))
147  {
148  mu[e] = 0;
149  }
150  else if(physfield[nvariables-1][e] >= (Phi0 + DeltaPhi))
151  {
152  mu[e] = m_mu0;
153  }
154  else if(abs(physfield[nvariables-1][e]-Phi0) < DeltaPhi)
155  {
156  mu[e] = m_mu0/2*(1+sin(M_PI*
157  (physfield[nvariables-1][e]-Phi0)/(2*DeltaPhi)));
158  }
159  }
160 }
double NekDouble
NekDouble m_FacL
Parameters.
VariableConverterSharedPtr m_varConv
Auxiliary object to convert variables.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.

Friends And Related Function Documentation

friend class MemoryManager< SmoothShockCapture >
friend

Definition at line 53 of file SmoothShockCapture.h.

Member Data Documentation

std::string Nektar::SmoothShockCapture::className
static
Initial value:
RegisterCreatorFunction("Smooth",
"Smooth artificial diffusion for shock capture.")

Name of the class.

Definition at line 67 of file SmoothShockCapture.h.