Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Static Public Member Functions | Static Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::SubSteppingExtrapolate Class Reference

#include <SubSteppingExtrapolate.h>

Inheritance diagram for Nektar::SubSteppingExtrapolate:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SubSteppingExtrapolate:
Collaboration graph
[legend]

Public Member Functions

 SubSteppingExtrapolate (const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
virtual ~SubSteppingExtrapolate ()
- Public Member Functions inherited from Nektar::Extrapolate
 Extrapolate (const LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields, MultiRegions::ExpListSharedPtr pPressure, const Array< OneD, int > pVel, const SolverUtils::AdvectionSharedPtr advObject)
virtual ~Extrapolate ()
void GenerateHOPBCMap ()
void SubSteppingTimeIntegration (const int intMethod, const LibUtilities::TimeIntegrationWrapperSharedPtr &IntegrationScheme)
void SubStepSaveFields (const int nstep)
void SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, const NekDouble Aii_DT, NekDouble kinvis)
void SubStepAdvance (const LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln, const int nstep, NekDouble time)
void MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)
void EvaluatePressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
Array< OneD, NekDouble > GetMaxStdVelocity (const Array< OneD, Array< OneD, NekDouble > > inarray)

Static Public Member Functions

static ExtrapolateSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, MultiRegions::ExpListSharedPtr &pPressure, const Array< OneD, int > &pVel, const SolverUtils::AdvectionSharedPtr &advObject)
 Creates an instance of this class.

Static Public Attributes

static std::string className
 Name of class.

Protected Member Functions

virtual void v_SubSteppingTimeIntegration (int intMethod, const LibUtilities::TimeIntegrationWrapperSharedPtr &IntegrationScheme)
virtual void v_SubStepSaveFields (int nstep)
virtual void v_SubStepSetPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, NekDouble Aii_Dt, NekDouble kinvis)
virtual void v_SubStepAdvance (const LibUtilities::TimeIntegrationSolutionSharedPtr &integrationSoln, int nstep, NekDouble time)
virtual void v_MountHOPBCs (int HBCdata, NekDouble kinvis, Array< OneD, NekDouble > &Q, Array< OneD, const NekDouble > &Advection)
void AddDuDt (const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble Aii_Dt)
void AddDuDt2D (const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble Aii_Dt)
void AddDuDt3D (const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble Aii_Dt)
void SubStepAdvection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void SubStepProjection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
void SubStepExtrapolateField (NekDouble toff, Array< OneD, Array< OneD, NekDouble > > &ExtVel)
void AddAdvectionPenaltyFlux (const Array< OneD, const Array< OneD, NekDouble > > &velfield, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &outarray)
NekDouble GetSubstepTimeStep ()
- Protected Member Functions inherited from Nektar::Extrapolate
void CalcNeumannPressureBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
void CalcOutflowBCs (const Array< OneD, const Array< OneD, NekDouble > > &fields, const Array< OneD, const Array< OneD, NekDouble > > &N, NekDouble kinvis)
void RollOver (Array< OneD, Array< OneD, NekDouble > > &input)
void CurlCurl (Array< OneD, Array< OneD, const NekDouble > > &Vel, Array< OneD, Array< OneD, NekDouble > > &Q, const int j)

Protected Attributes

LibUtilities::TimeIntegrationWrapperSharedPtr m_subStepIntegrationScheme
LibUtilities::TimeIntegrationSchemeOperators m_subStepIntegrationOps
Array< OneD, Array< OneD,
NekDouble > > 		m_previousVelFields
NekDouble m_cflSafetyFactor
int m_infosteps
int minsubsteps
- Protected Attributes inherited from Nektar::Extrapolate
LibUtilities::SessionReaderSharedPtr m_session
LibUtilities::CommSharedPtr m_comm
Array< OneD,
MultiRegions::ExpListSharedPtr > 		m_fields
MultiRegions::ExpListSharedPtr m_pressure
 Pointer to field holding pressure field.
Array< OneD, int > m_velocity
 int which identifies which components of m_fields contains the velocity (u,v,w);
SolverUtils::AdvectionSharedPtr m_advObject
Array< OneD, Array< OneD,
NekDouble > > 		m_previousVelFields
int m_curl_dim
 Curl-curl dimensionality.
int m_bnd_dim
 bounday dimensionality
Array< OneD, const
SpatialDomains::BoundaryConditionShPtr > 		m_PBndConds
 pressure boundary conditions container
Array< OneD,
MultiRegions::ExpListSharedPtr > 		m_PBndExp
 pressure boundary conditions expansion container
int m_pressureCalls
 number of times the high-order pressure BCs have been called
int m_pressureBCsMaxPts
 Maximum points used in pressure BC evaluation.
int m_pressureBCsElmtMaxPts
 Maximum points used in Element adjacent to pressure BC evaluation.
int m_intSteps
 Maximum points used in pressure BC evaluation.
