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

#include <EulerCFE.h>

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

Public Member Functions

virtual ~EulerCFE ()
 problem type selector More...
 
- Public Member Functions inherited from Nektar::CompressibleFlowSystem
virtual ~CompressibleFlowSystem ()
 Destructor for CompressibleFlowSystem class. More...
 
NekDouble GetStabilityLimit (int n)
 Function to calculate the stability limit for DG/CG. More...
 
Array< OneD, NekDoubleGetStabilityLimitVector (const Array< OneD, int > &ExpOrder)
 Function to calculate the stability limit for DG/CG (a vector of them). More...
 
- Public Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
virtual SOLVER_UTILS_EXPORT ~UnsteadySystem ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray)
 Calculate the larger time-step mantaining the problem stable. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::EquationSystem
virtual SOLVER_UTILS_EXPORT ~EquationSystem ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT void SetUpTraceNormals (void)
 
SOLVER_UTILS_EXPORT void InitObject ()
 Initialises the members of this object. More...
 
SOLVER_UTILS_EXPORT void DoInitialise ()
 Perform any initialisation necessary before solving the problem. More...
 
SOLVER_UTILS_EXPORT void DoSolve ()
 Solve the problem. More...
 
SOLVER_UTILS_EXPORT void TransCoeffToPhys ()
 Transform from coefficient to physical space. More...
 
SOLVER_UTILS_EXPORT void TransPhysToCoeff ()
 Transform from physical to coefficient space. More...
 
SOLVER_UTILS_EXPORT void Output ()
 Perform output operations after solve. More...
 
SOLVER_UTILS_EXPORT NekDouble LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Linf error computation. More...
 
SOLVER_UTILS_EXPORT std::string GetSessionName ()
 Get Session name. More...
 
template<class T >
boost::shared_ptr< T > as ()
 
SOLVER_UTILS_EXPORT void ResetSessionName (std::string newname)
 Reset Session name. More...
 
SOLVER_UTILS_EXPORT
LibUtilities::SessionReaderSharedPtr 
GetSession ()
 Get Session name. More...
 
SOLVER_UTILS_EXPORT
MultiRegions::ExpListSharedPtr 
GetPressure ()
 Get pressure field if available. More...
 
SOLVER_UTILS_EXPORT void PrintSummary (std::ostream &out)
 Print a summary of parameters and solver characteristics. More...
 
SOLVER_UTILS_EXPORT void SetLambda (NekDouble lambda)
 Set parameter m_lambda. More...
 
SOLVER_UTILS_EXPORT void EvaluateFunction (Array< OneD, Array< OneD, NekDouble > > &pArray, std::string pFunctionName, const NekDouble pTime=0.0, const int domain=0)
 Evaluates a function as specified in the session file. More...
 
SOLVER_UTILS_EXPORT void EvaluateFunction (std::vector< std::string > pFieldNames, Array< OneD, Array< OneD, NekDouble > > &pFields, const std::string &pName, const NekDouble &pTime=0.0, const int domain=0)
 Populate given fields with the function from session. More...
 
SOLVER_UTILS_EXPORT void EvaluateFunction (std::vector< std::string > pFieldNames, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const std::string &pName, const NekDouble &pTime=0.0, const int domain=0)
 Populate given fields with the function from session. More...
 
SOLVER_UTILS_EXPORT void EvaluateFunction (std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
 
SOLVER_UTILS_EXPORT std::string DescribeFunction (std::string pFieldName, const std::string &pFunctionName, const int domain)
 Provide a description of a function for a given field name. More...
 
SOLVER_UTILS_EXPORT void InitialiseBaseFlow (Array< OneD, Array< OneD, NekDouble > > &base)
 Perform initialisation of the base flow. More...
 
SOLVER_UTILS_EXPORT void SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 Initialise the data in the dependent fields. More...
 
SOLVER_UTILS_EXPORT void EvaluateExactSolution (int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 Evaluates an exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln, bool Normalised=false)
 Compute the L2 error between fields and a given exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, bool Normalised=false)
 Compute the L2 error of the fields. More...
 
SOLVER_UTILS_EXPORT Array
< OneD, NekDouble
ErrorExtraPoints (unsigned int field)
 Compute error (L2 and L_inf) over an larger set of quadrature points return [L2 Linf]. More...
 
SOLVER_UTILS_EXPORT void WeakAdvectionGreensDivergenceForm (const Array< OneD, Array< OneD, NekDouble > > &F, Array< OneD, NekDouble > &outarray)
 Compute the inner product $ (\nabla \phi \cdot F) $. More...
 
SOLVER_UTILS_EXPORT void WeakAdvectionDivergenceForm (const Array< OneD, Array< OneD, NekDouble > > &F, Array< OneD, NekDouble > &outarray)
 Compute the inner product $ (\phi, \nabla \cdot F) $. More...
 
SOLVER_UTILS_EXPORT void WeakAdvectionNonConservativeForm (const Array< OneD, Array< OneD, NekDouble > > &V, const Array< OneD, const NekDouble > &u, Array< OneD, NekDouble > &outarray, bool UseContCoeffs=false)
 Compute the inner product $ (\phi, V\cdot \nabla u) $. More...
 
f SOLVER_UTILS_EXPORT void AdvectionNonConservativeForm (const Array< OneD, Array< OneD, NekDouble > > &V, const Array< OneD, const NekDouble > &u, Array< OneD, NekDouble > &outarray, Array< OneD, NekDouble > &wk=NullNekDouble1DArray)
 Compute the non-conservative advection. More...
 
SOLVER_UTILS_EXPORT void WeakDGAdvection (const Array< OneD, Array< OneD, NekDouble > > &InField, Array< OneD, Array< OneD, NekDouble > > &OutField, bool NumericalFluxIncludesNormal=true, bool InFieldIsInPhysSpace=false, int nvariables=0)
 Calculate the weak discontinuous Galerkin advection. More...
 
SOLVER_UTILS_EXPORT void WeakDGDiffusion (const Array< OneD, Array< OneD, NekDouble > > &InField, Array< OneD, Array< OneD, NekDouble > > &OutField, bool NumericalFluxIncludesNormal=true, bool InFieldIsInPhysSpace=false)
 Calculate weak DG Diffusion in the LDG form. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n)
 Write checkpoint file of m_fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 Write checkpoint file of custom data fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_BaseFlow (const int n)
 Write base flow file of m_fields. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname)
 Write field data to the given filename. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 Write input fields to the given filename. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Input field data from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFldToMultiDomains (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int ndomains)
 Input field data from the given file to multiple domains. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, std::vector< std::string > &fieldStr, Array< OneD, Array< OneD, NekDouble > > &coeffs)
 Output a field. Input field data into array from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, MultiRegions::ExpListSharedPtr &pField, std::string &pFieldName)
 Output a field. Input field data into ExpList from the given file. More...
 
SOLVER_UTILS_EXPORT void ScanForHistoryPoints ()
 Builds map of which element holds each history point. More...
 
SOLVER_UTILS_EXPORT void WriteHistoryData (std::ostream &out)
 Probe each history point and write to file. More...
 
SOLVER_UTILS_EXPORT void SessionSummary (SummaryList &vSummary)
 Write out a session summary. More...
 
SOLVER_UTILS_EXPORT Array
< OneD,
MultiRegions::ExpListSharedPtr > & 
UpdateFields ()
 
SOLVER_UTILS_EXPORT
LibUtilities::FieldMetaDataMap
UpdateFieldMetaDataMap ()
 Get hold of FieldInfoMap so it can be updated. More...
 
SOLVER_UTILS_EXPORT NekDouble GetFinalTime ()
 Return final time. More...
 
SOLVER_UTILS_EXPORT int GetNcoeffs ()
 
SOLVER_UTILS_EXPORT int GetNcoeffs (const int eid)
 
SOLVER_UTILS_EXPORT int GetNumExpModes ()
 
SOLVER_UTILS_EXPORT const
Array< OneD, int > 
GetNumExpModesPerExp ()
 
SOLVER_UTILS_EXPORT int GetNvariables ()
 
SOLVER_UTILS_EXPORT const
std::string 
GetVariable (unsigned int i)
 
SOLVER_UTILS_EXPORT int GetTraceTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTraceNpoints ()
 
SOLVER_UTILS_EXPORT int GetExpSize ()
 
SOLVER_UTILS_EXPORT int GetPhys_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetCoeff_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTotPoints (int n)
 
SOLVER_UTILS_EXPORT int GetNpoints ()
 
SOLVER_UTILS_EXPORT int GetNumElmVelocity ()
 
SOLVER_UTILS_EXPORT int GetSteps ()
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep ()
 
SOLVER_UTILS_EXPORT void CopyFromPhysField (const int i, Array< OneD, NekDouble > &output)
 
SOLVER_UTILS_EXPORT void CopyToPhysField (const int i, Array< OneD, NekDouble > &output)
 
SOLVER_UTILS_EXPORT void SetSteps (const int steps)
 
SOLVER_UTILS_EXPORT void ZeroPhysFields ()
 
SOLVER_UTILS_EXPORT void FwdTransFields ()
 
