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

Outflow characteristic boundary conditions for compressible flow problems. More...

#include <EnforceEntropyPressure.h>

Inheritance diagram for Nektar::EnforceEntropyPressure:
[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 Array< OneD, Array< OneD, NekDouble > > &pGridVelocity, 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

void v_Apply (Array< OneD, Array< OneD, NekDouble > > &Fwd, Array< OneD, Array< OneD, NekDouble > > &physarray, const NekDouble &time) override
 
void GenerateRotationMatrices (const Array< OneD, const Array< OneD, NekDouble > > &normals)
 
void FromToRotation (Array< OneD, const NekDouble > &from, Array< OneD, const NekDouble > &to, NekDouble *mat)
 
SOLVER_UTILS_EXPORT void rotateToNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
SOLVER_UTILS_EXPORT void rotateFromNormal (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const Array< OneD, const Array< OneD, NekDouble > > &normals, const Array< OneD, const Array< OneD, NekDouble > > &vecLocs, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
- 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 Array< OneD, Array< OneD, NekDouble > > &pGridVelocity, const int pSpaceDim, const int bcRegion, const int cnt)
 Constructor. More...
 
virtual ~CFSBndCond ()=default
 
virtual void v_Apply (Array< OneD, Array< OneD, NekDouble > > &Fwd, Array< OneD, Array< OneD, NekDouble > > &physarray, const NekDouble &time)=0
 
virtual void v_ApplyBwdWeight ()
 

Protected Attributes

int m_npts
 
Array< OneD, NekDoublem_rhoBC
 
Array< OneD, Array< OneD, NekDouble > > m_velBC
 
Array< OneD, NekDoublem_pBC
 
Array< OneD, NekDoublem_VnInf
 Reference normal velocity. More...
 
Array< OneD, int > m_bndToTraceMap
 
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...
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
 Grid Velocity. 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
 
NekDouble m_angVel
 
int m_bcRegion
 Id of the boundary region. More...
 
int m_offset
 Offset. More...
 

Private Member Functions

 EnforceEntropyPressure (const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const Array< OneD, Array< OneD, NekDouble > > &pGridVelocity, const int pSpaceDim, const int bcRegion, const int cnt)
 
 ~EnforceEntropyPressure (void) override=default
 

Friends

class MemoryManager< EnforceEntropyPressure >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::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

Outflow characteristic boundary conditions for compressible flow problems.

Definition at line 51 of file EnforceEntropyPressure.h.

Constructor & Destructor Documentation

◆ EnforceEntropyPressure()

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

Definition at line 49 of file EnforceEntropyPressure.cpp.

55 : CFSBndCond(pSession, pFields, pTraceNormals, pGridVelocity, pSpaceDim,
56 bcRegion, cnt)
57{
58
60 m_fields[0]->GetBndCondExpansions()[m_bcRegion];
61
62 //-> Gather a list of index from trace to this boundary
63 m_npts = bndexp->GetTotPoints();
64
65 m_bndToTraceMap = Array<OneD, int>(m_npts, -1);
66
67 const Array<OneD, const int> &traceBndMap = m_fields[0]->GetTraceBndMap();
68
69 // Construct a map for the boundary to trace map for easy acess to
70 // phys space points
71 int cnt1 = 0;
72 for (int e = 0; e < bndexp->GetNumElmts(); ++e)
73 {
74 int nTracePts = bndexp->GetExp(e)->GetTotPoints();
75
76 int id =
77 m_fields[0]->GetTrace()->GetPhys_Offset(traceBndMap[m_offset + e]);
78
79 // Loop on the points of the m_bcRegion
80 for (int i = 0; i < nTracePts; i++)
81 {
82 // the ith point in region e
83 m_bndToTraceMap[cnt1++] = id + i;
84 }
85 }
86
87 Array<OneD, Array<OneD, NekDouble>> BCvals(m_fields.size());
88 m_bndPhys = Array<OneD, Array<OneD, NekDouble>>(m_fields.size());
89
90 for (int i = 0; i < m_fields.size(); ++i)
91 {
92 m_bndPhys[i] =
93 m_fields[i]->GetBndCondExpansions()[m_bcRegion]->UpdatePhys();
94
95 BCvals[i] = Array<OneD, NekDouble>(m_npts);
96 Vmath::Vcopy(m_npts, m_bndPhys[i], 1, BCvals[i], 1);
97 }
98
99 // Set up boudnary required BCs
100 m_rhoBC = Array<OneD, NekDouble>(m_npts);
101 Vmath::Vcopy(m_npts, BCvals[0], 1, m_rhoBC, 1);
102
103 m_velBC = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
104 // Evaluate velocity on boundary
105 for (int i = 0; i < m_spacedim; ++i)
106 {
107 m_velBC[i] = Array<OneD, NekDouble>(m_npts);
108 Vmath::Vcopy(m_npts, BCvals[i + 1], 1, m_velBC[i], 1);
109 Vmath::Vdiv(m_npts, m_velBC[i], 1, m_rhoBC, 1, m_velBC[i], 1);
110 }
111 m_pBC = Array<OneD, NekDouble>(m_npts);
112 m_varConv->GetPressure(BCvals, m_pBC);
113
114 // Computing the normal velocity for characteristics coming
115 // from outside the computational domain
116 m_VnInf = Array<OneD, NekDouble>(m_npts, 0.0);
117 for (int i = 0; i < m_spacedim; i++)
118 {
119 for (int j = 0; j < m_npts; ++j)
120 {
121 m_VnInf[j] += m_traceNormals[i][m_bndToTraceMap[j]] * m_velBC[i][j];
122 }
123 }
124}
int m_spacedim
Space dimension.
Definition: CFSBndCond.h:98
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Trace normals.
Definition: CFSBndCond.h:94
int m_bcRegion
Id of the boundary region.
Definition: CFSBndCond.h:113
VariableConverterSharedPtr m_varConv
Auxiliary object to convert variables.
Definition: CFSBndCond.h:100
CFSBndCond(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const Array< OneD, Array< OneD, NekDouble > > &pGridVelocity, const int pSpaceDim, const int bcRegion, const int cnt)
Constructor.
Definition: CFSBndCond.cpp:46
int m_offset
Offset.
Definition: CFSBndCond.h:115
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array of fields.
Definition: CFSBndCond.h:92
Array< OneD, Array< OneD, NekDouble > > m_velBC
Array< OneD, NekDouble > m_VnInf
Reference normal velocity.
Array< OneD, NekDouble > m_rhoBC
Array< OneD, Array< OneD, NekDouble > > m_bndPhys
std::shared_ptr< ExpList > ExpListSharedPtr
Shared pointer to an ExpList object.
void Vdiv(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x/y.
Definition: Vmath.hpp:126
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References Nektar::CFSBndCond::m_bcRegion, m_bndPhys, m_bndToTraceMap, Nektar::CFSBndCond::m_fields, m_npts, Nektar::CFSBndCond::m_offset, m_pBC, m_rhoBC, Nektar::CFSBndCond::m_spacedim, Nektar::CFSBndCond::m_traceNormals, Nektar::CFSBndCond::m_varConv, m_velBC, m_VnInf, Vmath::Vcopy(), and Vmath::Vdiv().

◆ ~EnforceEntropyPressure()

Nektar::EnforceEntropyPressure::~EnforceEntropyPressure ( void  )
overrideprivatedefault

Member Function Documentation

◆ create()

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

Creates an instance of this class.

Definition at line 57 of file EnforceEntropyPressure.h.

63 {
66 pSession, pFields, pTraceNormals, pGridVelocity, pSpaceDim,
67 bcRegion, cnt);
68 return p;
69 }
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:52

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

◆ FromToRotation()

void Nektar::EnforceEntropyPressure::FromToRotation ( Array< OneD, const NekDouble > &  from,
Array< OneD, const NekDouble > &  to,
NekDouble mat 
)
protected

◆ GenerateRotationMatrices()

void Nektar::EnforceEntropyPressure::GenerateRotationMatrices ( const Array< OneD, const Array< OneD, NekDouble > > &  normals)
protected

◆ rotateFromNormal()

SOLVER_UTILS_EXPORT void Nektar::EnforceEntropyPressure::rotateFromNormal ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
const Array< OneD, const Array< OneD, NekDouble > > &  normals,
const Array< OneD, const Array< OneD, NekDouble > > &  vecLocs,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protected

◆ rotateToNormal()

SOLVER_UTILS_EXPORT void Nektar::EnforceEntropyPressure::rotateToNormal ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
const Array< OneD, const Array< OneD, NekDouble > > &  normals,
const Array< OneD, const Array< OneD, NekDouble > > &  vecLocs,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protected

◆ v_Apply()

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

Implements Nektar::CFSBndCond.

Definition at line 126 of file EnforceEntropyPressure.cpp.

130{
131 int i, j;
132 int nDimensions = m_spacedim;
133
134 Array<OneD, Array<OneD, NekDouble>> FwdBnd(Fwd.size());
135 Array<OneD, Array<OneD, NekDouble>> bndPhys(Fwd.size());
136
137 // make a local copy of Fwd along boundary of interest
138 for (i = 0; i < Fwd.size(); ++i)
139 {
140 FwdBnd[i] = Array<OneD, NekDouble>(m_npts);
141 for (j = 0; j < m_npts; ++j)
142 {
143 FwdBnd[i][j] = Fwd[i][m_bndToTraceMap[j]];
144 }
145 }
146
147 // Computing the normal velocity for characteristics coming
148 // from inside the computational domain
149 Array<OneD, NekDouble> Vn(m_npts, 0.0);
150
151 for (i = 0; i < nDimensions; ++i)
152 {
153 for (j = 0; j < m_npts; ++j)
154 {
155 Vn[j] += m_traceNormals[i][m_bndToTraceMap[j]] * FwdBnd[i + 1][j];
156 }
157 }
158 // divide by density.
159 Vmath::Vdiv(m_npts, Vn, 1, FwdBnd[0], 1, Vn, 1);
160
161 // Get speed of sound
162 Array<OneD, NekDouble> pressure(m_npts);
163 Array<OneD, NekDouble> soundSpeed(m_npts);
164
165 m_varConv->GetPressure(FwdBnd, pressure);
166 m_varConv->GetSoundSpeed(FwdBnd, soundSpeed);
167
168 // Get Mach. Note: it is computed by Vn/c
169 Array<OneD, NekDouble> Mach(m_npts, 0.0);
170 Vmath::Vdiv(m_npts, Vn, 1, soundSpeed, 1, Mach, 1);
171 Vmath::Vabs(m_npts, Mach, 1, Mach, 1);
172
173 // Auxiliary variables
174 Array<OneD, NekDouble> velBC(nDimensions, 0.0);
175
176 // L represents properties outside boundary
177 // R represents properties inside boundary (numerical state)
178 NekDouble rhoL, uL, pL;
179 NekDouble EBC, rR, cstar, pstar, rhostar, ustar; // vn
180
181 NekDouble gamMinOne = m_gamma - 1.0;
182 NekDouble twoOverGamMinOne = 2.0 / gamMinOne;
183 NekDouble gamInv = 1.0 / m_gamma;
184
185 // Loop on m_bcRegions
186 for (int pnt = 0; pnt < m_npts; ++pnt)
187 {
188 // Impose inflow Riemann invariant
189 if (Vn[pnt] <= 0.0)
190 {
191 // Subsonic flows
192 if (Mach[pnt] < 1.00)
193 {
194 // right characteristic
195 rR = -Vn[pnt] - sqrt(m_gamma * pressure[pnt] / FwdBnd[0][pnt]) *
196 twoOverGamMinOne;
197 // vn = -m_VnInf[pnt]; //vn BC
198
199 // fix rhostar and pstar to be the input values
200 // compute ustar using left-pointing characteristic line IR^-
201 pstar = m_pBC[pnt];
202 rhostar = m_rhoBC[pnt];
203 cstar = sqrt(m_gamma * pstar / rhostar);
204 ustar = rR + cstar * twoOverGamMinOne;
205
206 // add supplement equation that rhoL=rhostar
207 // then pL=pstar, according to IL^0
208 // and uL=ustar, according to IL^+
209 rhoL = rhostar;
210 pL = pstar;
211 uL = ustar;
212
213 // std subsnoic inflow
214 // rhoL = m_rhoBC[pnt];
215 // uL = -m_VnInf[pnt];
216 // pL = m_pBC[pnt];
217 }
218 else // Supersonic inflow
219 {
220 // all characteristics are from left so just impose
221 // star state to left values
222 // Note: m_vnInf is the negative of the normal velocity
223 // across boundary
224 rhoL = m_rhoBC[pnt];
225 uL = -m_VnInf[pnt];
226 pL = m_pBC[pnt];
227 }
228
229 // Boundary energy
230 EBC = pL * twoOverGamMinOne * 0.5;
231
232 // evaluate the different between the left state normal
233 // velocity and that from the desired condition (note
234 // m_VnInf is using an outwards normal definition.
235 NekDouble VnDiff = uL + m_VnInf[pnt];
236
237 // Boundary velocities & Kinite energy
238 // Note: normals are negated since they point outwards in
239 // the domain
240
241 // Note: Can just use the BC values directly!!
242 for (j = 0; j < nDimensions; ++j)
243 {
244 // Set velocity to the desired conditions modified to
245 // take account of the normal state for Riemann
246 // problem. (Negative accounts for outwards normal definition)
247 // velBC[j] = m_velBC[j][pnt];
248 velBC[j] = m_velBC[j][pnt] -
249 VnDiff * m_traceNormals[j][m_bndToTraceMap[pnt]];
250
251 EBC += 0.5 * rhoL * velBC[j] * velBC[j];
252 }
253
254 //-------------------------------------------------------------------------
255 // Impose Left hand Riemann Invariant boundary conditions
256 m_bndPhys[0][pnt] = rhoL;
257 for (j = 0; j < nDimensions; ++j)
258 {
259 m_bndPhys[j + 1][pnt] = rhoL * velBC[j];
260 }
261 m_bndPhys[nDimensions + 1][pnt] = EBC;
262 }
263 else // Outflow
264 {
265
266 // Note: Allowing the switch can cause worse convergence in this
267 // type BC.
268 // So improve it later.
269 if (Mach[pnt] < 1.00)
270 {
271 // subsonic outflow: fix pstar
272 rR = -Vn[pnt] - sqrt(m_gamma * pressure[pnt] / FwdBnd[0][pnt]) *
273 twoOverGamMinOne;
274 // vn = -m_VnInf[pnt];
275
276 pstar = m_pBC[pnt];
277 rhostar = FwdBnd[0][pnt] * pow((pstar / pressure[pnt]), gamInv);
278 cstar = sqrt(m_gamma * pstar / rhostar);
279 ustar = rR + cstar * twoOverGamMinOne;
280
281 rhoL = rhostar;
282 uL = ustar;
283 pL = pstar;
284 }
285 else
286 {
287 // supersonic outflow
288 // Just set to imposed state and let Riemann BC dictate values
289 rhoL = m_rhoBC[pnt];
290 uL = -m_VnInf[pnt];
291 pL = m_pBC[pnt];
292 }
293
294 // Boundary energy
295 EBC = pL * twoOverGamMinOne * 0.5;
296
297 // Boundary velocities & Kinite energy
298 // Note: normals are negated since they point outwards in
299 // the domain
300 for (j = 0; j < nDimensions; ++j)
301 {
302 velBC[j] = -1.0 * uL * m_traceNormals[j][m_bndToTraceMap[pnt]];
303 EBC += 0.5 * rhoL * velBC[j] * velBC[j];
304 }
305
306 // Impose Left hand Riemann Invariant boundary conditions
307 m_bndPhys[0][pnt] = rhoL;
308 for (j = 0; j < nDimensions; ++j)
309 {
310 m_bndPhys[j + 1][pnt] = rhoL * velBC[j];
311 }
312 m_bndPhys[nDimensions + 1][pnt] = EBC;
313 }
314 }
315}
NekDouble m_gamma
Parameters of the flow.
Definition: CFSBndCond.h:105
double NekDouble
void Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
Definition: Vmath.hpp:352
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References m_bndPhys, m_bndToTraceMap, Nektar::CFSBndCond::m_gamma, m_npts, m_pBC, m_rhoBC, Nektar::CFSBndCond::m_spacedim, Nektar::CFSBndCond::m_traceNormals, Nektar::CFSBndCond::m_varConv, m_velBC, m_VnInf, CG_Iterations::pressure, tinysimd::sqrt(), Vmath::Vabs(), and Vmath::Vdiv().

Friends And Related Function Documentation

◆ MemoryManager< EnforceEntropyPressure >

friend class MemoryManager< EnforceEntropyPressure >
friend

Definition at line 1 of file EnforceEntropyPressure.h.

Member Data Documentation

◆ className

std::string Nektar::EnforceEntropyPressure::className
static
Initial value:
=
"EnforceEntropyPressure", EnforceEntropyPressure::create,
"Riemann boundary condition enforcing entropy and pressure.")
static CFSBndCondSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const Array< OneD, Array< OneD, NekDouble > > &pTraceNormals, const Array< OneD, Array< OneD, NekDouble > > &pGridVelocity, const int pSpaceDim, const int bcRegion, const int cnt)
Creates an instance of this class.
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
CFSBndCondFactory & GetCFSBndCondFactory()
Declaration of the boundary condition factory singleton.
Definition: CFSBndCond.cpp:40

Name of the class.

Definition at line 72 of file EnforceEntropyPressure.h.

◆ m_bndPhys

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

Definition at line 89 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_bndToTraceMap

Array<OneD, int> Nektar::EnforceEntropyPressure::m_bndToTraceMap
protected

Definition at line 85 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_npts

int Nektar::EnforceEntropyPressure::m_npts
protected

Definition at line 75 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_pBC

Array<OneD, NekDouble> Nektar::EnforceEntropyPressure::m_pBC
protected

Definition at line 81 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_rhoBC

Array<OneD, NekDouble> Nektar::EnforceEntropyPressure::m_rhoBC
protected

Definition at line 77 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_velBC

Array<OneD, Array<OneD, NekDouble> > Nektar::EnforceEntropyPressure::m_velBC
protected

Definition at line 79 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().

◆ m_VnInf

Array<OneD, NekDouble> Nektar::EnforceEntropyPressure::m_VnInf
protected

Reference normal velocity.

Definition at line 83 of file EnforceEntropyPressure.h.

Referenced by EnforceEntropyPressure(), and v_Apply().