NekDouble m_timestep
bool m_SingleMode
 Flag to determine if single homogeneous mode is used.
bool m_HalfMode
 Flag to determine if half homogeneous mode is used.
bool m_MultipleModes
 Flag to determine if use multiple homogenenous modes are used.
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous)
int m_npointsX
 number of points in X direction (if homogeneous)
int m_npointsY
 number of points in Y direction (if homogeneous)
int m_npointsZ
 number of points in Z direction (if homogeneous)
Array< OneD, int > m_pressureBCtoElmtID
 Id of element to which pressure boundary condition belongs.
Array< OneD, int > m_pressureBCtoTraceID
 Id of edge (2D) or face (3D) to which pressure boundary condition belongs.
Array< OneD, Array< OneD,
NekDouble > > 		m_pressureHBCs
 Storage for current and previous levels of high order pressure boundary conditions.
Array< OneD, Array< OneD,
NekDouble > > 		m_acceleration
 Storage for current and previous levels of the acceleration term.
Array< OneD, HBCInfo > m_HBCdata
 data structure to old all the information regarding High order pressure BCs
Array< OneD, NekDouble > m_wavenumber
 wave number 2 pi k /Lz
Array< OneD, NekDouble > m_negWavenumberSq
 minus Square of wavenumber
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 		m_outflowVel
 Storage for current and previous velocity fields at the otuflow for high order outflow BCs.

Detailed Description

Definition at line 58 of file SubSteppingExtrapolate.h.

Constructor & Destructor Documentation

Nektar::SubSteppingExtrapolate::SubSteppingExtrapolate ( const LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr >  pFields,
MultiRegions::ExpListSharedPtr  pPressure,
const Array< OneD, int >  pVel,
const SolverUtils::AdvectionSharedPtr  advObject 
)

Definition at line 49 of file SubSteppingExtrapolate.cpp.

References m_cflSafetyFactor, m_infosteps, Nektar::Extrapolate::m_session, and minsubsteps.

: Extrapolate(pSession,pFields,pPressure,pVel,advObject)
{
m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
m_session->LoadParameter("CFL", m_cflSafetyFactor, 0.5);
m_session->LoadParameter("MinSubSteps", minsubsteps,0);
}
Nektar::SubSteppingExtrapolate::~SubSteppingExtrapolate ( )
virtual

Definition at line 62 of file SubSteppingExtrapolate.cpp.

{
}

Member Function Documentation

