Nektar++
Static Public Member Functions | Static Public Attributes | Protected Member Functions | Private Member Functions | Friends | List of all members
Nektar::WallBC Class Reference

Wall boundary conditions for compressible flow problems. More...

#include <WallBC.h>

Inheritance diagram for Nektar::WallBC:
[legend]

Static Public Member Functions

static CFSBndCondSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const int pSpaceDim, const int bcRegion, const int cnt)
 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_Apply (Array< OneD, Array< OneD, NekDouble > > &Fwd, Array< OneD, Array< OneD, NekDouble > > &physarray, const NekDouble &time)
 
- Protected Member Functions inherited from Nektar::CFSBndCond
 CFSBndCond (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const int pSpaceDim, const int bcRegion, const int cnt)
 Constructor. More...
 
virtual void v_ApplyBwdWeight ()
 

Private Member Functions

 WallBC (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const int pSpaceDim, const int bcRegion, const int cnt)
 
virtual ~WallBC (void)
 

Friends

class MemoryManager< WallBC >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::CFSBndCond
virtual ~CFSBndCond ()
 
void Apply (Array< OneD, Array< OneD, NekDouble > > &Fwd, Array< OneD, Array< OneD, NekDouble > > &physarray, const NekDouble &time=0)
 Apply the boundary condition. More...
 
void ApplyBwdWeight ()
 Apply the Weight of boundary condition. More...
 
- Protected Attributes inherited from Nektar::CFSBndCond
LibUtilities::SessionReaderSharedPtr m_session
 Session reader. More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 Array of fields. More...
 
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
 Trace normals. More...
 
int m_spacedim
 Space dimension. More...
 
VariableConverterSharedPtr m_varConv
 Auxiliary object to convert variables. More...
 
NekDouble m_diffusionAveWeight
 Weight for average calculation of diffusion term. More...
 
NekDouble m_gamma
 Parameters of the flow. More...
 
NekDouble m_rhoInf
 
NekDouble m_pInf
 
NekDouble m_pOut
 
Array< OneD, NekDoublem_velInf
 
int m_bcRegion
 Id of the boundary region. More...
 
int m_offset
 Offset. More...
 

Detailed Description

Wall boundary conditions for compressible flow problems.

Definition at line 47 of file WallBC.h.

Constructor & Destructor Documentation

◆ WallBC()

Nektar::WallBC::WallBC ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const Array< OneD, Array< OneD, NekDouble > > &  pTraceNormals,
const int  pSpaceDim,
const int  bcRegion,
const int  cnt 
)
private

Definition at line 49 of file WallBC.cpp.

55  : CFSBndCond(pSession, pFields, pTraceNormals, pSpaceDim, bcRegion, cnt)
56 {
58 }
CFSBndCond(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const int pSpaceDim, const int bcRegion, const int cnt)
Constructor.
Definition: CFSBndCond.cpp:47
NekDouble m_diffusionAveWeight
Weight for average calculation of diffusion term.
Definition: CFSBndCond.h:99

References Nektar::CFSBndCond::m_diffusionAveWeight.

◆ ~WallBC()

virtual Nektar::WallBC::~WallBC ( void  )
inlineprivatevirtual

Definition at line 84 of file WallBC.h.

84 {};

Member Function Documentation

◆ create()

static CFSBndCondSharedPtr Nektar::WallBC::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
const Array< OneD, Array< OneD, NekDouble > > &  pTraceNormals,
const int  pSpaceDim,
const int  bcRegion,
const int  cnt 
)
inlinestatic

Creates an instance of this class.

Definition at line 54 of file WallBC.h.

