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::NonSmoothShockCapture Class Reference

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

#include <NonSmoothShockCapture.h>

Inheritance diagram for Nektar::NonSmoothShockCapture:
Inheritance graph
[legend]
Collaboration diagram for Nektar::NonSmoothShockCapture:
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_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...
 
virtual void v_DoArtificialDiffusion (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 

Private Member Functions

 NonSmoothShockCapture (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
 
virtual ~NonSmoothShockCapture (void)
 

Friends

class MemoryManager< NonSmoothShockCapture >
 

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

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

Definition at line 49 of file NonSmoothShockCapture.h.

Constructor & Destructor Documentation

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

Definition at line 48 of file NonSmoothShockCapture.cpp.

52  : ArtificialDiffusion(pSession, pFields, spacedim)
53 {
54 }
ArtificialDiffusion(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int spacedim)
Constructor.
virtual Nektar::NonSmoothShockCapture::~NonSmoothShockCapture ( void  )
inlineprivatevirtual

Definition at line 81 of file NonSmoothShockCapture.h.

81 {};

Member Function Documentation

static ArtificialDiffusionSharedPtr Nektar::NonSmoothShockCapture::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 NonSmoothShockCapture.h.

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

60  {
63  AllocateSharedPtr(pSession, pFields, spacedim);
64  return p;
65  }
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::NonSmoothShockCapture::v_GetArtificialViscosity ( const Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, NekDouble > &  mu 
)
protectedvirtual

Implements Nektar::ArtificialDiffusion.

Definition at line 56 of file NonSmoothShockCapture.cpp.

References Nektar::ArtificialDiffusion::m_fields, Nektar::ArtificialDiffusion::m_Kappa, Nektar::ArtificialDiffusion::m_mu0, Nektar::ArtificialDiffusion::m_Skappa, and Nektar::ArtificialDiffusion::m_varConv.

59 {
60  const int nElements = m_fields[0]->GetExpSize();
61  int PointCount = 0;
62  int nTotQuadPoints = m_fields[0]->GetTotPoints();
63 
64  Array<OneD, NekDouble> Sensor (nTotQuadPoints, 0.0);
65  Array<OneD, NekDouble> SensorKappa(nTotQuadPoints, 0.0);
66 
67  m_varConv->GetSensor(m_fields[0], physfield, Sensor, SensorKappa);
68 
69  for (int e = 0; e < nElements; e++)
70  {
71  // Threshold value specified in C. Biottos thesis. Based on a 1D
72  // shock tube problem S_k = log10(1/p^4). See G.E. Barter and
73  // D.L. Darmofal. Shock Capturing with PDE-based artificial
74  // diffusion for DGFEM: Part 1 Formulation, Journal of Computational
75  // Physics 229 (2010) 1810-1827 for further reference
76 
77  // Adjustable depending on the coarsness of the mesh. Might want to
78  // move this variable into the session file
79 
80  int nQuadPointsElement = m_fields[0]->GetExp(e)->GetTotPoints();
81 
82  for (int n = 0; n < nQuadPointsElement; n++)
83  {
84  NekDouble mu_0 = m_mu0;
85 
86  if (Sensor[n+PointCount] < (m_Skappa-m_Kappa))
87  {
88  mu[n+PointCount] = 0;
89  }
90  else if (Sensor[n+PointCount] >= (m_Skappa-m_Kappa)
91  && Sensor[n+PointCount] <= (m_Skappa+m_Kappa))
92  {
93  mu[n+PointCount] = mu_0 * (0.5 * (1 + sin(
94  M_PI * (Sensor[n+PointCount] -
95  m_Skappa - m_Kappa) /
96  (2*m_Kappa))));
97  }
98  else if (Sensor[n+PointCount] > (m_Skappa+m_Kappa))
99  {
100  mu[n+PointCount] = mu_0;
101  }
102  }
103 
104  PointCount += nQuadPointsElement;
105  }
106 }
double NekDouble
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< NonSmoothShockCapture >
friend

Definition at line 53 of file NonSmoothShockCapture.h.

Member Data Documentation

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

Name of the class.

Definition at line 68 of file NonSmoothShockCapture.h.