void Nektar::SubSteppingExtrapolate::AddAdvectionPenaltyFlux ( const Array< OneD, const Array< OneD, NekDouble > > &  velfield,
const Array< OneD, const Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protected

Number of trace points

Number of spatial dimensions

Forward state array

Backward state array

upwind numerical flux state array

Normal velocity array

Extract forwards/backwards trace spaces Note: Needs to have correct i value to get boundary conditions

Upwind between elements

Construct difference between numflux and Fwd,Bwd

Calculate the numerical fluxes multipling Fwd, Bwd and numflux by the normal advection velocity

Definition at line 390 of file SubSteppingExtrapolate.cpp.

References ASSERTL1, Nektar::Extrapolate::m_curl_dim, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_velocity, Vmath::Vmul(), Vmath::Vsub(), and Vmath::Vvtvp().

Referenced by SubStepAdvection().

{
ASSERTL1(
physfield.num_elements() == Outarray.num_elements(),
"Physfield and outarray are of different dimensions");
int i;
/// Number of trace points
int nTracePts = m_fields[0]->GetTrace()->GetNpoints();
/// Number of spatial dimensions
int nDimensions = m_curl_dim;
Array<OneD, Array<OneD, NekDouble> > traceNormals(m_curl_dim);
for(i = 0; i < nDimensions; ++i)
{
traceNormals[i] = Array<OneD, NekDouble> (nTracePts);
}
//trace normals
m_fields[0]->GetTrace()->GetNormals(traceNormals);
/// Forward state array
Array<OneD, NekDouble> Fwd(3*nTracePts);
/// Backward state array
Array<OneD, NekDouble> Bwd = Fwd + nTracePts;
/// upwind numerical flux state array
Array<OneD, NekDouble> numflux = Bwd + nTracePts;
/// Normal velocity array
Array<OneD, NekDouble> Vn (nTracePts, 0.0);
// Extract velocity field along the trace space and multiply by trace normals
for(i = 0; i < nDimensions; ++i)
{
m_fields[0]->ExtractTracePhys(m_fields[m_velocity[i]]->GetPhys(), Fwd);
Vmath::Vvtvp(nTracePts, traceNormals[i], 1, Fwd, 1, Vn, 1, Vn, 1);
}
for(i = 0; i < physfield.num_elements(); ++i)
{
/// Extract forwards/backwards trace spaces
/// Note: Needs to have correct i value to get boundary conditions
m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
/// Upwind between elements
m_fields[0]->GetTrace()->Upwind(Vn, Fwd, Bwd, numflux);
/// Construct difference between numflux and Fwd,Bwd
Vmath::Vsub(nTracePts, numflux, 1, Fwd, 1, Fwd, 1);
Vmath::Vsub(nTracePts, numflux, 1, Bwd, 1, Bwd, 1);
/// Calculate the numerical fluxes multipling Fwd, Bwd and
/// numflux by the normal advection velocity
Vmath::Vmul(nTracePts, Fwd, 1, Vn, 1, Fwd, 1);
Vmath::Vmul(nTracePts, Bwd, 1, Vn, 1, Bwd, 1);
m_fields[0]->AddFwdBwdTraceIntegral(Fwd,Bwd,Outarray[i]);
}
}
void Nektar::SubSteppingExtrapolate::AddDuDt ( const Array< OneD, const Array< OneD, NekDouble > > &  N,
NekDouble  Aii_Dt 
)
protected

Definition at line 515 of file SubSteppingExtrapolate.cpp.

References AddDuDt2D(), AddDuDt3D(), ASSERTL0, and Nektar::Extrapolate::m_velocity.

Referenced by v_SubStepSetPressureBCs().

{
switch(m_velocity.num_elements())
{
case 1:
ASSERTL0(
false,
"Velocity correction scheme not designed to have just one velocity component");
break;
case 2:
AddDuDt2D(N,Aii_Dt);
break;
case 3:
AddDuDt3D(N,Aii_Dt);
break;
}
}
void Nektar::SubSteppingExtrapolate::AddDuDt2D ( const Array< OneD, const Array< OneD, NekDouble > > &  N,
NekDouble  Aii_Dt 
)
protected

Definition at line 538 of file SubSteppingExtrapolate.cpp.

References ASSERTL0, Nektar::SpatialDomains::eHigh, Nektar::SpatialDomains::eNoUserDefined, Nektar::SpatialDomains::eTimeDependent, Nektar::StdRegions::StdExpansion::GetTotPoints(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_pressureBCsMaxPts, Nektar::Extrapolate::m_pressureBCtoElmtID, Nektar::Extrapolate::m_pressureBCtoTraceID, Nektar::Extrapolate::m_velocity, and Vmath::Vsub().

Referenced by AddDuDt().

{
int i,n;
ASSERTL0(m_velocity.num_elements() == 2," Routine currently only set up for 2D");
int pindex=2;
Array<OneD, const SpatialDomains::BoundaryConditionShPtr > PBndConds;
Array<OneD, MultiRegions::ExpListSharedPtr> PBndExp,UBndExp,VBndExp;
PBndConds = m_fields[pindex]->GetBndConditions();
PBndExp = m_fields[pindex]->GetBndCondExpansions();
UBndExp = m_fields[m_velocity[0]]->GetBndCondExpansions();
VBndExp = m_fields[m_velocity[1]]->GetBndCondExpansions();
StdRegions::StdExpansionSharedPtr elmt;
StdRegions::StdExpansion1DSharedPtr Pbc;
int cnt,elmtid,nq,offset, boundary,ncoeffs;
Array<OneD, NekDouble> N1(m_pressureBCsMaxPts), N2(m_pressureBCsMaxPts);
Array<OneD, NekDouble> ubc(m_pressureBCsMaxPts),vbc(m_pressureBCsMaxPts);
Array<OneD, NekDouble> Pvals(m_pressureBCsMaxPts);
Array<OneD, NekDouble> Nu,Nv,Ptmp;
for(cnt = n = 0; n < PBndConds.num_elements(); ++n)
{
SpatialDomains::BndUserDefinedType type = PBndConds[n]->GetUserDefined();
if(type == SpatialDomains::eHigh)
{
for(i = 0; i < PBndExp[n]->GetExpSize(); ++i,cnt++)
{
// find element and edge of this expansion.
elmtid = m_pressureBCtoElmtID[cnt];
elmt = m_fields[0]->GetExp(elmtid);
offset = m_fields[0]->GetPhys_Offset(elmtid);
Nu = N[0] + offset;
Nv = N[1] + offset;
Pbc = boost::dynamic_pointer_cast<StdRegions::StdExpansion1D> (PBndExp[n]->GetExp(i));
nq = Pbc->GetTotPoints();
ncoeffs = Pbc->GetNcoeffs();
boundary = m_pressureBCtoTraceID[cnt];
// Get velocity bc
UBndExp[n]->GetExp(i)->BwdTrans(UBndExp[n]->GetCoeffs() + UBndExp[n]->GetCoeff_Offset(i),ubc);
VBndExp[n]->GetExp(i)->BwdTrans(VBndExp[n]->GetCoeffs() + VBndExp[n]->GetCoeff_Offset(i),vbc);
// Get edge values and put into Nu,Nv
elmt->GetEdgePhysVals(boundary,Pbc,Nu,N1);
elmt->GetEdgePhysVals(boundary,Pbc,Nv,N2);
// Take different as Forward Euler but N1,N2
// actually contain the integration of the
// previous steps from the time integration
// scheme.
Vmath::Vsub(nq,ubc,1,N1,1,ubc,1);
Vmath::Vsub(nq,vbc,1,N2,1,vbc,1);
// Divide by aii_Dt to get correct Du/Dt. This is
// because all coefficients in the integration
// scheme are normalised so u^{n+1} has unit
// coefficient and N is already multiplied by
// local coefficient when taken from integration
// scheme
Blas::Dscal(nq,1.0/Aii_Dt,&ubc[0],1);
Blas::Dscal(nq,1.0/Aii_Dt,&vbc[0],1);
// subtrace off du/dt derivative
Pbc->NormVectorIProductWRTBase(ubc,vbc,Pvals);
Vmath::Vsub(ncoeffs,Ptmp = PBndExp[n]->UpdateCoeffs()
+PBndExp[n]->GetCoeff_Offset(i),1, Pvals,1,
Ptmp = PBndExp[n]->UpdateCoeffs()+PBndExp[n]->GetCoeff_Offset(i),1);
}
}
// setting if just standard BC no High order
else if(type == SpatialDomains::eNoUserDefined || type == SpatialDomains::eTimeDependent)
{
cnt += PBndExp[n]->GetExpSize();
}
else
{
ASSERTL0(false,"Unknown USERDEFINEDTYPE in pressure boundary condition");
}
}
}
void Nektar::SubSteppingExtrapolate::AddDuDt3D ( const Array< OneD, const Array< OneD, NekDouble > > &  N,
NekDouble  Aii_Dt 
)
protected

Definition at line 638 of file SubSteppingExtrapolate.cpp.

References ASSERTL0, Nektar::SpatialDomains::eHigh, Nektar::SpatialDomains::eNoUserDefined, Nektar::SpatialDomains::eTimeDependent, Nektar::StdRegions::StdExpansion::GetTotPoints(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_pressureBCsMaxPts, Nektar::Extrapolate::m_pressureBCtoElmtID, Nektar::Extrapolate::m_pressureBCtoTraceID, Nektar::Extrapolate::m_velocity, and Vmath::Vsub().

Referenced by AddDuDt().

{
int i,n;
ASSERTL0(m_velocity.num_elements() == 3," Routine currently only set up for 3D");
int pindex=3;
Array<OneD, const SpatialDomains::BoundaryConditionShPtr > PBndConds;
Array<OneD, MultiRegions::ExpListSharedPtr> PBndExp,UBndExp,VBndExp,WBndExp;
PBndConds = m_fields[pindex]->GetBndConditions();
PBndExp = m_fields[pindex]->GetBndCondExpansions();
UBndExp = m_fields[m_velocity[0]]->GetBndCondExpansions();
VBndExp = m_fields[m_velocity[1]]->GetBndCondExpansions();
WBndExp = m_fields[m_velocity[2]]->GetBndCondExpansions();
StdRegions::StdExpansionSharedPtr elmt;
StdRegions::StdExpansion2DSharedPtr Pbc;
int cnt,elmtid,nq,offset, boundary,ncoeffs;
Array<OneD, NekDouble> N1(m_pressureBCsMaxPts), N2(m_pressureBCsMaxPts),
N3(m_pressureBCsMaxPts);
Array<OneD, NekDouble> ubc(m_pressureBCsMaxPts),vbc(m_pressureBCsMaxPts),
wbc(m_pressureBCsMaxPts);
Array<OneD, NekDouble> Pvals(m_pressureBCsMaxPts);
Array<OneD, NekDouble> Nu,Nv,Nw,Ptmp;
for(cnt = n = 0; n < PBndConds.num_elements(); ++n)
{
SpatialDomains::BndUserDefinedType type = PBndConds[n]->GetUserDefined();
if(type == SpatialDomains::eHigh)
{
for(i = 0; i < PBndExp[n]->GetExpSize(); ++i,cnt++)
{
// find element and face of this expansion.
elmtid = m_pressureBCtoElmtID[cnt];
elmt = m_fields[0]->GetExp(elmtid);
offset = m_fields[0]->GetPhys_Offset(elmtid);
Nu = N[0] + offset;
Nv = N[1] + offset;
Nw = N[2] + offset;
Pbc = boost::dynamic_pointer_cast<StdRegions::StdExpansion2D> (PBndExp[n]->GetExp(i));
nq = Pbc->GetTotPoints();
ncoeffs = Pbc->GetNcoeffs();
boundary = m_pressureBCtoTraceID[cnt];
// Get velocity bc
UBndExp[n]->GetExp(i)->BwdTrans(UBndExp[n]->GetCoeffs() +
UBndExp[n]->GetCoeff_Offset(i),ubc);
VBndExp[n]->GetExp(i)->BwdTrans(VBndExp[n]->GetCoeffs() +
VBndExp[n]->GetCoeff_Offset(i),vbc);
WBndExp[n]->GetExp(i)->BwdTrans(WBndExp[n]->GetCoeffs() +
WBndExp[n]->GetCoeff_Offset(i),wbc);
// Get edge values and put into N1,N2,N3
elmt->GetFacePhysVals(boundary,Pbc,Nu,N1);
elmt->GetFacePhysVals(boundary,Pbc,Nv,N2);
elmt->GetFacePhysVals(boundary,Pbc,Nw,N3);
// Take different as Forward Euler but N1,N2,N3
// actually contain the integration of the
// previous steps from the time integration
// scheme.
Vmath::Vsub(nq,ubc,1,N1,1,ubc,1);
Vmath::Vsub(nq,vbc,1,N2,1,vbc,1);
Vmath::Vsub(nq,wbc,1,N3,1,wbc,1);
// Divide by aii_Dt to get correct Du/Dt. This is
// because all coefficients in the integration
// scheme are normalised so u^{n+1} has unit
// coefficient and N is already multiplied by
// local coefficient when taken from integration
// scheme
Blas::Dscal(nq,1.0/Aii_Dt,&ubc[0],1);
Blas::Dscal(nq,1.0/Aii_Dt,&vbc[0],1);
Blas::Dscal(nq,1.0/Aii_Dt,&wbc[0],1);
// subtrace off du/dt derivative
Pbc->NormVectorIProductWRTBase(ubc,vbc,wbc,Pvals);
Vmath::Vsub(
ncoeffs,
Ptmp = PBndExp[n]->UpdateCoeffs()+PBndExp[n]->GetCoeff_Offset(i),
1,Pvals,1,
Ptmp = PBndExp[n]->UpdateCoeffs()+PBndExp[n]->GetCoeff_Offset(i),
1);
}
}
// setting if just standard BC no High order
else if(type == SpatialDomains::eNoUserDefined || type == SpatialDomains::eTimeDependent)
{
cnt += PBndExp[n]->GetExpSize();
}
else
{
ASSERTL0(false,"Unknown USERDEFINEDTYPE in pressure boundary condition");
}
}
}
static ExtrapolateSharedPtr Nektar::SubSteppingExtrapolate::create ( const LibUtilities::SessionReaderSharedPtr &  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr > &  pFields,
MultiRegions::ExpListSharedPtr &  pPressure,
const Array< OneD, int > &  pVel,
const SolverUtils::AdvectionSharedPtr &  advObject 
)
inlinestatic

Creates an instance of this class.

Definition at line 63 of file SubSteppingExtrapolate.h.

{
ExtrapolateSharedPtr p = MemoryManager<SubSteppingExtrapolate>
::AllocateSharedPtr(pSession,pFields,pPressure,pVel,advObject);
return p;
}
NekDouble Nektar::SubSteppingExtrapolate::GetSubstepTimeStep ( )
protected

Definition at line 355 of file SubSteppingExtrapolate.cpp.

References Nektar::Extrapolate::GetMaxStdVelocity(), m_cflSafetyFactor, Nektar::Extrapolate::m_comm, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_velocity, Nektar::LibUtilities::ReduceMin, and Vmath::Vmin().

Referenced by v_SubStepAdvance().

{
int n_element = m_fields[0]->GetExpSize();
const Array<OneD, int> ExpOrder=m_fields[0]->EvalBasisNumModesMaxPerExp();
Array<OneD, int> ExpOrderList (n_element, ExpOrder);
const NekDouble cLambda = 0.2; // Spencer book pag. 317
Array<OneD, NekDouble> tstep (n_element, 0.0);
Array<OneD, NekDouble> stdVelocity(n_element, 0.0);
Array<OneD, Array<OneD, NekDouble> > velfields(m_velocity.num_elements());
for(int i = 0; i < m_velocity.num_elements(); ++i)
{
velfields[i] = m_fields[m_velocity[i]]->UpdatePhys();
}
stdVelocity = GetMaxStdVelocity(velfields);
for(int el = 0; el < n_element; ++el)
{
tstep[el] = m_cflSafetyFactor /
(stdVelocity[el] * cLambda *
(ExpOrder[el]-1) * (ExpOrder[el]-1));
}
NekDouble TimeStep = Vmath::Vmin(n_element, tstep, 1);
m_comm->AllReduce(TimeStep,LibUtilities::ReduceMin);
return TimeStep;
}
void Nektar::SubSteppingExtrapolate::SubStepAdvection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Explicit Advection terms used by SubStepAdvance time integration

Get the number of coefficients

Define an auxiliary variable to compute the RHS

Operations to compute the RHS

Multiply the flux by the inverse of the mass matrix

Store in outarray the physical values of the RHS

Definition at line 125 of file SubSteppingExtrapolate.cpp.

References AddAdvectionPenaltyFlux(), Nektar::Extrapolate::m_advObject, Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, Vmath::Neg(), and SubStepExtrapolateField().

Referenced by v_SubSteppingTimeIntegration().

{
int i;
int nVariables = inarray.num_elements();
int nQuadraturePts = inarray[0].num_elements();
/// Get the number of coefficients
int ncoeffs = m_fields[0]->GetNcoeffs();
/// Define an auxiliary variable to compute the RHS
Array<OneD, Array<OneD, NekDouble> > WeakAdv(nVariables);
WeakAdv[0] = Array<OneD, NekDouble> (ncoeffs*nVariables);
for(i = 1; i < nVariables; ++i)
{
WeakAdv[i] = WeakAdv[i-1] + ncoeffs;
}
Array<OneD, Array<OneD, NekDouble> > Velfields(m_velocity.num_elements());
Array<OneD, int> VelIds(m_velocity.num_elements());
Velfields[0] = Array<OneD, NekDouble> (nQuadraturePts*m_velocity.num_elements());
for(i = 1; i < m_velocity.num_elements(); ++i)
{
Velfields[i] = Velfields[i-1] + nQuadraturePts;
}
SubStepExtrapolateField(fmod(time,m_timestep), Velfields);
m_advObject->Advect(m_velocity.num_elements(), m_fields, Velfields, inarray, outarray, time);
for(i = 0; i < nVariables; ++i)
{
m_fields[i]->IProductWRTBase(outarray[i],WeakAdv[i]);
// negation requried due to sign of DoAdvection term to be consistent
Vmath::Neg(ncoeffs, WeakAdv[i], 1);
}
AddAdvectionPenaltyFlux(Velfields, inarray, WeakAdv);
/// Operations to compute the RHS
for(i = 0; i < nVariables; ++i)
{
// Negate the RHS
Vmath::Neg(ncoeffs, WeakAdv[i], 1);
/// Multiply the flux by the inverse of the mass matrix
m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
/// Store in outarray the physical values of the RHS
m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
}
//add the force
/*if(m_session->DefinesFunction("BodyForce"))
{
if(m_SingleMode || m_HalfMode)
{
for(int i = 0; i < m_nConvectiveFields; ++i)
{
m_forces[i]->SetWaveSpace(true);
m_forces[i]->BwdTrans(m_forces[i]->GetCoeffs(),
m_forces[i]->UpdatePhys());
}
}
int nqtot = m_fields[0]->GetTotPoints();
for(int i = 0; i < m_nConvectiveFields; ++i)
{
Vmath::Vadd(nqtot,outarray[i],1,(m_forces[i]->GetPhys()),1,outarray[i],1);
}
}*/
}
void Nektar::SubSteppingExtrapolate::SubStepExtrapolateField ( NekDouble  toff,
Array< OneD, Array< OneD, NekDouble > > &  ExtVel 
)
protected

Extrapolate field using equally time spaced field un,un-1,un-2, (at dt intervals) to time n+t at order Ord

Definition at line 462 of file SubSteppingExtrapolate.cpp.

References Vmath::Fill(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_intSteps, m_previousVelFields, Nektar::Extrapolate::m_timestep, Nektar::Extrapolate::m_velocity, npts, and Vmath::Smul().

Referenced by SubStepAdvection().

{
int npts = m_fields[0]->GetTotPoints();
int nvel = m_velocity.num_elements();
int i,j;
Array<OneD, NekDouble> l(4);
int ord = m_intSteps;
// calculate Lagrange interpolants
Vmath::Fill(4,1.0,l,1);
for(i = 0; i <= ord; ++i)
{
for(j = 0; j <= ord; ++j)
{
if(i != j)
{
l[i] *= (j*m_timestep+toff);
l[i] /= (j*m_timestep-i*m_timestep);
}
}
}
for(i = 0; i < nvel; ++i)
{
Vmath::Smul(npts,l[0],m_previousVelFields[i],1,ExtVel[i],1);
for(j = 1; j <= ord; ++j)
{
Blas::Daxpy(npts,l[j],m_previousVelFields[j*nvel+i],1,
ExtVel[i],1);
}
}
}
void Nektar::SubSteppingExtrapolate::SubStepProjection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Projection used by SubStepAdvance time integration

Definition at line 205 of file SubSteppingExtrapolate.cpp.

References ASSERTL1, and Vmath::Vcopy().

Referenced by v_SubSteppingTimeIntegration().

{
ASSERTL1(inarray.num_elements() == outarray.num_elements(),"Inarray and outarray of different sizes ");
for(int i = 0; i < inarray.num_elements(); ++i)
{
Vmath::Vcopy(inarray[i].num_elements(),inarray[i],1,outarray[i],1);
}
}
void Nektar::SubSteppingExtrapolate::v_MountHOPBCs ( int  HBCdata,
NekDouble  kinvis,
Array< OneD, NekDouble > &  Q,
Array< OneD, const NekDouble > &  Advection 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 503 of file SubSteppingExtrapolate.cpp.

References Vmath::Smul().

{
Vmath::Smul(HBCdata,-kinvis,Q,1,Q,1);
}
void Nektar::SubSteppingExtrapolate::v_SubStepAdvance ( const LibUtilities::TimeIntegrationSolutionSharedPtr &  integrationSoln,
int  nstep,
NekDouble  time 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 292 of file SubSteppingExtrapolate.cpp.

References GetSubstepTimeStep(), m_cflSafetyFactor, Nektar::Extrapolate::m_comm, Nektar::Extrapolate::m_fields, m_infosteps, Nektar::Extrapolate::m_intSteps, m_subStepIntegrationOps, m_subStepIntegrationScheme, Nektar::Extrapolate::m_timestep, and minsubsteps.

{
int n;
int nsubsteps;
NekDouble dt;
Array<OneD, Array<OneD, NekDouble> > fields, velfields;
static int ncalls = 1;
int nint = min(ncalls++, m_intSteps);
Array<OneD, NekDouble> CFL(m_fields[0]->GetExpSize(),
m_cflSafetyFactor);
//this needs to change
m_comm = m_fields[0]->GetComm()->GetRowComm();
// Get the proper time step with CFL control
dt = GetSubstepTimeStep();
nsubsteps = (m_timestep > dt)? ((int)(m_timestep/dt)+1):1;
nsubsteps = max(minsubsteps, nsubsteps);
dt = m_timestep/nsubsteps;
if (m_infosteps && !((nstep+1)%m_infosteps) && m_comm->GetRank() == 0)
{
cout << "Sub-integrating using "<< nsubsteps
<< " steps over Dt = " << m_timestep
<< " (SubStep CFL=" << m_cflSafetyFactor << ")"<< endl;
}
for (int m = 0; m < nint; ++m)
{
// We need to update the fields held by the m_integrationSoln
fields = integrationSoln->UpdateSolutionVector()[m];
// Initialise NS solver which is set up to use a GLM method
// with calls to EvaluateAdvection_SetPressureBCs and
// SolveUnsteadyStokesSystem
LibUtilities::TimeIntegrationSolutionSharedPtr
SubIntegrationSoln = m_subStepIntegrationScheme->
InitializeScheme(dt, fields, time, m_subStepIntegrationOps);
for(n = 0; n < nsubsteps; ++n)
{
fields = m_subStepIntegrationScheme->TimeIntegrate(
n, dt, SubIntegrationSoln,
m_subStepIntegrationOps);
}
// Reset time integrated solution in m_integrationSoln
integrationSoln->SetSolVector(m,fields);
}
}
void Nektar::SubSteppingExtrapolate::v_SubSteppingTimeIntegration ( int  intMethod,
const LibUtilities::TimeIntegrationWrapperSharedPtr &  IntegrationScheme 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 66 of file SubSteppingExtrapolate.cpp.

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, >::CreateInstance(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineProjection(), Nektar::LibUtilities::eBackwardEuler, Nektar::LibUtilities::eBDFImplicitOrder1, Nektar::LibUtilities::eBDFImplicitOrder2, Nektar::LibUtilities::GetTimeIntegrationWrapperFactory(), Nektar::Extrapolate::m_fields, Nektar::Extrapolate::m_intSteps, m_previousVelFields, m_subStepIntegrationOps, m_subStepIntegrationScheme, Nektar::Extrapolate::m_velocity, SubStepAdvection(), and SubStepProjection().

{
int i;
// Set to 1 for first step and it will then be increased in
// time advance routines
switch(intMethod)
{
case LibUtilities::eBackwardEuler:
case LibUtilities::eBDFImplicitOrder1:
{
m_subStepIntegrationScheme = LibUtilities::GetTimeIntegrationWrapperFactory().CreateInstance("ForwardEuler");
// Fields for linear interpolation
m_previousVelFields = Array<OneD, Array<OneD, NekDouble> >(2*m_fields.num_elements());
int ntotpts = m_fields[0]->GetTotPoints();
m_previousVelFields[0] = Array<OneD, NekDouble>(2*m_fields.num_elements()*ntotpts);
for(i = 1; i < 2*m_fields.num_elements(); ++i)
{
m_previousVelFields[i] = m_previousVelFields[i-1] + ntotpts;
}
}
break;
case LibUtilities::eBDFImplicitOrder2:
{
m_subStepIntegrationScheme = LibUtilities::GetTimeIntegrationWrapperFactory().CreateInstance("RungeKutta2_ImprovedEuler");
int nvel = m_velocity.num_elements();
// Fields for quadratic interpolation
m_previousVelFields = Array<OneD, Array<OneD, NekDouble> >(3*nvel);
int ntotpts = m_fields[0]->GetTotPoints();
m_previousVelFields[0] = Array<OneD, NekDouble>(3*nvel*ntotpts);
for(i = 1; i < 3*nvel; ++i)
{
m_previousVelFields[i] = m_previousVelFields[i-1] + ntotpts;
}
}
break;
default:
ASSERTL0(0,"Integration method not suitable: Options include BackwardEuler or BDFImplicitOrder1");
break;
}
m_intSteps = IntegrationScheme->GetIntegrationSteps();
// set explicit time-integration class operators
m_subStepIntegrationOps.DefineOdeRhs(&SubSteppingExtrapolate::SubStepAdvection, this);
m_subStepIntegrationOps.DefineProjection(&SubSteppingExtrapolate::SubStepProjection, this);
}
void Nektar::SubSteppingExtrapolate::v_SubStepSaveFields ( int  nstep)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 243 of file SubSteppingExtrapolate.cpp.

References Nektar::Extrapolate::m_fields, m_previousVelFields, Nektar::Extrapolate::m_velocity, npts, and Vmath::Vcopy().

{
int i,n;
int nvel = m_velocity.num_elements();
int npts = m_fields[0]->GetTotPoints();
// rotate fields
int nblocks = m_previousVelFields.num_elements()/nvel;
Array<OneD, NekDouble> save;
for(n = 0; n < nvel; ++n)
{
save = m_previousVelFields[(nblocks-1)*nvel+n];
for(i = nblocks-1; i > 0; --i)
{
m_previousVelFields[i*nvel+n] = m_previousVelFields[(i-1)*nvel+n];
}
m_previousVelFields[n] = save;
}
// Put previous field
for(i = 0; i < nvel; ++i)
{
m_fields[m_velocity[i]]->BwdTrans(m_fields[m_velocity[i]]->GetCoeffs(),
m_fields[m_velocity[i]]->UpdatePhys());
Vmath::Vcopy(npts,m_fields[m_velocity[i]]->GetPhys(),1,
m_previousVelFields[i],1);
}
if(nstep == 0)// initialise all levels with first field
{
for(n = 0; n < nvel; ++n)
{
for(i = 1; i < nblocks; ++i)
{
Vmath::Vcopy(npts,m_fields[m_velocity[n]]->GetPhys(),1,
m_previousVelFields[i*nvel+n],1);
}
}
}
}
void Nektar::SubSteppingExtrapolate::v_SubStepSetPressureBCs ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
NekDouble  Aii_Dt,
NekDouble  kinvis 
)
protectedvirtual

Implements Nektar::Extrapolate.

Definition at line 222 of file SubSteppingExtrapolate.cpp.

References AddDuDt(), Nektar::Extrapolate::EvaluatePressureBCs(), Nektar::Extrapolate::m_fields, and Nektar::Extrapolate::m_velocity.

{
int nConvectiveFields =m_fields.num_elements()-1;
Array<OneD, Array<OneD, NekDouble> > velfields(nConvectiveFields);
for(int i = 0; i < nConvectiveFields; ++i)
{
velfields[i] = m_fields[m_velocity[i]]->GetPhys();
}
EvaluatePressureBCs(velfields,inarray,kinvis);
AddDuDt(inarray,Aii_Dt);
}

Member Data Documentation

std::string Nektar::SubSteppingExtrapolate::className
static
Initial value:
GetExtrapolateFactory().RegisterCreatorFunction(
"SubStepping",
SubSteppingExtrapolate::create,
"SubStepping")

Name of class.

Registers the class with the Factory.

Definition at line 76 of file SubSteppingExtrapolate.h.

NekDouble Nektar::SubSteppingExtrapolate::m_cflSafetyFactor
protected

Definition at line 149 of file SubSteppingExtrapolate.h.

Referenced by GetSubstepTimeStep(), SubSteppingExtrapolate(), and v_SubStepAdvance().

int Nektar::SubSteppingExtrapolate::m_infosteps
protected

Definition at line 150 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

Array<OneD, Array<OneD, NekDouble> > Nektar::SubSteppingExtrapolate::m_previousVelFields
protected

Definition at line 147 of file SubSteppingExtrapolate.h.

Referenced by SubStepExtrapolateField(), v_SubSteppingTimeIntegration(), and v_SubStepSaveFields().

LibUtilities::TimeIntegrationSchemeOperators Nektar::SubSteppingExtrapolate::m_subStepIntegrationOps
protected

Definition at line 145 of file SubSteppingExtrapolate.h.

Referenced by v_SubStepAdvance(), and v_SubSteppingTimeIntegration().

LibUtilities::TimeIntegrationWrapperSharedPtr Nektar::SubSteppingExtrapolate::m_subStepIntegrationScheme
protected

Definition at line 144 of file SubSteppingExtrapolate.h.

Referenced by v_SubStepAdvance(), and v_SubSteppingTimeIntegration().

int Nektar::SubSteppingExtrapolate::minsubsteps
protected

Definition at line 151 of file SubSteppingExtrapolate.h.

Referenced by SubSteppingExtrapolate(), and v_SubStepAdvance().

Generated on Wed Mar 4 2015 16:33:34 for Nektar++ by   doxygen 1.8.1.2