59  {
61  AllocateSharedPtr(pSession, pFields,
62  pTraceNormals, pSpaceDim, bcRegion, cnt);
63  return p;
64  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< CFSBndCond > CFSBndCondSharedPtr
A shared pointer to a boundary condition object.
Definition: CFSBndCond.h:48

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

◆ v_Apply()

void Nektar::WallBC::v_Apply ( Array< OneD, Array< OneD, NekDouble > > &  Fwd,
Array< OneD, Array< OneD, NekDouble > > &  physarray,
const NekDouble time 
)
protectedvirtual

Implements Nektar::CFSBndCond.

Definition at line 60 of file WallBC.cpp.

64 {
65  boost::ignore_unused(time);
66 
67  int i;
68  int nVariables = physarray.size();
69 
70  const Array<OneD, const int> &traceBndMap
71  = m_fields[0]->GetTraceBndMap();
72 
73  // Adjust the physical values of the trace to take
74  // user defined boundaries into account
75  int e, id1, id2, nBCEdgePts, eMax;
76 
77  eMax = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->GetExpSize();
78 
79  for (e = 0; e < eMax; ++e)
80  {
81  nBCEdgePts = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->
82  GetExp(e)->GetTotPoints();
83  id1 = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->
84  GetPhys_Offset(e);
85  id2 = m_fields[0]->GetTrace()->GetPhys_Offset(traceBndMap[m_offset+e]);
86 
87  // Boundary condition for epsilon term.
88  if (nVariables == m_spacedim+3)
89  {
90  Vmath::Zero(nBCEdgePts, &Fwd[nVariables-1][id2], 1);
91  }
92 
93  // For 2D/3D, define: v* = v - 2(v.n)n
94  Array<OneD, NekDouble> tmp(nBCEdgePts, 0.0);
95 
96  // Calculate (v.n)
97  for (i = 0; i < m_spacedim; ++i)
98  {
99  Vmath::Vvtvp(nBCEdgePts,
100  &Fwd[1+i][id2], 1,
101  &m_traceNormals[i][id2], 1,
102  &tmp[0], 1,
103  &tmp[0], 1);
104  }
105 
106  // Calculate 2.0(v.n)
107  Vmath::Smul(nBCEdgePts, -2.0, &tmp[0], 1, &tmp[0], 1);
108 
109  // Calculate v* = v - 2.0(v.n)n
110  for (i = 0; i < m_spacedim; ++i)
111  {
112  Vmath::Vvtvp(nBCEdgePts,
113  &tmp[0], 1,
114  &m_traceNormals[i][id2], 1,
115  &Fwd[1+i][id2], 1,
116  &Fwd[1+i][id2], 1);
117  }
118 
119  // Copy boundary adjusted values into the boundary expansion
120  for (i = 0; i < nVariables; ++i)
121  {
122  Vmath::Vcopy(nBCEdgePts, &Fwd[i][id2], 1,
123  &(m_fields[i]->GetBndCondExpansions()[m_bcRegion]->
124  UpdatePhys())[id1], 1);
125  }
126  }
127 }
int m_spacedim
Space dimension.
Definition: CFSBndCond.h:95
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Trace normals.
Definition: CFSBndCond.h:93
int m_bcRegion
Id of the boundary region.
Definition: CFSBndCond.h:109
int m_offset
Offset.
Definition: CFSBndCond.h:111
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
Definition: CFSBndCond.h:91
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:513
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.cpp:225
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:436
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1199

References Nektar::CFSBndCond::m_bcRegion, Nektar::CFSBndCond::m_fields, Nektar::CFSBndCond::m_offset, Nektar::CFSBndCond::m_spacedim, Nektar::CFSBndCond::m_traceNormals, Vmath::Smul(), Vmath::Vcopy(), Vmath::Vvtvp(), and Vmath::Zero().

Friends And Related Function Documentation

◆ MemoryManager< WallBC >

friend class MemoryManager< WallBC >
friend

Definition at line 1 of file WallBC.h.

Member Data Documentation

◆ className

std::string Nektar::WallBC::className
static
Initial value:
RegisterCreatorFunction("Wall",
"Slip wall boundary condition.")
static CFSBndCondSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const int pSpaceDim, const int bcRegion, const int cnt)
Creates an instance of this class.
Definition: WallBC.h:54
CFSBndCondFactory & GetCFSBndCondFactory()
Declaration of the boundary condition factory singleton.
Definition: CFSBndCond.cpp:41

Name of the class.

Definition at line 67 of file WallBC.h.