Commit 8de2f756 authored by Dave Moxey's avatar Dave Moxey

Tidy up, fix another couple of missing BCs in Euler solver

(cherry picked from commit dce42b36)
parent 11a435c1
......@@ -250,7 +250,7 @@ namespace Nektar
fields[0]->GetTrace()->Upwind(Vn, Fwd[i], Bwd[i], numflux[i]);
}
// Extract internal values of the scalar variables for Neumann bcs
// Extract internal values of the scalar variables for Neumann bcs
Array< OneD, Array<OneD, NekDouble > > uplus(nScalars);
for (i = 0; i < nScalars; ++i)
......@@ -258,7 +258,7 @@ namespace Nektar
uplus[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
fields[i]->ExtractTracePhys(inarray[i], uplus[i]);
}
// Modify the values in case of boundary interfaces
if (fields[0]->GetBndCondExpansions().num_elements())
{
......@@ -300,17 +300,13 @@ namespace Nektar
Array<OneD, NekDouble> Tw(nTracePts, m_Twall);
Array< OneD, Array<OneD, NekDouble > > scalarVariables(nScalars);
//Array< OneD, Array<OneD, NekDouble > > uplus(nScalars);
// Extract internal values of the scalar variables for Neumann bcs
for (i = 0; i < nScalars; ++i)
{
scalarVariables[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
// uplus[i] = Array<OneD, NekDouble>(nTracePts, 0.0);
// fields[i]->ExtractTracePhys(inarray[i], uplus[i]);
}
// Compute boundary conditions for velocities
for (i = 0; i < nScalars-1; ++i)
{
......@@ -513,9 +509,8 @@ namespace Nektar
Vn, 1, Vn, 1);
}
Array<OneD, NekDouble > qtemp(nTracePts);
Array<OneD, NekDouble > qtemp(nTracePts);
// Evaulate Riemann flux
// qflux = \hat{q} \cdot u = q \cdot n
// Notice: i = 1 (first row of the viscous tensor is zero)
......@@ -534,9 +529,9 @@ namespace Nektar
Vmath::Vmul(nTracePts, m_traceNormals[j], 1, qfluxtemp, 1,
qfluxtemp, 1);
// Extract the physical values of the solution at the boundaries
// Extract the physical values of the solution at the boundaries
fields[i]->ExtractTracePhys(qfield[j][i], qtemp);
// Impose weak boundary condition with flux
if (fields[0]->GetBndCondExpansions().num_elements())
{
......@@ -560,7 +555,7 @@ namespace Nektar
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
const Array<OneD, const NekDouble> &qtemp,
const Array<OneD, const NekDouble> &qtemp,
Array<OneD, NekDouble> &penaltyflux)
{
int cnt = 0;
......@@ -572,11 +567,7 @@ namespace Nektar
int nBndRegions = fields[var]->GetBndCondExpansions().num_elements();
Array<OneD, NekDouble > uterm(nTracePts);
//Array<OneD, NekDouble > qtemp(nTracePts);
// Extract the physical values of the solution at the boundaries
//fields[var]->ExtractTracePhys(qfield, qtemp);
// Loop on the boundary regions to apply appropriate bcs
for (i = 0; i < nBndRegions; ++i)
{
......
......@@ -93,7 +93,7 @@ namespace Nektar
virtual void v_WeakPenaltyO1(
const Array<OneD, MultiRegions::ExpListSharedPtr> &fields,
const Array<OneD, Array<OneD, NekDouble> > &inarray,
const Array<OneD, Array<OneD, NekDouble> > &uplus,
const Array<OneD, Array<OneD, NekDouble> > &uplus,
Array<OneD, Array<OneD, NekDouble> > &penaltyfluxO1);
virtual void v_NumericalFluxO2(
......@@ -107,7 +107,7 @@ namespace Nektar
const int var,
const int dir,
const Array<OneD, const NekDouble> &qfield,
const Array<OneD, const NekDouble> &qtemp,
const Array<OneD, const NekDouble> &qtemp,
Array<OneD, NekDouble> &penaltyflux);
virtual void v_SetHomoDerivs(
......
......@@ -406,7 +406,6 @@ namespace Nektar
Array<OneD, Array<OneD, NekDouble> > &physarray)
{
int i;
int nTracePts = GetTraceTotPoints();
int nVariables = physarray.num_elements();
const Array<OneD, const int> &traceBndMap
......@@ -496,7 +495,6 @@ namespace Nektar
Array<OneD, Array<OneD, NekDouble> > &physarray)
{
int i;
int nTracePts = GetTraceTotPoints();
int nVariables = physarray.num_elements();
const Array<OneD, const int> &traceBndMap
......@@ -564,7 +562,6 @@ namespace Nektar
Array<OneD, Array<OneD, NekDouble> > &physarray)
{
int i;
int nTracePts = GetTraceTotPoints();
int nVariables = physarray.num_elements();
const Array<OneD, const int> &traceBndMap
......@@ -654,7 +651,6 @@ namespace Nektar
{
int i, j;
int nTracePts = GetTraceTotPoints();
int nVariables = physarray.num_elements();
int nDimensions = m_spacedim;
const Array<OneD, const int> &traceBndMap
......@@ -1523,7 +1519,6 @@ namespace Nektar
int i, j;
int e, pnt;
int id1, id2, nBCEdgePts;
int nTracePts = GetTraceTotPoints();
int nVariables = physarray.num_elements();
int nDimensions = m_spacedim;
......
......@@ -166,51 +166,50 @@ namespace Nektar
int &cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &inarray);
void WallBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void WallViscousBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void SymmetryBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void RiemannInvariantBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void PressureOutflowNonReflectiveBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void PressureOutflowBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void PressureOutflowFileBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void PressureInflowFileBC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void ExtrapOrder0BC(
int bcRegion,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
void GetVelocityVector(
const Array<OneD, Array<OneD, NekDouble> > &physfield,
......
......@@ -225,14 +225,14 @@ namespace Nektar
int cnt = 0;
int nvariables = inarray.num_elements();
int nTracePts = GetTraceTotPoints();
Array<OneD, Array<OneD, NekDouble> > Fwd(nvariables);
for (int i = 0; i < nvariables; ++i)
{
Fwd[i] = Array<OneD, NekDouble>(nTracePts);
m_fields[i]->ExtractTracePhys(inarray[i], Fwd[i]);
}
std::string userDefStr;
int nreg = m_fields[0]->GetBndConditions().num_elements();
// Loop over Boundary Regions
......@@ -249,11 +249,11 @@ namespace Nektar
else if(boost::iequals(userDefStr, "IsentropicVortex"))
{
// Isentropic Vortex Boundary Condition
SetBoundaryIsentropicVortex(n, time, cnt, inarray);
SetBoundaryIsentropicVortex(n, time, cnt, Fwd, inarray);
}
else if (boost::iequals(userDefStr,"RinglebFlow"))
{
SetBoundaryRinglebFlow(n, time, cnt, inarray);
SetBoundaryRinglebFlow(n, time, cnt, Fwd, inarray);
}
else
{
......@@ -399,21 +399,13 @@ namespace Nektar
void EulerCFE::SetBoundaryIsentropicVortex(
int bcRegion,
NekDouble time,
int cnt,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray)
{
int nvariables = physarray.num_elements();
int nTraceNumPoints = GetTraceTotPoints();
Array<OneD, Array<OneD, NekDouble> > Fwd(nvariables);
const Array<OneD, const int> &bndTraceMap = m_fields[0]->GetTraceBndMap();
// Get physical values of the forward trace (from exp to phys)
for (int i = 0; i < nvariables; ++i)
{
Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
m_fields[i]->ExtractTracePhys(physarray[i], Fwd[i]);
}
int id2, e_max;
e_max = m_fields[0]->GetBndCondExpansions()[bcRegion]->GetExpSize();
......@@ -845,11 +837,10 @@ namespace Nektar
int bcRegion,
NekDouble time,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray)
{
int nvariables = physarray.num_elements();
int nTraceNumPoints = GetTraceTotPoints();
Array<OneD, Array<OneD, NekDouble> > Fwd(nvariables);
// For 3DHomogenoeus1D
int n_planes = 1;
......@@ -858,14 +849,6 @@ namespace Nektar
int nPointsTot = m_fields[0]->GetTotPoints();
int nPointsTot_plane = m_fields[0]->GetPlane(0)->GetTotPoints();
n_planes = nPointsTot/nPointsTot_plane;
nTraceNumPoints = nTraceNumPoints * n_planes;
}
// Get physical values of the forward trace (from exp to phys)
for (int i = 0; i < nvariables; ++i)
{
Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
m_fields[i]->ExtractTracePhys(physarray[i], Fwd[i]);
}
int id2, id2_plane, e_max;
......
......@@ -123,7 +123,9 @@ namespace Nektar
void SetBoundaryIsentropicVortex(
int bcRegion,
NekDouble time,
int cnt, Array<OneD, Array<OneD, NekDouble> > &physarray);
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
/// Ringleb Flow Test Case.
void GetExactRinglebFlow(
......@@ -135,6 +137,7 @@ namespace Nektar
int bcRegion,
NekDouble time,
int cnt,
Array<OneD, Array<OneD, NekDouble> > &Fwd,
Array<OneD, Array<OneD, NekDouble> > &physarray);
};
}
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment