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

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

#include <WallViscousBC.h>

Inheritance diagram for Nektar::WallViscousBC:
[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 classNameViscous
 Name of the class. More...
 
static std::string classNameAdiabatic
 

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 ()
 

Protected Attributes

Array< OneD, Array< OneD, NekDouble > > m_bndPhys
 
- 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...
 

Private Member Functions

 WallViscousBC (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 ~WallViscousBC (void)
 

Friends

class MemoryManager< WallViscousBC >
 

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...
 

Detailed Description

Wall boundary conditions for viscous compressible flow problems.

Definition at line 47 of file WallViscousBC.h.

Constructor & Destructor Documentation

◆ WallViscousBC()

Nektar::WallViscousBC::WallViscousBC ( 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 54 of file WallViscousBC.cpp.

60  : CFSBndCond(pSession, pFields, pTraceNormals, pSpaceDim, bcRegion, cnt)
61 {
63 
64  m_bndPhys = Array<OneD, Array<OneD, NekDouble> > (m_fields.size());
65 }
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
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
Definition: CFSBndCond.h:91
Array< OneD, Array< OneD, NekDouble > > m_bndPhys
Definition: WallViscousBC.h:74

References m_bndPhys, Nektar::CFSBndCond::m_diffusionAveWeight, and Nektar::CFSBndCond::m_fields.

◆ ~WallViscousBC()

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

Definition at line 89 of file WallViscousBC.h.

89 {};

Member Function Documentation

◆ create()

static CFSBndCondSharedPtr Nektar::WallViscousBC::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 WallViscousBC.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::WallViscousBC::v_Apply ( Array< OneD, Array< OneD, NekDouble > > &  Fwd,
Array< OneD, Array< OneD, NekDouble > > &  physarray,
const NekDouble time 
)
protectedvirtual

Implements Nektar::CFSBndCond.

Definition at line 67 of file WallViscousBC.cpp.

71 {
72  boost::ignore_unused(time);
73 
74  int i;
75  int nVariables = physarray.size();
76 
77  // Find the fields whose WallViscous/Adiabatic-BC is time-dependent
78  // Update variables on the boundaries of these fields
79  // Get the updated variables on the WallViscous/Adiabatic boundary
80  //
81  // Maybe the EvaluateBoundaryConditions() should be put upstream to
82  // CompressibleFlowSystem::NumCalRiemFluxJac(), So that the BCs will not
83  // be repeatedly updated when there are more than one time-dependent BC.
84  std::string varName;
85  for (i = 0; i < nVariables; ++i)
86  {
87  if (m_fields[i]->GetBndConditions()[m_bcRegion]->IsTimeDependent())
88  {
89  varName = m_session->GetVariable(i);
90  m_fields[i]->EvaluateBoundaryConditions(time, varName);
91 
92  m_bndPhys[i] = m_fields[i]->GetBndCondExpansions()[m_bcRegion]
93  ->UpdatePhys();
94  }
95  }
96 
97 
98  const Array<OneD, const int> &traceBndMap
99  = m_fields[0]->GetTraceBndMap();
100 
101  // Take into account that for PDE based shock capturing, eps = 0 at the
102  // wall. Adjust the physical values of the trace to take user defined
103  // boundaries into account
104  int e, id1, id2, nBCEdgePts, eMax;
105 
106  eMax = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->GetExpSize();
107 
108  for (e = 0; e < eMax; ++e)
109  {
110  nBCEdgePts = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->
111  GetExp(e)->GetTotPoints();
112  id1 = m_fields[0]->GetBndCondExpansions()[m_bcRegion]->
113  GetPhys_Offset(e);
114  id2 = m_fields[0]->GetTrace()->GetPhys_Offset(traceBndMap[m_offset+e]);
115 
116  // Boundary condition for epsilon term.
117  if (nVariables == m_spacedim+3)
118  {
119  Vmath::Zero(nBCEdgePts, &Fwd[nVariables-1][id2], 1);
120  }
121 
122  // V = - Vin
123  for (i = 0; i < m_spacedim; i++)
124  {
125  Vmath::Neg(nBCEdgePts, &Fwd[i+1][id2], 1);
126  }
127 
128  // Superimpose the perturbation
129  for (i = 0; i < nVariables; ++i)
130  {
131  if (m_fields[i]->GetBndConditions()[m_bcRegion]->IsTimeDependent())
132  {
133  Vmath::Vadd(nBCEdgePts, &m_bndPhys[i][id1], 1,
134  &Fwd[i][id2], 1, &Fwd[i][id2], 1);
135  }
136  }
137 
138  // Copy boundary adjusted values into the boundary expansion
139  for (i = 0; i < nVariables; ++i)
140  {
141  Vmath::Vcopy(nBCEdgePts, &Fwd[i][id2], 1,
142  &(m_fields[i]->GetBndCondExpansions()[m_bcRegion]->
143  UpdatePhys())[id1], 1);
144  }
145  }
146 }
LibUtilities::SessionReaderSharedPtr m_session
Session reader.
Definition: CFSBndCond.h:89
int m_spacedim
Space dimension.
Definition: CFSBndCond.h:95
int m_bcRegion
Id of the boundary region.
Definition: CFSBndCond.h:109
int m_offset
Offset.
Definition: CFSBndCond.h:111
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:461
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:322
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, m_bndPhys, Nektar::CFSBndCond::m_fields, Nektar::CFSBndCond::m_offset, Nektar::CFSBndCond::m_session, Nektar::CFSBndCond::m_spacedim, Vmath::Neg(), Vmath::Vadd(), Vmath::Vcopy(), and Vmath::Zero().

Friends And Related Function Documentation

◆ MemoryManager< WallViscousBC >

friend class MemoryManager< WallViscousBC >
friend

Definition at line 1 of file WallViscousBC.h.

Member Data Documentation

◆ classNameAdiabatic

std::string Nektar::WallViscousBC::classNameAdiabatic
static
Initial value:
RegisterCreatorFunction("WallAdiabatic",
"Adiabatic 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: WallViscousBC.h:54
CFSBndCondFactory & GetCFSBndCondFactory()
Declaration of the boundary condition factory singleton.
Definition: CFSBndCond.cpp:41

Definition at line 68 of file WallViscousBC.h.

◆ classNameViscous

std::string Nektar::WallViscousBC::classNameViscous
static
Initial value:
RegisterCreatorFunction("WallViscous",
"No-slip (viscous) wall boundary condition.")

Name of the class.

Definition at line 67 of file WallViscousBC.h.

◆ m_bndPhys

Array<OneD, Array<OneD, NekDouble> > Nektar::WallViscousBC::m_bndPhys
protected

Definition at line 74 of file WallViscousBC.h.

Referenced by v_Apply(), and WallViscousBC().