SOLVER_UTILS_EXPORT void GetFluxVector (const int i, Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
 
SOLVER_UTILS_EXPORT void GetFluxVector (const int i, Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &fluxX, Array< OneD, Array< OneD, NekDouble > > &fluxY)
 
SOLVER_UTILS_EXPORT void GetFluxVector (const int i, const int j, Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
 
SOLVER_UTILS_EXPORT void NumericalFlux (Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numflux)
 
SOLVER_UTILS_EXPORT void NumericalFlux (Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numfluxX, Array< OneD, Array< OneD, NekDouble > > &numfluxY)
 
SOLVER_UTILS_EXPORT void NumFluxforScalar (const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &uflux)
 
SOLVER_UTILS_EXPORT void NumFluxforVector (const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
 
SOLVER_UTILS_EXPORT void SetModifiedBasis (const bool modbasis)
 
SOLVER_UTILS_EXPORT int NoCaseStringCompare (const std::string &s1, const std::string &s2)
 Perform a case-insensitive string comparison. More...
 
SOLVER_UTILS_EXPORT int GetCheckpointNumber ()
 
SOLVER_UTILS_EXPORT void SetCheckpointNumber (int num)
 
SOLVER_UTILS_EXPORT int GetCheckpointSteps ()
 
SOLVER_UTILS_EXPORT void SetCheckpointSteps (int num)
 
SOLVER_UTILS_EXPORT void SetTime (const NekDouble time)
 
SOLVER_UTILS_EXPORT void SetInitialStep (const int step)
 
SOLVER_UTILS_EXPORT void SetBoundaryConditions (NekDouble time)
 Evaluates the boundary conditions at the given time. More...
 
virtual SOLVER_UTILS_EXPORT bool v_NegatedOp ()
 Virtual function to identify if operator is negated in DoSolve. More...
 

Static Public Member Functions

static
SolverUtils::EquationSystemSharedPtr 
create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 
- Static Public Member Functions inherited from Nektar::CompressibleFlowSystem
static
SolverUtils::EquationSystemSharedPtr 
create (const LibUtilities::SessionReaderSharedPtr &pSession)
 Creates an instance of this class. More...
 

Public Attributes

ProblemType m_problemType
 
- Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem
NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 

Static Public Attributes

static std::string className
 Name of class. More...
 
- Static Public Attributes inherited from Nektar::CompressibleFlowSystem
static std::string className
 Name of class. More...
 

Protected Member Functions

 EulerCFE (const LibUtilities::SessionReaderSharedPtr &pSession)
 
virtual void v_InitObject ()
 Initialization object for CompressibleFlowSystem class. More...
 
virtual void v_GenerateSummary (SolverUtils::SummaryList &s)
 Print a summary of time stepping parameters. More...
 
void DoOdeRhs (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 Compute the right-hand side. More...
 
void DoOdeProjection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 Compute the projection and call the method for imposing the boundary conditions in case of discontinuous projection. More...
 
virtual void v_SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 Set the initial conditions. More...
 
virtual void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time=0.0)
 Get the exact solutions for isentropic vortex and Ringleb flow problems. More...
 
- Protected Member Functions inherited from Nektar::CompressibleFlowSystem
 CompressibleFlowSystem (const LibUtilities::SessionReaderSharedPtr &pSession)
 
void GetFluxVector (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &flux)
 Return the flux vector for the compressible Euler equations. More...
 
void GetFluxVectorDeAlias (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &flux)
 Return the flux vector for the compressible Euler equations by using the de-aliasing technique. More...
 
void GetViscousFluxVector (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &derivatives, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &viscousTensor)
 Return the flux vector for the LDG diffusion problem. More...
 
void GetFluxVectorPDESC (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &flux)
 
void GetViscousFluxVectorDeAlias (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &derivatives, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &viscousTensor)
 Return the flux vector for the LDG diffusion problem. More...
 
void SetCommonBC (const std::string &userDefStr, const int n, const NekDouble time, int &cnt, Array< OneD, Array< OneD, NekDouble > > &inarray)
 Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem level since they are compressible solver specific; b) Time dependent BCs. More...
 
void WallBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Wall boundary conditions for compressible flow problems. More...
 
void WallViscousBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Wall boundary conditions for viscous compressible flow problems. More...
 
void SymmetryBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Symmetry boundary conditions for compressible flow problems. More...
 
void RiemannInvariantBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Outflow characteristic boundary conditions for compressible flow problems. More...
 
void PressureOutflowNonReflectiveBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Pressure outflow non-reflective boundary conditions for compressible flow problems. More...
 
void PressureOutflowBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Pressure outflow boundary conditions for compressible flow problems. More...
 
void PressureOutflowFileBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Pressure outflow boundary conditions for compressible flow problems. More...
 
void PressureInflowFileBC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Pressure inflow boundary conditions for compressible flow problems where either the density and the velocities are assigned from a file or the full state is assigned from a file (depending on the problem type, either subsonic or supersonic). More...
 
void ExtrapOrder0BC (int bcRegion, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Extrapolation of order 0 for all the variables such that, at the boundaries, a trivial Riemann problem is solved. More...
 
void GetVelocityVector (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &velocity)
 
void GetSoundSpeed (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure, Array< OneD, NekDouble > &soundspeed)
 
void GetMach (Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &soundspeed, Array< OneD, NekDouble > &mach)
 
void GetTemperature (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure, Array< OneD, NekDouble > &temperature)
 
void GetPressure (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure)
 Calculate the pressure field $ p = (\gamma-1)(E-\frac{1}{2}\rho\| \mathbf{v} \|^2) $ assuming an ideal gas law. More...
 
void GetPressure (const Array< OneD, const Array< OneD, NekDouble > > &physfield, const Array< OneD, const Array< OneD, NekDouble > > &velocity, Array< OneD, NekDouble > &pressure)
 
void GetEnthalpy (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure, Array< OneD, NekDouble > &enthalpy)
 
void GetEntropy (const Array< OneD, const Array< OneD, NekDouble > > &physfield, const Array< OneD, const NekDouble > &pressure, const Array< OneD, const NekDouble > &temperature, Array< OneD, NekDouble > &entropy)
 
void GetSmoothArtificialViscosity (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &eps_bar)
 
void GetDynamicViscosity (const Array< OneD, const NekDouble > &temperature, Array< OneD, NekDouble > &mu)
 
void GetStdVelocity (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &stdV)
 
virtual bool v_PostIntegrate (int step)
 
bool CalcSteadyState (bool output)
 
void GetSensor (const Array< OneD, const Array< OneD, NekDouble > > &physarray, Array< OneD, NekDouble > &Sensor, Array< OneD, NekDouble > &SensorKappa)
 
void GetElementDimensions (Array< OneD, Array< OneD, NekDouble > > &outarray, Array< OneD, NekDouble > &hmin)
 
void GetAbsoluteVelocity (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &Vtot)
 
void GetArtificialDynamicViscosity (const Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &mu_var)
 
void SetVarPOrderElmt (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &PolyOrder)
 
void GetForcingTerm (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > outarrayForcing)
 
virtual NekDouble v_GetTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray)
 Return the timestep to be used for the next step in the time-marching loop. More...
 
NekDouble GetGasConstant ()
 
NekDouble GetGamma ()
 
const Array< OneD, const Array
< OneD, NekDouble > > & 
GetVecLocs ()
 
const Array< OneD, const Array
< OneD, NekDouble > > & 
GetNormals ()
 
virtual void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 
- Protected Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession)
 Initialises UnsteadySystem class members. More...
 
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator ()
 Get the maximum timestep estimator for cfl control. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoSolve ()
 Solves an unsteady problem. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoInitialise ()
 Sets up initial conditions. More...
 
virtual SOLVER_UTILS_EXPORT void v_AppendOutput1D (Array< OneD, Array< OneD, NekDouble > > &solution1D)
 Print the solution at each solution point in a txt file. More...
 
virtual SOLVER_UTILS_EXPORT void v_NumericalFlux (Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numflux)
 
virtual SOLVER_UTILS_EXPORT void v_NumericalFlux (Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numfluxX, Array< OneD, Array< OneD, NekDouble > > &numfluxY)
 
virtual SOLVER_UTILS_EXPORT void v_NumFluxforScalar (const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &uflux)
 
virtual SOLVER_UTILS_EXPORT void v_NumFluxforVector (const Array< OneD, Array< OneD, NekDouble > > &ufield, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &qfield, Array< OneD, Array< OneD, NekDouble > > &qflux)
 
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_SteadyStateCheck (int step)
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time)
 
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff (const Array< OneD, Array< OneD, NekDouble > > vel, StdRegions::VarCoeffMap &varCoeffMap)
 Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::EquationSystem
SOLVER_UTILS_EXPORT EquationSystem (const LibUtilities::SessionReaderSharedPtr &pSession)
 Initialises EquationSystem class members. More...
 
int nocase_cmp (const std::string &s1, const std::string &s2)
 
virtual SOLVER_UTILS_EXPORT
NekDouble 
v_LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Virtual function for the L_inf error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT
NekDouble 
v_L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the L_2 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
SOLVER_UTILS_EXPORT void SetUpBaseFields (SpatialDomains::MeshGraphSharedPtr &mesh)
 
SOLVER_UTILS_EXPORT void ImportFldBase (std::string pInfile, SpatialDomains::MeshGraphSharedPtr pGraph)
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT
MultiRegions::ExpListSharedPtr 
v_GetPressure (void)
 

Private Member Functions

void SetBoundaryConditions (Array< OneD, Array< OneD, NekDouble > > &physarray, NekDouble time)
 Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem level since they are compressible solver specific; b) Isentropic vortex and Ringleb flow BCs implemented at EulerCFE level since they are Euler solver specific; c) Time dependent BCs. More...
 
void EvaluateIsentropicVortex (const Array< OneD, NekDouble > &x, const Array< OneD, NekDouble > &y, const Array< OneD, NekDouble > &z, Array< OneD, Array< OneD, NekDouble > > &u, NekDouble time, const int o=0)
 Isentropic Vortex Test Case. More...
 
void GetExactIsentropicVortex (int field, Array< OneD, NekDouble > &outarray, NekDouble time)
 Compute the exact solution for the isentropic vortex problem. More...
 
void SetInitialIsentropicVortex (NekDouble initialtime)
 Set the initial condition for the isentropic vortex problem. More...
 
void SetBoundaryIsentropicVortex (int bcRegion, NekDouble time, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Set the boundary conditions for the isentropic vortex problem. More...
 
void GetExactRinglebFlow (int field, Array< OneD, NekDouble > &outarray)
 Ringleb Flow Test Case. More...
 
void SetInitialRinglebFlow (void)
 Set the initial condition for the Ringleb flow problem. More...
 
void SetBoundaryRinglebFlow (int bcRegion, NekDouble time, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
 Set the boundary conditions for the Ringleb flow problem. More...
 

Friends

class MemoryManager< EulerCFE >
 

Additional Inherited Members

- Protected Types inherited from Nektar::SolverUtils::EquationSystem
enum  HomogeneousType { eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, eNotHomogeneous }
 Parameter for homogeneous expansions. More...
 
- Protected Attributes inherited from Nektar::CompressibleFlowSystem
SolverUtils::RiemannSolverSharedPtr m_riemannSolver
 
SolverUtils::RiemannSolverSharedPtr m_riemannSolverLDG
 
SolverUtils::AdvectionSharedPtr m_advection
 
SolverUtils::DiffusionSharedPtr m_diffusion
 
Array< OneD, Array< OneD,
NekDouble > > 
m_vecLocs
 
NekDouble m_gamma
 
NekDouble m_pInf
 
NekDouble m_rhoInf
 
NekDouble m_uInf
 
NekDouble m_vInf
 
NekDouble m_wInf
 
NekDouble m_UInf
 
NekDouble m_gasConstant
 
NekDouble m_Twall
 
std::string m_ViscosityType
 
std::string m_shockCaptureType
 
std::string m_EqTypeStr
 
NekDouble m_mu
 
NekDouble m_Skappa
 
NekDouble m_Kappa
 
NekDouble m_mu0
 
NekDouble m_FacL
 
NekDouble m_FacH
 
NekDouble m_eps_max
 
NekDouble m_thermalConductivity
 
NekDouble m_Cp
 
NekDouble m_C1
 
NekDouble m_C2
 
NekDouble m_hFactor
 
NekDouble m_Prandtl
 
NekDouble m_amplitude
 
NekDouble m_omega
 
std::ofstream m_errFile
 
int m_steadyStateSteps
 
NekDouble m_steadyStateTol
 
std::vector
< SolverUtils::ForcingSharedPtr
m_forcing
 
StdRegions::StdQuadExpSharedPtr m_OrthoQuadExp
 
StdRegions::StdHexExpSharedPtr m_OrthoHexExp
 
bool m_smoothDiffusion
 
Array< OneD, NekDoublem_pressureStorage
 
Array< OneD, Array< OneD,
NekDouble > > 
m_fieldStorage
 
Array< OneD, Array< OneD,
NekDouble > > 
m_un
 
- Protected Attributes inherited from Nektar::SolverUtils::UnsteadySystem
int m_infosteps
 Number of time steps between outputting status information. More...
 
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
 Wrapper to the time integration scheme. More...
 
LibUtilities::TimeIntegrationSchemeOperators m_ode
 The time integration scheme operators to use. More...
 
LibUtilities::TimeIntegrationSolutionSharedPtr m_intSoln
 
NekDouble m_epsilon
 
bool m_explicitDiffusion
 Indicates if explicit or implicit treatment of diffusion is used. More...
 
bool m_explicitAdvection
 Indicates if explicit or implicit treatment of advection is used. More...
 
bool m_explicitReaction
 Indicates if explicit or implicit treatment of reaction is used. More...
 
bool m_homoInitialFwd
 Flag to determine if simulation should start in homogeneous forward transformed state. More...
 
std::vector< int > m_intVariables
 
std::vector< FilterSharedPtrm_filters
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 The session reader. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 Field input/output. More...
 
std::map< std::string, Array
< OneD, Array< OneD, float > > > 
m_interpWeights
 Map of the interpolation weights for a specific filename. More...
 
std::map< std::string, Array
< OneD, Array< OneD, unsigned
int > > > 
m_interpInds
 Map of the interpolation indices for a specific filename. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_fields
 Array holding all dependent variables. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_base
 Base fields. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_derivedfields
 Array holding all dependent variables. More...
 
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions
 Pointer to boundary conditions object. More...
 
SpatialDomains::MeshGraphSharedPtr m_graph
 Pointer to graph defining mesh. More...
 
std::string m_sessionName
 Name of the session. More...
 
NekDouble m_time
 Current time of simulation. More...
 
int m_initialStep
 Number of the step where the simulation should begin. More...
 
NekDouble m_fintime
 Finish time of the simulation. More...
 
NekDouble m_timestep
 Time step size. More...
 
NekDouble m_lambda
 Lambda constant in real system if one required. More...
 
std::set< std::string > m_loadedFields
 
NekDouble m_checktime
 Time between checkpoints. More...
 
int m_nchk
 Number of checkpoints written so far. More...
 
int m_steps
 Number of steps to take. More...
 
int m_checksteps
 Number of steps between checkpoints. More...
 
int m_spacedim
 Spatial dimension (>= expansion dim). More...
 
int m_expdim
 Expansion dimension. More...
 
bool m_singleMode
 Flag to determine if single homogeneous mode is used. More...
 
bool m_halfMode
 Flag to determine if half homogeneous mode is used. More...
 
bool m_multipleModes
 Flag to determine if use multiple homogenenous modes are used. More...
 
bool m_useFFT
 Flag to determine if FFT is used for homogeneous transform. More...
 
bool m_homogen_dealiasing
 Flag to determine if dealiasing is used for homogeneous simulations. More...
 
bool m_specHP_dealiasing
 Flag to determine if dealisising is usde for the Spectral/hp element discretisation. More...
 
enum MultiRegions::ProjectionType m_projectionType
 Type of projection; e.g continuous or discontinuous. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_traceNormals
 Array holding trace normals for DG simulations in the forwards direction. More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_gradtan
 1 x nvariable x nq More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_tanbasis
 2 x m_spacedim x nq More...
 
Array< OneD, bool > m_checkIfSystemSingular
 Flag to indicate if the fields should be checked for singularity. More...
 
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
 Map to identify relevant solver info to dump in output fields. More...
 
int m_NumQuadPointsError
 Number of Quadrature points used to work out the error. More...
 
enum HomogeneousType m_HomogeneousType
 
NekDouble m_LhomX
 physical length in X direction (if homogeneous) More...
 
NekDouble m_LhomY
 physical length in Y direction (if homogeneous) More...
 
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous) More...
 
int m_npointsX
 number of points in X direction (if homogeneous) More...
 
int m_npointsY
 number of points in Y direction (if homogeneous) More...
 
int m_npointsZ
 number of points in Z direction (if homogeneous) More...
 
int m_HomoDirec
 number of homogenous directions More...
 

Detailed Description

Definition at line 58 of file EulerCFE.h.

Constructor & Destructor Documentation

Nektar::EulerCFE::~EulerCFE ( )
virtual

problem type selector

Destructor for EulerCFE class.

Definition at line 97 of file EulerCFE.cpp.

98  {
99  }
Nektar::EulerCFE::EulerCFE ( const LibUtilities::SessionReaderSharedPtr pSession)
protected

Definition at line 54 of file EulerCFE.cpp.

56  : CompressibleFlowSystem(pSession)
57  {
58  }
CompressibleFlowSystem(const LibUtilities::SessionReaderSharedPtr &pSession)

Member Function Documentation

static SolverUtils::EquationSystemSharedPtr Nektar::EulerCFE::create ( const LibUtilities::SessionReaderSharedPtr pSession)
inlinestatic

Creates an instance of this class.

Definition at line 64 of file EulerCFE.h.

References Nektar::MemoryManager< DataType >::AllocateSharedPtr().

66  {
68  p->InitObject();
69  return p;
70  }
static boost::shared_ptr< DataType > AllocateSharedPtr()
Allocate a shared pointer from the memory pool.
boost::shared_ptr< EquationSystem > EquationSystemSharedPtr
A shared pointer to an EquationSystem object.
void Nektar::EulerCFE::DoOdeProjection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Compute the projection and call the method for imposing the boundary conditions in case of discontinuous projection.

Definition at line 175 of file EulerCFE.cpp.

References ASSERTL0, Nektar::MultiRegions::eDiscontinuous, Nektar::MultiRegions::eGalerkin, Nektar::MultiRegions::eMixed_CG_Discontinuous, Nektar::SolverUtils::EquationSystem::GetNpoints(), Nektar::SolverUtils::EquationSystem::m_projectionType, SetBoundaryConditions(), and Vmath::Vcopy().

Referenced by v_InitObject().

179  {
180  int i;
181  int nvariables = inarray.num_elements();
182 
183  switch (m_projectionType)
184  {
186  {
187  // Just copy over array
188  int npoints = GetNpoints();
189 
190  for(i = 0; i < nvariables; ++i)
191  {
192  Vmath::Vcopy(npoints, inarray[i], 1, outarray[i], 1);
193  }
194  SetBoundaryConditions(outarray, time);
195  break;
196  }
199  {
200  ASSERTL0(false, "No Continuous Galerkin for Euler equations");
201  break;
202  }
203  default:
204  ASSERTL0(false, "Unknown projection scheme");
205  break;
206  }
207  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
SOLVER_UTILS_EXPORT int GetNpoints()
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void SetBoundaryConditions(Array< OneD, Array< OneD, NekDouble > > &physarray, NekDouble time)
Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem ...
Definition: EulerCFE.cpp:220
void Nektar::EulerCFE::DoOdeRhs ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Compute the right-hand side.

Definition at line 144 of file EulerCFE.cpp.

References Nektar::SolverUtils::EquationSystem::GetNpoints(), Nektar::CompressibleFlowSystem::m_advection, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::CompressibleFlowSystem::m_forcing, Nektar::SolverUtils::EquationSystem::m_spacedim, and Vmath::Neg().

Referenced by v_InitObject().

148  {
149  int i;
150  int nvariables = inarray.num_elements();
151  int npoints = GetNpoints();
152 
153  Array<OneD, Array<OneD, NekDouble> > advVel(m_spacedim);
154 
155  m_advection->Advect(nvariables, m_fields, advVel, inarray,
156  outarray, time);
157 
158  for (i = 0; i < nvariables; ++i)
159  {
160  Vmath::Neg(npoints, outarray[i], 1);
161  }
162 
163  // Add sponge layer if defined in the session file
164  std::vector<SolverUtils::ForcingSharedPtr>::const_iterator x;
165  for (x = m_forcing.begin(); x != m_forcing.end(); ++x)
166  {
167  (*x)->Apply(m_fields, inarray, outarray, time);
168  }
169  }
int m_spacedim
Spatial dimension (>= expansion dim).
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:382
SolverUtils::AdvectionSharedPtr m_advection
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
std::vector< SolverUtils::ForcingSharedPtr > m_forcing
void Nektar::EulerCFE::EvaluateIsentropicVortex ( const Array< OneD, NekDouble > &  x,
const Array< OneD, NekDouble > &  y,
const Array< OneD, NekDouble > &  z,
Array< OneD, Array< OneD, NekDouble > > &  u,
NekDouble  time,
const int  o = 0 
)
private

Isentropic Vortex Test Case.

Definition at line 290 of file EulerCFE.cpp.

References Nektar::CompressibleFlowSystem::m_gamma, Nektar::SolverUtils::EquationSystem::m_spacedim, and Vmath::Zero().

Referenced by GetExactIsentropicVortex(), SetBoundaryIsentropicVortex(), and SetInitialIsentropicVortex().

297  {
298  int nq = x.num_elements();
299 
300  // Flow parameters
301  const NekDouble x0 = 5.0;
302  const NekDouble y0 = 0.0;
303  const NekDouble beta = 5.0;
304  const NekDouble u0 = 1.0;
305  const NekDouble v0 = 0.5;
306  const NekDouble gamma = m_gamma;
307  NekDouble r, xbar, ybar, tmp;
308  NekDouble fac = 1.0/(16.0*gamma*M_PI*M_PI);
309 
310  // In 3D zero rhow field.
311  if (m_spacedim == 3)
312  {
313  Vmath::Zero(nq, &u[3][o], 1);
314  }
315 
316  // Fill storage
317  for (int i = 0; i < nq; ++i)
318  {
319  xbar = x[i] - u0*time - x0;
320  ybar = y[i] - v0*time - y0;
321  r = sqrt(xbar*xbar + ybar*ybar);
322  tmp = beta*exp(1-r*r);
323  u[0][i+o] = pow(1.0 - (gamma-1.0)*tmp*tmp*fac, 1.0/(gamma-1.0));
324  u[1][i+o] = u[0][i+o]*(u0 - tmp*ybar/(2*M_PI));
325  u[2][i+o] = u[0][i+o]*(v0 + tmp*xbar/(2*M_PI));
326  u[m_spacedim+1][i+o] = pow(u[0][i+o], gamma)/(gamma-1.0) +
327  0.5*(u[1][i+o]*u[1][i+o] + u[2][i+o]*u[2][i+o]) / u[0][i+o];
328  }
329  }
int m_spacedim
Spatial dimension (>= expansion dim).
double NekDouble
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:359
void Nektar::EulerCFE::GetExactIsentropicVortex ( int  field,
Array< OneD, NekDouble > &  outarray,
NekDouble  time 
)
private

Compute the exact solution for the isentropic vortex problem.

Definition at line 334 of file EulerCFE.cpp.

References EvaluateIsentropicVortex(), Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_spacedim, and Vmath::Vcopy().

Referenced by v_EvaluateExactSolution().

338  {
339  int nTotQuadPoints = GetTotPoints();
340  Array<OneD, NekDouble> x(nTotQuadPoints);
341  Array<OneD, NekDouble> y(nTotQuadPoints);
342  Array<OneD, NekDouble> z(nTotQuadPoints);
343  Array<OneD, Array<OneD, NekDouble> > u(m_spacedim+2);
344 
345  m_fields[0]->GetCoords(x, y, z);
346 
347  for (int i = 0; i < m_spacedim + 2; ++i)
348  {
349  u[i] = Array<OneD, NekDouble>(nTotQuadPoints);
350  }
351 
352  EvaluateIsentropicVortex(x, y, z, u, time);
353 
354  Vmath::Vcopy(nTotQuadPoints, u[field], 1, outarray, 1);
355  }
void EvaluateIsentropicVortex(const Array< OneD, NekDouble > &x, const Array< OneD, NekDouble > &y, const Array< OneD, NekDouble > &z, Array< OneD, Array< OneD, NekDouble > > &u, NekDouble time, const int o=0)
Isentropic Vortex Test Case.
Definition: EulerCFE.cpp:290
SOLVER_UTILS_EXPORT int GetTotPoints()
int m_spacedim
Spatial dimension (>= expansion dim).
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Nektar::EulerCFE::GetExactRinglebFlow ( int  field,
Array< OneD, NekDouble > &  outarray 
)
private

Ringleb Flow Test Case.

Compute the exact solution for the Ringleb flow problem.

Definition at line 440 of file EulerCFE.cpp.

References ASSERTL0, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, and Nektar::CompressibleFlowSystem::m_gamma.

Referenced by v_EvaluateExactSolution().

443  {
444  int nTotQuadPoints = GetTotPoints();
445 
446  Array<OneD, NekDouble> rho(nTotQuadPoints, 100.0);
447  Array<OneD, NekDouble> rhou(nTotQuadPoints);
448  Array<OneD, NekDouble> rhov(nTotQuadPoints);
449  Array<OneD, NekDouble> E(nTotQuadPoints);
450  Array<OneD, NekDouble> x(nTotQuadPoints);
451  Array<OneD, NekDouble> y(nTotQuadPoints);
452  Array<OneD, NekDouble> z(nTotQuadPoints);
453 
454  m_fields[0]->GetCoords(x, y, z);
455 
456  // Flow parameters
457  NekDouble c, k, phi, r, J, VV, pp, sint, P, ss;
458  NekDouble J11, J12, J21, J22, det;
459  NekDouble Fx, Fy;
460  NekDouble xi, yi;
461  NekDouble dV;
462  NekDouble dtheta;
463  NekDouble par1;
464  NekDouble theta = M_PI / 4.0;
465  NekDouble kExt = 0.7;
466  NekDouble V = kExt * sin(theta);
467  NekDouble toll = 1.0e-8;
468  NekDouble errV = 1.0;
469  NekDouble errTheta = 1.0;
470  NekDouble gamma = m_gamma;
471  NekDouble gamma_1_2 = (gamma - 1.0) / 2.0;
472 
473  for (int i = 0; i < nTotQuadPoints; ++i)
474  {
475  while ((abs(errV) > toll) || (abs(errTheta) > toll))
476  {
477  VV = V * V;
478  sint = sin(theta);
479  c = sqrt(1.0 - gamma_1_2 * VV);
480  k = V / sint;
481  phi = 1.0 / k;
482  pp = phi * phi;
483  J = 1.0 / c + 1.0 / (3.0 * c * c * c) +
484  1.0 / (5.0 * c * c * c * c * c) -
485  0.5 * log((1.0 + c) / (1.0 - c));
486 
487  r = pow(c, 1.0 / gamma_1_2);
488  xi = 1.0 / (2.0 * r) * (1.0 / VV - 2.0 * pp) + J / 2.0;
489  yi = phi / (r * V) * sqrt(1.0 - VV * pp);
490  par1 = 25.0 - 5.0 * VV;
491  ss = sint * sint;
492 
493  Fx = xi - x[i];
494  Fy = yi - y[i];
495 
496  J11 = 39062.5 / pow(par1, 3.5) * (1.0 / VV - 2.0 / VV * ss) *
497  V + 1562.5 / pow(par1, 2.5) * (-2.0 / (VV * V) + 4.0 /
498  (VV * V) * ss) + 12.5 / pow(par1, 1.5) * V + 312.5 /
499  pow(par1, 2.5) * V + 7812.5 / pow(par1, 3.5) * V -
500  0.25 * (-1.0 / pow(par1, 0.5) * V/(1.0 - 0.2 *
501  pow(par1, 0.5)) - (1.0 + 0.2 * pow(par1, 0.5)) /
502  pow((1.0 - 0.2 * pow(par1, 0.5)), 2.0) /
503  pow(par1, 0.5) * V) / (1.0 + 0.2 * pow(par1, 0.5)) *
504  (1.0 - 0.2 * pow(par1, 0.5));
505 
506  J12 = -6250.0 / pow(par1, 2.5) / VV * sint * cos(theta);
507  J21 = -6250.0 / (VV * V) * sint /
508  pow(par1, 2.5) * pow((1.0 - ss), 0.5) +
509  78125.0 / V * sint / pow(par1, 3.5) *
510  pow((1.0 - ss), 0.5);
511 
512  // the matrix is singular when theta = pi/2
513  if(abs(y[i])<toll && abs(cos(theta))<toll)
514  {
515  J22 = -39062.5 / pow(par1, 3.5) / V + 3125 /
516  pow(par1, 2.5) / (VV * V) + 12.5 / pow(par1, 1.5) *
517  V + 312.5 / pow(par1, 2.5) * V + 7812.5 /
518  pow(par1, 3.5) * V - 0.25 * (-1.0 / pow(par1, 0.5) *
519  V / (1.0 - 0.2 * pow(par1, 0.5)) - (1.0 + 0.2 *
520  pow(par1, 0.5)) / pow((1.0 - 0.2 *
521  pow(par1, 0.5)), 2.0) / pow(par1, 0.5) * V) /
522  (1.0 + 0.2 * pow(par1, 0.5)) * (1.0 - 0.2 *
523  pow(par1,0.5));
524 
525  // dV = -dV/dx * Fx
526  dV = -1.0 / J22 * Fx;
527  dtheta = 0.0;
528  theta = M_PI / 2.0;
529  }
530  else
531  {
532  J22 = 3125.0 / VV * cos(theta) / pow(par1, 2.5) *
533  pow((1.0 - ss), 0.5) - 3125.0 / VV * ss /
534  pow(par1, 2.5) / pow((1.0 - ss), 0.5) * cos(theta);
535 
536  det = -1.0 / (J11 * J22 - J12 * J21);
537 
538  // [dV dtheta]' = -[invJ]*[Fx Fy]'
539  dV = det * ( J22 * Fx - J12 * Fy);
540  dtheta = det * (-J21 * Fx + J11 * Fy);
541  }
542 
543  V = V + dV;
544  theta = theta + dtheta;
545 
546  errV = abs(dV);
547  errTheta = abs(dtheta);
548 
549  }
550 
551  c = sqrt(1.0 - gamma_1_2 * VV);
552  r = pow(c, 1.0 / gamma_1_2);
553 
554  rho[i] = r;
555  rhou[i] = rho[i] * V * cos(theta);
556  rhov[i] = rho[i] * V * sin(theta);
557  P = (c * c) * rho[i] / gamma;
558  E[i] = P / (gamma - 1.0) + 0.5 *
559  (rhou[i] * rhou[i] / rho[i] + rhov[i] * rhov[i] / rho[i]);
560 
561  // Resetting the guess value
562  errV = 1.0;
563  errTheta = 1.0;
564  theta = M_PI/4.0;
565  V = kExt*sin(theta);
566  }
567 
568  switch (field)
569  {
570  case 0:
571  outarray = rho;
572  break;
573  case 1:
574  outarray = rhou;
575  break;
576  case 2:
577  outarray = rhov;
578  break;
579  case 3:
580  outarray = E;
581  break;
582  default:
583  ASSERTL0(false, "Error in variable number!");
584  break;
585  }
586  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
SOLVER_UTILS_EXPORT int GetTotPoints()
double NekDouble
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Nektar::EulerCFE::SetBoundaryConditions ( Array< OneD, Array< OneD, NekDouble > > &  inarray,
NekDouble  time 
)
private

Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem level since they are compressible solver specific; b) Isentropic vortex and Ringleb flow BCs implemented at EulerCFE level since they are Euler solver specific; c) Time dependent BCs.

Parameters
inarrayfields array.
timetime.

Definition at line 220 of file EulerCFE.cpp.

References ASSERTL0, Nektar::SolverUtils::EquationSystem::m_fields, SetBoundaryIsentropicVortex(), SetBoundaryRinglebFlow(), and Nektar::CompressibleFlowSystem::SetCommonBC().

Referenced by DoOdeProjection().

223  {
224  std::string varName;
225  int cnt = 0;
226 
227  std::string userDefStr;
228  int nreg = m_fields[0]->GetBndConditions().num_elements();
229  // Loop over Boundary Regions
230  for (int n = 0; n < nreg; ++n)
231  {
232  userDefStr = m_fields[0]->GetBndConditions()[n]->GetUserDefined();
233  if(!userDefStr.empty())
234  {
235  if(boost::iequals(userDefStr,"WallViscous"))
236  {
237  ASSERTL0(false, "WallViscous is a wrong bc for the "
238  "Euler equations");
239  }
240  else if(boost::iequals(userDefStr, "IsentropicVortex"))
241  {
242  // Isentropic Vortex Boundary Condition
243  SetBoundaryIsentropicVortex(n, time, cnt, inarray);
244  }
245  else if (boost::iequals(userDefStr,"RinglebFlow"))
246  {
247  SetBoundaryRinglebFlow(n, time, cnt, inarray);
248  }
249  else
250  {
251  // set up userdefined BC common to all solvers
252  SetCommonBC(userDefStr,n,time, cnt,inarray);
253  }
254  }
255 
256  // no User Defined conditions provided so skip cnt
257  cnt += m_fields[0]->GetBndCondExpansions()[n]->GetExpSize();
258  }
259  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
void SetCommonBC(const std::string &userDefStr, const int n, const NekDouble time, int &cnt, Array< OneD, Array< OneD, NekDouble > > &inarray)
Set boundary conditions which can be: a) Wall and Symmerty BCs implemented at CompressibleFlowSystem ...
void SetBoundaryRinglebFlow(int bcRegion, NekDouble time, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
Set the boundary conditions for the Ringleb flow problem.
Definition: EulerCFE.cpp:835
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void SetBoundaryIsentropicVortex(int bcRegion, NekDouble time, int cnt, Array< OneD, Array< OneD, NekDouble > > &physarray)
Set the boundary conditions for the isentropic vortex problem.
Definition: EulerCFE.cpp:390
void Nektar::EulerCFE::SetBoundaryIsentropicVortex ( int  bcRegion,
NekDouble  time,
int  cnt,
Array< OneD, Array< OneD, NekDouble > > &  physarray 
)
private

Set the boundary conditions for the isentropic vortex problem.

Definition at line 390 of file EulerCFE.cpp.

References EvaluateIsentropicVortex(), Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, and Vmath::Vcopy().

Referenced by SetBoundaryConditions().

395  {
396  int nvariables = physarray.num_elements();
397  int nTraceNumPoints = GetTraceTotPoints();
398  Array<OneD, Array<OneD, NekDouble> > Fwd(nvariables);
399  const Array<OneD, const int> &bndTraceMap = m_fields[0]->GetTraceBndMap();
400 
401  // Get physical values of the forward trace (from exp to phys)
402  for (int i = 0; i < nvariables; ++i)
403  {
404  Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
405  m_fields[i]->ExtractTracePhys(physarray[i], Fwd[i]);
406  }
407 
408  int id2, e_max;
409  e_max = m_fields[0]->GetBndCondExpansions()[bcRegion]->GetExpSize();
410 
411  for(int e = 0; e < e_max; ++e)
412  {
413  int npoints = m_fields[0]->
414  GetBndCondExpansions()[bcRegion]->GetExp(e)->GetTotPoints();
415  int id1 = m_fields[0]->
416  GetBndCondExpansions()[bcRegion]->GetPhys_Offset(e);
417  id2 = m_fields[0]->GetTrace()->GetPhys_Offset(bndTraceMap[cnt++]);
418 
419  Array<OneD,NekDouble> x(npoints, 0.0);
420  Array<OneD,NekDouble> y(npoints, 0.0);
421  Array<OneD,NekDouble> z(npoints, 0.0);
422 
423  m_fields[0]->GetBndCondExpansions()[bcRegion]->
424  GetExp(e)->GetCoords(x, y, z);
425 
426  EvaluateIsentropicVortex(x, y, z, Fwd, time, id2);
427 
428  for (int i = 0; i < nvariables; ++i)
429  {
430  Vmath::Vcopy(npoints, &Fwd[i][id2], 1,
431  &(m_fields[i]->GetBndCondExpansions()[bcRegion]->
432  UpdatePhys())[id1], 1);
433  }
434  }
435  }
void EvaluateIsentropicVortex(const Array< OneD, NekDouble > &x, const Array< OneD, NekDouble > &y, const Array< OneD, NekDouble > &z, Array< OneD, Array< OneD, NekDouble > > &u, NekDouble time, const int o=0)
Isentropic Vortex Test Case.
Definition: EulerCFE.cpp:290
SOLVER_UTILS_EXPORT int GetTraceTotPoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Nektar::EulerCFE::SetBoundaryRinglebFlow ( int  bcRegion,
NekDouble  time,
int  cnt,
Array< OneD, Array< OneD, NekDouble > > &  physarray 
)
private

Set the boundary conditions for the Ringleb flow problem.

Definition at line 835 of file EulerCFE.cpp.

References Nektar::LibUtilities::eFunctionTypeFile, Nektar::SolverUtils::EquationSystem::eHomogeneous1D, Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_expdim, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::CompressibleFlowSystem::m_gamma, Nektar::SolverUtils::EquationSystem::m_HomogeneousType, Nektar::SolverUtils::EquationSystem::m_session, and Vmath::Vcopy().

Referenced by SetBoundaryConditions().

840  {
841  int nvariables = physarray.num_elements();
842  int nTraceNumPoints = GetTraceTotPoints();
843  Array<OneD, Array<OneD, NekDouble> > Fwd(nvariables);
844 
845  // For 3DHomogenoeus1D
846  int n_planes = 1;
848  {
849  int nPointsTot = m_fields[0]->GetTotPoints();
850  int nPointsTot_plane = m_fields[0]->GetPlane(0)->GetTotPoints();
851  n_planes = nPointsTot/nPointsTot_plane;
852  nTraceNumPoints = nTraceNumPoints * n_planes;
853  }
854 
855  // Get physical values of the forward trace (from exp to phys)
856  for (int i = 0; i < nvariables; ++i)
857  {
858  Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
859  m_fields[i]->ExtractTracePhys(physarray[i], Fwd[i]);
860  }
861 
862  int id2, id2_plane, e_max;
863 
864  e_max = m_fields[0]->GetBndCondExpansions()[bcRegion]->GetExpSize();
865 
866  for(int e = 0; e < e_max; ++e)
867  {
868  int npoints = m_fields[0]->
869  GetBndCondExpansions()[bcRegion]->GetExp(e)->GetTotPoints();
870  int id1 = m_fields[0]->
871  GetBndCondExpansions()[bcRegion]->GetPhys_Offset(e);
872 
873  // For 3DHomogenoeus1D
875  {
876  int cnt_plane = cnt/n_planes;
877  int e_plane;
878  int e_max_plane = e_max/n_planes;
879  int nTracePts_plane = GetTraceTotPoints();
880 
881  int planeID = floor((e + 0.5 )/ e_max_plane );
882  e_plane = e - e_max_plane*planeID;
883 
884  id2_plane = m_fields[0]->GetTrace()->GetPhys_Offset(
885  m_fields[0]->GetTraceMap()->
886  GetBndCondCoeffsToGlobalCoeffsMap(
887  cnt_plane + e_plane));
888  id2 = id2_plane + planeID*nTracePts_plane;
889  }
890  else // For general case
891  {
892  id2 = m_fields[0]->
893  GetTrace()->GetPhys_Offset(m_fields[0]->GetTraceMap()->
894  GetBndCondTraceToGlobalTraceMap(cnt++));
895  }
896 
897  Array<OneD,NekDouble> x0(npoints, 0.0);
898  Array<OneD,NekDouble> x1(npoints, 0.0);
899  Array<OneD,NekDouble> x2(npoints, 0.0);
900 
901  m_fields[0]->GetBndCondExpansions()[bcRegion]->
902  GetExp(e)->GetCoords(x0, x1, x2);
903 
904  // Flow parameters
905  NekDouble c, k, phi, r, J, VV, pp, sint, P, ss;
906  NekDouble J11, J12, J21, J22, det;
907  NekDouble Fx, Fy;
908  NekDouble xi, yi;
909  NekDouble dV;
910  NekDouble dtheta;
911  NekDouble par1;
912  NekDouble theta = M_PI / 4.0;
913  NekDouble kExt = 0.7;
914  NekDouble V = kExt * sin(theta);
915  NekDouble toll = 1.0e-8;
916  NekDouble errV = 1.0;
917  NekDouble errTheta = 1.0;
918  NekDouble gamma = m_gamma;
919  NekDouble gamma_1_2 = (gamma - 1.0) / 2.0;
920 
921  // Loop on all the points of that edge
922  for (int j = 0; j < npoints; j++)
923  {
924 
925  while ((abs(errV) > toll) || (abs(errTheta) > toll))
926  {
927  VV = V * V;
928  sint = sin(theta);
929  c = sqrt(1.0 - gamma_1_2 * VV);
930  k = V / sint;
931  phi = 1.0 / k;
932  pp = phi * phi;
933  J = 1.0 / c + 1.0 / (3.0 * c * c * c) +
934  1.0 / (5.0 * c * c * c * c * c) -
935  0.5 * log((1.0 + c) / (1.0 - c));
936 
937  r = pow(c, 1.0 / gamma_1_2);
938  xi = 1.0 / (2.0 * r) * (1.0 / VV - 2.0 * pp) + J / 2.0;
939  yi = phi / (r * V) * sqrt(1.0 - VV * pp);
940  par1 = 25.0 - 5.0 * VV;
941  ss = sint * sint;
942 
943  Fx = xi - x0[j];
944  Fy = yi - x1[j];
945 
946  J11 = 39062.5 / pow(par1, 3.5) *
947  (1.0 / VV - 2.0 / VV * ss) * V + 1562.5 /
948  pow(par1, 2.5) * (-2.0 / (VV * V) + 4.0 /
949  (VV * V) * ss) + 12.5 / pow(par1, 1.5) * V +
950  312.5 / pow(par1, 2.5) * V + 7812.5 /
951  pow(par1, 3.5) * V - 0.25 *
952  (-1.0 / pow(par1, 0.5) * V / (1.0 - 0.2 *
953  pow(par1, 0.5)) - (1.0 + 0.2 * pow(par1, 0.5)) /
954  pow((1.0 - 0.2 * pow(par1, 0.5)), 2.0) /
955  pow(par1, 0.5) * V) / (1.0 + 0.2 * pow(par1, 0.5)) *
956  (1.0 - 0.2 * pow(par1, 0.5));
957 
958  J12 = -6250.0 / pow(par1, 2.5) / VV * sint * cos(theta);
959  J21 = -6250.0 / (VV * V) * sint / pow(par1, 2.5) *
960  pow((1.0 - ss), 0.5) + 78125.0 / V * sint /
961  pow(par1, 3.5) * pow((1.0 - ss), 0.5);
962 
963  // the matrix is singular when theta = pi/2
964  if (abs(x1[j]) < toll && abs(cos(theta)) < toll)
965  {
966  J22 = -39062.5 / pow(par1, 3.5) / V + 3125 /
967  pow(par1, 2.5) / (VV * V) + 12.5 /
968  pow(par1, 1.5) * V + 312.5 / pow(par1, 2.5) *
969  V + 7812.5 / pow(par1, 3.5) * V - 0.25 *
970  (-1.0 / pow(par1, 0.5) * V / (1.0 - 0.2 *
971  pow(par1, 0.5)) - (1.0 + 0.2 * pow(par1, 0.5)) /
972  pow((1.0 - 0.2 * pow(par1, 0.5)), 2.0) /
973  pow(par1, 0.5) * V) / (1.0 + 0.2 *
974  pow(par1, 0.5)) * (1.0 - 0.2 * pow(par1, 0.5));
975 
976  // dV = -dV/dx * Fx
977  dV = -1.0 / J22 * Fx;
978  dtheta = 0.0;
979  theta = M_PI / 2.0;
980  }
981  else
982  {
983  J22 = 3125.0 / VV * cos(theta) / pow(par1, 2.5) *
984  pow((1.0 - ss), 0.5) - 3125.0 / VV * ss /
985  pow(par1, 2.5) / pow((1.0 - ss), 0.5) *
986  cos(theta);
987 
988  det = -1.0 / (J11 * J22 - J12 * J21);
989 
990  // [dV dtheta]' = -[invJ]*[Fx Fy]'
991  dV = det * ( J22 * Fx - J12 * Fy);
992  dtheta = det * (-J21 * Fx + J11 * Fy);
993  }
994 
995  V = V + dV;
996  theta = theta + dtheta;
997 
998  errV = abs(dV);
999  errTheta = abs(dtheta);
1000  }
1001 
1002  c = sqrt(1.0 - gamma_1_2 * VV);
1003  int kk = id2 + j;
1004  NekDouble timeramp = 200.0;
1005  std::string restartstr = "RESTART";
1006  if (time<timeramp &&
1007  !(m_session->DefinesFunction("InitialConditions") &&
1008  m_session->GetFunctionType("InitialConditions", 0) ==
1010  {
1011  Fwd[0][kk] = pow(c, 1.0 / gamma_1_2) *
1012  exp(-1.0 + time /timeramp);
1013 
1014  Fwd[1][kk] = Fwd[0][kk] * V * cos(theta) *
1015  exp(-1 + time / timeramp);
1016 
1017  Fwd[2][kk] = Fwd[0][kk] * V * sin(theta) *
1018  exp(-1 + time / timeramp);
1019  }
1020  else
1021  {
1022  Fwd[0][kk] = pow(c, 1.0 / gamma_1_2);
1023  Fwd[1][kk] = Fwd[0][kk] * V * cos(theta);
1024  Fwd[2][kk] = Fwd[0][kk] * V * sin(theta);
1025  }
1026 
1027  P = (c * c) * Fwd[0][kk] / gamma;
1028  Fwd[3][kk] = P / (gamma - 1.0) + 0.5 *
1029  (Fwd[1][kk] * Fwd[1][kk] / Fwd[0][kk] +
1030  Fwd[2][kk] * Fwd[2][kk] / Fwd[0][kk]);
1031 
1032  errV = 1.0;
1033  errTheta = 1.0;
1034  theta = M_PI / 4.0;
1035  V = kExt * sin(theta);
1036  }
1037 
1038  for (int i = 0; i < nvariables; ++i)
1039  {
1040  Vmath::Vcopy(npoints, &Fwd[i][id2], 1,
1041  &(m_fields[i]->GetBndCondExpansions()[bcRegion]->
1042  UpdatePhys())[id1],1);
1043  }
1044  }
1045  }
int m_expdim
Expansion dimension.
double NekDouble
SOLVER_UTILS_EXPORT int GetTraceTotPoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
enum HomogeneousType m_HomogeneousType
void Nektar::EulerCFE::SetInitialIsentropicVortex ( NekDouble  initialtime)
private

Set the initial condition for the isentropic vortex problem.

Definition at line 360 of file EulerCFE.cpp.

References EvaluateIsentropicVortex(), Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_spacedim, and Vmath::Vcopy().

Referenced by v_SetInitialConditions().

361  {
362  int nTotQuadPoints = GetTotPoints();
363  Array<OneD, NekDouble> x(nTotQuadPoints);
364  Array<OneD, NekDouble> y(nTotQuadPoints);
365  Array<OneD, NekDouble> z(nTotQuadPoints);
366  Array<OneD, Array<OneD, NekDouble> > u(m_spacedim+2);
367 
368  m_fields[0]->GetCoords(x, y, z);
369 
370  for (int i = 0; i < m_spacedim + 2; ++i)
371  {
372  u[i] = Array<OneD, NekDouble>(nTotQuadPoints);
373  }
374 
375  EvaluateIsentropicVortex(x, y, z, u, initialtime);
376 
377  // Forward transform to fill the coefficient space
378  for(int i = 0; i < m_fields.num_elements(); ++i)
379  {
380  Vmath::Vcopy(nTotQuadPoints, u[i], 1, m_fields[i]->UpdatePhys(), 1);
381  m_fields[i]->SetPhysState(true);
382  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
383  m_fields[i]->UpdateCoeffs());
384  }
385  }
void EvaluateIsentropicVortex(const Array< OneD, NekDouble > &x, const Array< OneD, NekDouble > &y, const Array< OneD, NekDouble > &z, Array< OneD, Array< OneD, NekDouble > > &u, NekDouble time, const int o=0)
Isentropic Vortex Test Case.
Definition: EulerCFE.cpp:290
SOLVER_UTILS_EXPORT int GetTotPoints()
int m_spacedim
Spatial dimension (>= expansion dim).
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1047
void Nektar::EulerCFE::SetInitialRinglebFlow ( void  )
private

Set the initial condition for the Ringleb flow problem.

Definition at line 591 of file EulerCFE.cpp.

References Nektar::SolverUtils::EquationSystem::m_fields, Nektar::CompressibleFlowSystem::m_gamma, Nektar::SolverUtils::EquationSystem::m_sessionName, and Nektar::SolverUtils::EquationSystem::WriteFld().

592  {
593  // Get number of different boundaries in the input file
594  int nbnd = m_fields[0]->GetBndConditions().num_elements();
595 
596  // Loop on all the edges of the input file
597  for(int bcRegion=0; bcRegion < nbnd; ++bcRegion)
598  {
599 
600  int npoints = m_fields[0]->
601  GetBndCondExpansions()[bcRegion]->GetNpoints();
602 
603  Array<OneD,NekDouble> x0(npoints, 0.0);
604  Array<OneD,NekDouble> x1(npoints, 0.0);
605  Array<OneD,NekDouble> x2(npoints, 0.0);
606 
607  Array<OneD, NekDouble> rho(npoints, 0.0);
608  Array<OneD, NekDouble> rhou(npoints, 0.0);
609  Array<OneD, NekDouble> rhov(npoints, 0.0);
610  Array<OneD, NekDouble> E(npoints, 0.0);
611 
612  m_fields[0]->GetBndCondExpansions()[bcRegion]->
613  GetCoords(x0, x1, x2);
614 
615  // Flow parameters
616  NekDouble c, k, phi, r, J, VV, pp, sint, P, ss;
617  NekDouble J11, J12, J21, J22, det;
618  NekDouble Fx, Fy;
619  NekDouble xi, yi;
620  NekDouble dV;
621  NekDouble dtheta;
622  NekDouble par1;
623  NekDouble theta = M_PI / 4.0;
624  NekDouble kExt = 0.7;
625  NekDouble V = kExt * sin(theta);
626  NekDouble toll = 1.0e-8;
627  NekDouble errV = 1.0;
628  NekDouble errTheta = 1.0;
629  NekDouble gamma = m_gamma;
630  NekDouble gamma_1_2 = (gamma - 1.0) / 2.0;
631 
632  // Loop on all the points of that edge
633  for (int j = 0; j < npoints; j++)
634  {
635  while ((abs(errV) > toll) || (abs(errTheta) > toll))
636  {
637 
638  VV = V * V;
639  sint = sin(theta);
640  c = sqrt(1.0 - gamma_1_2 * VV);
641  k = V / sint;
642  phi = 1.0 / k;
643  pp = phi * phi;
644  J = 1.0 / c + 1.0 / (3.0 * c * c * c) +
645  1.0 / (5.0 * c * c * c * c * c) -
646  0.5 * log((1.0 + c) / (1.0 - c));
647 
648  r = pow(c, 1.0 / gamma_1_2);
649  xi = 1.0 / (2.0 * r) * (1.0 / VV - 2.0 * pp) + J / 2.0;
650  yi = phi / (r * V) * sqrt(1.0 - VV * pp);
651  par1 = 25.0 - 5.0 * VV;
652  ss = sint * sint;
653 
654  Fx = xi - x0[j];
655  Fy = yi - x1[j];
656 
657  J11 = 39062.5 / pow(par1, 3.5) *
658  (1.0 / VV - 2.0 / VV * ss) * V +
659  1562.5 / pow(par1, 2.5) * (-2.0 /
660  (VV * V) + 4.0 / (VV * V) * ss) +
661  12.5 / pow(par1, 1.5) * V +
662  312.5 / pow(par1, 2.5) * V +
663  7812.5 / pow(par1, 3.5) * V -
664  0.25 * (-1.0 / pow(par1, 0.5) * V /
665  (1.0 - 0.2 * pow(par1, 0.5)) - (1.0 + 0.2 *
666  pow(par1, 0.5)) / pow((1.0 - 0.2 *
667  pow(par1, 0.5)), 2.0) /
668  pow(par1, 0.5) * V) /
669  (1.0 + 0.2 * pow(par1, 0.5)) *
670  (1.0 - 0.2 * pow(par1, 0.5));
671 
672  J12 = -6250.0 / pow(par1, 2.5) / VV * sint * cos(theta);
673  J21 = -6250.0 / (VV * V) * sint / pow(par1, 2.5) *
674  pow((1.0 - ss), 0.5) + 78125.0 / V * sint /
675  pow(par1, 3.5) * pow((1.0 - ss), 0.5);
676 
677  // the matrix is singular when theta = pi/2
678  if (abs(x1[j]) < toll && abs(cos(theta)) < toll)
679  {
680 
681  J22 = -39062.5 / pow(par1, 3.5) / V +
682  3125 / pow(par1, 2.5) / (VV * V) + 12.5 /
683  pow(par1, 1.5) * V + 312.5 / pow(par1, 2.5) * V +
684  7812.5 / pow(par1, 3.5) * V -
685  0.25 * (-1.0 / pow(par1, 0.5) * V /
686  (1.0 - 0.2 * pow(par1, 0.5)) -
687  (1.0 + 0.2 * pow(par1, 0.5)) /
688  pow((1.0 - 0.2* pow(par1, 0.5)), 2.0) /
689  pow(par1, 0.5) *V) /
690  (1.0 + 0.2 * pow(par1, 0.5)) *
691  (1.0 - 0.2 * pow(par1, 0.5));
692 
693  // dV = -dV/dx * Fx
694  dV = -1.0 / J22 * Fx;
695  dtheta = 0.0;
696  theta = M_PI / 2.0;
697  }
698  else
699  {
700 
701  J22 = 3125.0 / VV * cos(theta) / pow(par1, 2.5) *
702  pow((1.0 - ss), 0.5) - 3125.0 / VV * ss /
703  pow(par1, 2.5) / pow((1.0 - ss), 0.5) *
704  cos(theta);
705 
706  det = -1.0 / (J11 * J22 - J12 * J21);
707 
708  // [dV dtheta]' = -[invJ]*[Fx Fy]'
709  dV = det * ( J22 * Fx - J12 * Fy);
710  dtheta = det * (-J21 * Fx + J11 * Fy);
711  }
712 
713  V = V + dV;
714  theta = theta + dtheta;
715 
716  errV = abs(dV);
717  errTheta = abs(dtheta);
718 
719  }
720 
721  c = sqrt(1.0 - gamma_1_2 * VV);
722  rho[j] = pow(c, 1.0 / gamma_1_2) * exp(-1.0);
723  rhou[j] = rho[j] * V * cos(theta) * exp(-1.0);
724  rhov[j] = rho[j] * V * sin(theta) * exp(-1.0);
725  P = (c * c) * rho[j] / gamma;
726  E[j] = P / (gamma - 1.0) + 0.5 *
727  (rhou[j] * rhou[j] /
728  rho[j] + rhov[j] * rhov[j] / rho[j]);
729  errV = 1.0;
730  errTheta = 1.0;
731  theta = M_PI / 4.0;
732  V = kExt * sin(theta);
733 
734  }
735 
736  // Fill the physical space
737  m_fields[0]->GetBndCondExpansions()[bcRegion]->SetPhys(rho);
738  m_fields[1]->GetBndCondExpansions()[bcRegion]->SetPhys(rhou);
739  m_fields[2]->GetBndCondExpansions()[bcRegion]->SetPhys(rhov);
740  m_fields[3]->GetBndCondExpansions()[bcRegion]->SetPhys(E);
741 
742  // Forward transform to fill the coefficients space
743  for(int i = 0; i < m_fields.num_elements(); ++i)
744  {
745  m_fields[i]->GetBndCondExpansions()[bcRegion]->
746  FwdTrans_BndConstrained(
747  m_fields[i]->GetBndCondExpansions()[bcRegion]->
748  GetPhys(),
749  m_fields[i]->GetBndCondExpansions()[bcRegion]->
750  UpdateCoeffs());
751  }
752 
753  }
754 
755  // Calculation of the initial internal values as a weighted average
756  // over the distance from the boundaries
757  int nq = m_fields[0]->GetNpoints();
758 
759  Array<OneD,NekDouble> x0(nq);
760  Array<OneD,NekDouble> x1(nq);
761  Array<OneD,NekDouble> x2(nq);
762 
763  // Get the coordinates
764  // (assuming all fields have the same discretisation)
765  m_fields[0]->GetCoords(x0, x1, x2);
766 
767  for(int j = 0; j < nq; j++)
768  {
769  NekDouble Dist = 0.0;
770  NekDouble rho = 0.0;
771  NekDouble rhou = 0.0;
772  NekDouble rhov = 0.0;
773  NekDouble E = 0.0;
774  NekDouble SumDist = 0.0;
775 
776  // Calculation of all the distances
777  // Loop on all the edges of the input file
778  for (int bcRegion = 0; bcRegion < nbnd; ++bcRegion)
779  {
780  // Get quadrature points on the edge
781  int npoints = m_fields[0]->
782  GetBndCondExpansions()[bcRegion]->GetNpoints();
783 
784  Array<OneD,NekDouble> xb0(npoints, 0.0);
785  Array<OneD,NekDouble> xb1(npoints, 0.0);
786  Array<OneD,NekDouble> xb2(npoints, 0.0);
787 
788  m_fields[0]->GetBndCondExpansions()[bcRegion]->
789  GetCoords(xb0, xb1, xb2);
790 
791  for (int k = 0; k < npoints; k++)
792  {
793  Dist = sqrt((xb0[k] - x0[j]) * (xb0[k] - x0[j]) +
794  (xb1[k] - x1[j]) * (xb1[k] - x1[j]));
795 
796  SumDist += Dist;
797  rho += Dist * (m_fields[0]->
798  GetBndCondExpansions()[bcRegion]->GetPhys())[k];
799  rhou += Dist * (m_fields[1]->
800  GetBndCondExpansions()[bcRegion]->GetPhys())[k];
801  rhov += Dist * (m_fields[2]->
802  GetBndCondExpansions()[bcRegion]->GetPhys())[k];
803  E += Dist * (m_fields[3]->
804  GetBndCondExpansions()[bcRegion]->GetPhys())[k];
805  }
806  }
807 
808  rho = rho / SumDist;
809  rhou = rhou / SumDist;
810  rhov = rhov / SumDist;
811  E = E / SumDist;
812 
813  (m_fields[0]->UpdatePhys())[j] = rho;
814  (m_fields[1]->UpdatePhys())[j] = rhou;
815  (m_fields[2]->UpdatePhys())[j] = rhov;
816  (m_fields[3]->UpdatePhys())[j] = E;
817 
818  }
819 
820  for (int i = 0 ; i < m_fields.num_elements(); i++)
821  {
822  m_fields[i]->SetPhysState(true);
823  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
824  m_fields[i]->UpdateCoeffs());
825  }
826 
827  // Dump initial conditions to file
828  std::string outname = m_sessionName + "_initialRingleb.chk";
829  WriteFld(outname);
830  }
std::string m_sessionName
Name of the session.
double NekDouble
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Nektar::EulerCFE::v_EvaluateExactSolution ( unsigned int  field,
Array< OneD, NekDouble > &  outfield,
const NekDouble  time = 0.0 
)
protectedvirtual

Get the exact solutions for isentropic vortex and Ringleb flow problems.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 265 of file EulerCFE.cpp.

References Nektar::eIsentropicVortex, Nektar::eRinglebFlow, GetExactIsentropicVortex(), GetExactRinglebFlow(), m_problemType, and Nektar::SolverUtils::EquationSystem::v_EvaluateExactSolution().

269  {
270  switch(m_problemType)
271  {
272  case eIsentropicVortex:
273  {
274  GetExactIsentropicVortex(field, outfield, time);
275  break;
276  }
277  case eRinglebFlow:
278  {
279  GetExactRinglebFlow(field, outfield);
280  break;
281  }
282  default:
283  {
284  EquationSystem::v_EvaluateExactSolution(field, outfield, time);
285  break;
286  }
287  }
288  }
Ringleb Flow.
Definition: EulerCFE.h:47
void GetExactIsentropicVortex(int field, Array< OneD, NekDouble > &outarray, NekDouble time)
Compute the exact solution for the isentropic vortex problem.
Definition: EulerCFE.cpp:334
void GetExactRinglebFlow(int field, Array< OneD, NekDouble > &outarray)
Ringleb Flow Test Case.
Definition: EulerCFE.cpp:440
Isentropic Vortex.
Definition: EulerCFE.h:46
ProblemType m_problemType
Definition: EulerCFE.h:77
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
void Nektar::EulerCFE::v_GenerateSummary ( SolverUtils::SummaryList s)
protectedvirtual

Print a summary of time stepping parameters.

Print out a summary with some relevant information.

Reimplemented from Nektar::CompressibleFlowSystem.

Definition at line 104 of file EulerCFE.cpp.

References Nektar::SolverUtils::AddSummaryItem(), m_problemType, Nektar::ProblemTypeMap, and Nektar::CompressibleFlowSystem::v_GenerateSummary().

105  {
108  }
virtual void v_GenerateSummary(SolverUtils::SummaryList &s)
Print a summary of time stepping parameters.
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:50
ProblemType m_problemType
Definition: EulerCFE.h:77
const char *const ProblemTypeMap[]
Definition: EulerADCFE.h:50
void Nektar::EulerCFE::v_InitObject ( )
protectedvirtual

Initialization object for CompressibleFlowSystem class.

Reimplemented from Nektar::CompressibleFlowSystem.

Definition at line 60 of file EulerCFE.cpp.

References ASSERTL0, Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineProjection(), DoOdeProjection(), DoOdeRhs(), Nektar::SolverUtils::UnsteadySystem::m_explicitAdvection, Nektar::SolverUtils::UnsteadySystem::m_ode, m_problemType, Nektar::SolverUtils::EquationSystem::m_session, Nektar::ProblemTypeMap, Nektar::SIZE_ProblemType, and Nektar::CompressibleFlowSystem::v_InitObject().

61  {
63 
64  if(m_session->DefinesSolverInfo("PROBLEMTYPE"))
65  {
66  int i;
67  std::string ProblemTypeStr =
68  m_session->GetSolverInfo("PROBLEMTYPE");
69  for (i = 0; i < (int) SIZE_ProblemType; ++i)
70  {
71  if (boost::iequals(ProblemTypeMap[i], ProblemTypeStr))
72  {
74  break;
75  }
76  }
77  }
78  else
79  {
81  }
82 
84  {
87  }
88  else
89  {
90  ASSERTL0(false, "Implicit CFE not set up.");
91  }
92  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
virtual void v_InitObject()
Initialization object for CompressibleFlowSystem class.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DoOdeProjection(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the projection and call the method for imposing the boundary conditions in case of discontinu...
Definition: EulerCFE.cpp:175
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
Length of enum list.
Definition: EulerADCFE.h:47
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
Compute the right-hand side.
Definition: EulerCFE.cpp:144
ProblemType m_problemType
Definition: EulerCFE.h:77
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
ProblemType
Definition: EulerADCFE.h:44
const char *const ProblemTypeMap[]
Definition: EulerADCFE.h:50
void Nektar::EulerCFE::v_SetInitialConditions ( NekDouble  initialtime = 0.0,
bool  dumpInitialConditions = true,
const int  domain = 0 
)
protectedvirtual

Set the initial conditions.

Reimplemented from Nektar::CompressibleFlowSystem.

Definition at line 113 of file EulerCFE.cpp.

References Nektar::SolverUtils::EquationSystem::Checkpoint_Output(), Nektar::eIsentropicVortex, Nektar::SolverUtils::EquationSystem::m_checksteps, m_problemType, SetInitialIsentropicVortex(), Nektar::CompressibleFlowSystem::v_SetInitialConditions(), and Nektar::SolverUtils::EquationSystem::v_SetInitialConditions().

117  {
118  switch (m_problemType)
119  {
120  case eIsentropicVortex:
121  {
122  SetInitialIsentropicVortex(initialtime);
123  break;
124  }
125  default:
126  {
127  EquationSystem::v_SetInitialConditions(initialtime, false);
128  break;
129  }
130  }
131 
133 
134  if (dumpInitialConditions && m_checksteps)
135  {
136  // Dump initial conditions to file
138  }
139  }
void SetInitialIsentropicVortex(NekDouble initialtime)
Set the initial condition for the isentropic vortex problem.
Definition: EulerCFE.cpp:360
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
int m_checksteps
Number of steps between checkpoints.
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Isentropic Vortex.
Definition: EulerCFE.h:46
ProblemType m_problemType
Definition: EulerCFE.h:77
virtual void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)

Friends And Related Function Documentation

friend class MemoryManager< EulerCFE >
friend

Definition at line 61 of file EulerCFE.h.

Member Data Documentation

string Nektar::EulerCFE::className
static
Initial value:
=
"EulerCFE", EulerCFE::create,
"Euler equations in conservative variables without "
"artificial diffusion.")

Name of class.

Definition at line 72 of file EulerCFE.h.

ProblemType Nektar::EulerCFE::m_problemType