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

#include <NavierStokesCFE.h>

Inheritance diagram for Nektar::NavierStokesCFE:
[legend]

Public Member Functions

virtual ~NavierStokesCFE ()
 
- 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...
 
Array< OneD, NekDoubleGetElmtMinHP (void)
 Function to get estimate of min h/p factor per element. More...
 
virtual void GetPressure (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &pressure)
 Extract array with pressure from physfield. More...
 
virtual void GetDensity (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, NekDouble > &density)
 Extract array with density from physfield. More...
 
virtual bool HasConstantDensity ()
 
virtual void GetVelocity (const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &velocity)
 Extract array with velocity from physfield. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
SOLVER_UTILS_EXPORT AdvectionSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
virtual SOLVER_UTILS_EXPORT ~AdvectionSystem ()
 
SOLVER_UTILS_EXPORT AdvectionSharedPtr GetAdvObject ()
 Returns the advection object held by this instance. More...
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetElmtCFLVals (void)
 
SOLVER_UTILS_EXPORT NekDouble GetCFLEstimate (int &elmtid)
 
- 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 >
std::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 ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 
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 SessionFunctionSharedPtr GetFunction (std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
 Get a SessionFunction by name. 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, NekDoubleErrorExtraPoints (unsigned int field)
 Compute error (L2 and L_inf) over an larger set of quadrature points return [L2 Linf]. 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 SessionSummary (SummaryList &vSummary)
 Write out a session summary. More...
 
SOLVER_UTILS_EXPORT Array< OneD, MultiRegions::ExpListSharedPtr > & UpdateFields ()
 
SOLVER_UTILS_EXPORT LibUtilities::FieldMetaDataMapUpdateFieldMetaDataMap ()
 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 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 SetModifiedBasis (const bool modbasis)
 
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, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 

Static Public Attributes

static std::string className
 

Protected Member Functions

 NavierStokesCFE (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
virtual void v_InitObject ()
 Initialization object for CompressibleFlowSystem class. More...
 
virtual void v_DoDiffusion (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
 
virtual void v_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...
 
virtual void v_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...
 
virtual void v_GetFluxPenalty (const Array< OneD, Array< OneD, NekDouble > > &uFwd, const Array< OneD, Array< OneD, NekDouble > > &uBwd, Array< OneD, Array< OneD, NekDouble > > &penaltyCoeff)
 Return the penalty vector for the LDGNS diffusion problem. More...
 
void GetViscosityAndThermalCondFromTemp (const Array< OneD, NekDouble > &temperature, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &thermalCond)
 Update viscosity todo: add artificial viscosity here. More...
 
- Protected Member Functions inherited from Nektar::CompressibleFlowSystem
 CompressibleFlowSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
void InitialiseParameters ()
 Load CFS parameters from the session file. More...
 
void InitAdvection ()
 Create advection and diffusion objects for CFS. 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...
 
void DoAdvection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
 Compute the advection terms for the right-hand side. More...
 
void DoDiffusion (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const Array< OneD, Array< OneD, NekDouble > > &pFwd, const Array< OneD, Array< OneD, NekDouble > > &pBwd)
 Add the diffusions terms to the right-hand side. More...
 
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 SetBoundaryConditions (Array< OneD, Array< OneD, NekDouble > > &physarray, NekDouble time)
 
void GetElmtTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, NekDouble > &tstep)
 Calculate the maximum timestep on each element subject to CFL restrictions. More...
 
virtual NekDouble v_GetTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray)
 Calculate the maximum timestep subject to CFL restrictions. More...
 
virtual void v_SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 Set up logic for residual calculation. More...
 
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)
 
virtual Array< OneD, NekDoublev_GetMaxStdVelocity ()
 Compute the advection velocity in the standard space for each element of the expansion. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate (int step)
 
- Protected Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 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_GenerateSummary (SummaryList &s)
 Print a summary of time stepping parameters. 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 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...
 
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans ()
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time, int &nchk)
 
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, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises EquationSystem class members. More...
 
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...
 
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure (void)
 

Protected Attributes

std::string m_ViscosityType
 
NekDouble m_Cp
 
NekDouble m_Prandtl
 
NekDouble m_muRef
 
NekDouble m_thermalConductivityRef
 
Array< OneD, NekDoublem_mu
 
Array< OneD, NekDoublem_thermalConductivity
 
- Protected Attributes inherited from Nektar::CompressibleFlowSystem
SolverUtils::DiffusionSharedPtr m_diffusion
 
ArtificialDiffusionSharedPtr m_artificialDiffusion
 
Array< OneD, Array< OneD, NekDouble > > m_vecLocs
 
NekDouble m_gamma
 
std::string m_shockCaptureType
 
NekDouble m_filterAlpha
 
NekDouble m_filterExponent
 
NekDouble m_filterCutoff
 
bool m_useFiltering
 
bool m_useLocalTimeStep
 
VariableConverterSharedPtr m_varConv
 
std::vector< CFSBndCondSharedPtrm_bndConds
 
std::vector< SolverUtils::ForcingSharedPtrm_forcing
 
- Protected Attributes inherited from Nektar::SolverUtils::AdvectionSystem
SolverUtils::AdvectionSharedPtr m_advObject
 Advection term. More...
 
- Protected Attributes inherited from Nektar::SolverUtils::UnsteadySystem
int m_infosteps
 Number of time steps between outputting status information. More...
 
int m_abortSteps
 Number of steps between checks for abort conditions. More...
 
int m_filtersInfosteps
 Number of time steps between outputting filters information. More...
 
int m_nanSteps
 
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...
 
NekDouble m_steadyStateTol
 Tolerance to which steady state should be evaluated at. More...
 
int m_steadyStateSteps
 Check for steady state at step interval. More...
 
Array< OneD, Array< OneD, NekDouble > > m_previousSolution
 Storage for previous solution for steady-state check. More...
 
std::ofstream m_errFile
 
std::vector< int > m_intVariables
 
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
 
NekDouble m_filterTimeWarning
 Number of time steps between outputting status information. More...
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 The session reader. More...
 
std::map< std::string, SolverUtils::SessionFunctionSharedPtrm_sessionFunctions
 Map of known SessionFunctions. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 Field input/output. More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 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...
 
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, 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...
 

Friends

class MemoryManager< NavierStokesCFE >
 

Additional Inherited Members

- Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem
NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 
- Protected Types inherited from Nektar::SolverUtils::EquationSystem
enum  HomogeneousType { eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, eNotHomogeneous }
 Parameter for homogeneous expansions. More...
 
- Static Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
static std::string equationSystemTypeLookupIds []
 

Detailed Description

Definition at line 46 of file NavierStokesCFE.h.

Constructor & Destructor Documentation

◆ ~NavierStokesCFE()

Nektar::NavierStokesCFE::~NavierStokesCFE ( )
virtual

Definition at line 54 of file NavierStokesCFE.cpp.

55  {
56 
57  }

◆ NavierStokesCFE()

Nektar::NavierStokesCFE::NavierStokesCFE ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)
protected

Definition at line 46 of file NavierStokesCFE.cpp.

49  : UnsteadySystem(pSession, pGraph),
50  CompressibleFlowSystem(pSession, pGraph)
51  {
52  }
CompressibleFlowSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.

Member Function Documentation

◆ create()

static SolverUtils::EquationSystemSharedPtr Nektar::NavierStokesCFE::create ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)
inlinestatic

Definition at line 52 of file NavierStokesCFE.h.

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

55  {
58  p->InitObject();
59  return p;
60  }
static std::shared_ptr< DataType > AllocateSharedPtr(const Args &...args)
Allocate a shared pointer from the memory pool.
std::shared_ptr< EquationSystem > EquationSystemSharedPtr
A shared pointer to an EquationSystem object.

◆ GetViscosityAndThermalCondFromTemp()

void Nektar::NavierStokesCFE::GetViscosityAndThermalCondFromTemp ( const Array< OneD, NekDouble > &  temperature,
Array< OneD, NekDouble > &  mu,
Array< OneD, NekDouble > &  thermalCond 
)
protected

Update viscosity todo: add artificial viscosity here.

Definition at line 458 of file NavierStokesCFE.cpp.

References Vmath::Fill(), m_Cp, m_muRef, m_Prandtl, Nektar::CompressibleFlowSystem::m_varConv, m_ViscosityType, and Vmath::Smul().

Referenced by v_GetFluxPenalty(), Nektar::NavierStokesCFEAxisym::v_GetViscousFluxVector(), v_GetViscousFluxVector(), and v_GetViscousFluxVectorDeAlias().

462  {
463  int nPts = temperature.num_elements();
464 
465  // Variable viscosity through the Sutherland's law
466  if (m_ViscosityType == "Variable")
467  {
468  m_varConv->GetDynamicViscosity(temperature, mu);
469  }
470  else
471  {
472  Vmath::Fill(nPts, m_muRef, mu, 1);
473  }
474  NekDouble tRa = m_Cp / m_Prandtl;
475  Vmath::Smul(nPts, tRa, mu, 1, thermalCond, 1);
476 
477  }
VariableConverterSharedPtr m_varConv
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:45
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:216
double NekDouble

◆ v_DoDiffusion()

void Nektar::NavierStokesCFE::v_DoDiffusion ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const Array< OneD, Array< OneD, NekDouble > > &  pFwd,
const Array< OneD, Array< OneD, NekDouble > > &  pBwd 
)
protectedvirtual

Reimplemented from Nektar::CompressibleFlowSystem.

Reimplemented in Nektar::NavierStokesCFEAxisym.

Definition at line 121 of file NavierStokesCFE.cpp.

References Nektar::SolverUtils::EquationSystem::GetNpoints(), Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::CompressibleFlowSystem::m_diffusion, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::CompressibleFlowSystem::m_varConv, Nektar::NullNekDoubleArrayofArray, and Vmath::Vadd().

Referenced by Nektar::NavierStokesCFEAxisym::v_DoDiffusion().

126  {
127  int i;
128  int nvariables = inarray.num_elements();
129  int npoints = GetNpoints();
130  int nTracePts = GetTraceTotPoints();
131 
132  Array<OneD, Array<OneD, NekDouble> > outarrayDiff(nvariables);
133 
134  Array<OneD, Array<OneD, NekDouble> > inarrayDiff(nvariables-1);
135  Array<OneD, Array<OneD, NekDouble> > inFwd(nvariables-1);
136  Array<OneD, Array<OneD, NekDouble> > inBwd(nvariables-1);
137 
138  for (i = 0; i < nvariables; ++i)
139  {
140  outarrayDiff[i] = Array<OneD, NekDouble>(npoints);
141  }
142 
143  for (i = 0; i < nvariables-1; ++i)
144  {
145  inarrayDiff[i] = Array<OneD, NekDouble>(npoints);
146  inFwd[i] = Array<OneD, NekDouble>(nTracePts);
147  inBwd[i] = Array<OneD, NekDouble>(nTracePts);
148  }
149 
150  // Extract pressure
151  // (use inarrayDiff[0] as a temporary storage for the pressure)
152  m_varConv->GetPressure(inarray, inarrayDiff[0]);
153 
154  // Extract temperature
155  m_varConv->GetTemperature(inarray, inarrayDiff[nvariables-2]);
156 
157  // Extract velocities
158  m_varConv->GetVelocityVector(inarray, inarrayDiff);
159 
160  // Repeat calculation for trace space
161  if (pFwd == NullNekDoubleArrayofArray ||
163  {
166  }
167  else
168  {
169  m_varConv->GetPressure(pFwd, inFwd[0]);
170  m_varConv->GetPressure(pBwd, inBwd[0]);
171 
172  m_varConv->GetTemperature(pFwd, inFwd[nvariables-2]);
173  m_varConv->GetTemperature(pBwd, inBwd[nvariables-2]);
174 
175  m_varConv->GetVelocityVector(pFwd, inFwd);
176  m_varConv->GetVelocityVector(pBwd, inBwd);
177  }
178 
179  // Diffusion term in physical rhs form
180  m_diffusion->Diffuse(nvariables, m_fields, inarrayDiff, outarrayDiff,
181  inFwd, inBwd);
182 
183  for (i = 0; i < nvariables; ++i)
184  {
185  Vmath::Vadd(npoints,
186  outarrayDiff[i], 1,
187  outarray[i], 1,
188  outarray[i], 1);
189  }
190  }
VariableConverterSharedPtr m_varConv
SolverUtils::DiffusionSharedPtr m_diffusion
SOLVER_UTILS_EXPORT int GetTraceTotPoints()
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
static Array< OneD, Array< OneD, NekDouble > > NullNekDoubleArrayofArray
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:302

◆ v_GetFluxPenalty()

void Nektar::NavierStokesCFE::v_GetFluxPenalty ( const Array< OneD, Array< OneD, NekDouble > > &  uFwd,
const Array< OneD, Array< OneD, NekDouble > > &  uBwd,
Array< OneD, Array< OneD, NekDouble > > &  penaltyCoeff 
)
protectedvirtual

Return the penalty vector for the LDGNS diffusion problem.

Definition at line 415 of file NavierStokesCFE.cpp.

References GetViscosityAndThermalCondFromTemp(), Vmath::Svtsvtp(), Vmath::Vmul(), and Vmath::Vsub().

Referenced by v_InitObject().

419  {
420  unsigned int nTracePts = uFwd[0].num_elements();
421 
422  // Compute average temperature
423  unsigned int nVariables = uFwd.num_elements();
424  Array<OneD, NekDouble> tAve{nTracePts, 0.0};
425  Vmath::Svtsvtp(nTracePts, 0.5, uFwd[nVariables-1], 1,
426  0.5, uBwd[nVariables-1], 1, tAve, 1);
427 
428  // Get average viscosity and thermal conductivity
429  Array<OneD, NekDouble> muAve{nTracePts, 0.0};
430  Array<OneD, NekDouble> tcAve{nTracePts, 0.0};
431 
432  GetViscosityAndThermalCondFromTemp(tAve, muAve, tcAve);
433 
434  // Compute penalty term
435  for (int i = 0; i < nVariables; ++i)
436  {
437  // Get jump of u variables
438  Vmath::Vsub(nTracePts, uFwd[i], 1, uBwd[i], 1, penaltyCoeff[i], 1);
439  // Multiply by variable coefficient = {coeff} ( u^+ - u^- )
440  if ( i < nVariables-1 )
441  {
442  Vmath::Vmul(nTracePts, muAve, 1, penaltyCoeff[i], 1,
443  penaltyCoeff[i], 1);
444  }
445  else
446  {
447  Vmath::Vmul(nTracePts, tcAve, 1, penaltyCoeff[i], 1,
448  penaltyCoeff[i], 1);
449  }
450  }
451  }
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.cpp:346
void GetViscosityAndThermalCondFromTemp(const Array< OneD, NekDouble > &temperature, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &thermalCond)
Update viscosity todo: add artificial viscosity here.
void Svtsvtp(int n, const T alpha, const T *x, int incx, const T beta, const T *y, int incy, T *z, int incz)
vvtvvtp (scalar times vector plus scalar times vector):
Definition: Vmath.cpp:594
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:186

◆ v_GetViscousFluxVector()

void Nektar::NavierStokesCFE::v_GetViscousFluxVector ( const Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  derivativesO1,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  viscousTensor 
)
protectedvirtual

Return the flux vector for the LDG diffusion problem.

Todo:
Complete the viscous flux vector

Reimplemented in Nektar::NavierStokesCFEAxisym.

Definition at line 196 of file NavierStokesCFE.cpp.

References GetViscosityAndThermalCondFromTemp(), m_mu, Nektar::SolverUtils::EquationSystem::m_spacedim, m_thermalConductivity, Vmath::Smul(), Vmath::Vadd(), Vmath::Vcopy(), Vmath::Vmul(), Vmath::Vvtvp(), and Vmath::Zero().

Referenced by v_InitObject().

200  {
201  // Auxiliary variables
202  int nScalar = physfield.num_elements();
203  int nPts = physfield[0].num_elements();
204  Array<OneD, NekDouble > divVel (nPts, 0.0);
205 
206  // Stokes hypothesis
207  const NekDouble lambda = -2.0/3.0;
208 
209  // Update viscosity and thermal conductivity
210  GetViscosityAndThermalCondFromTemp(physfield[nScalar-1], m_mu,
212 
213  // Velocity divergence
214  for (int j = 0; j < m_spacedim; ++j)
215  {
216  Vmath::Vadd(nPts, divVel, 1, derivativesO1[j][j], 1,
217  divVel, 1);
218  }
219 
220  // Velocity divergence scaled by lambda * mu
221  Vmath::Smul(nPts, lambda, divVel, 1, divVel, 1);
222  Vmath::Vmul(nPts, m_mu, 1, divVel, 1, divVel, 1);
223 
224  // Viscous flux vector for the rho equation = 0
225  for (int i = 0; i < m_spacedim; ++i)
226  {
227  Vmath::Zero(nPts, viscousTensor[i][0], 1);
228  }
229 
230  // Viscous stress tensor (for the momentum equations)
231  for (int i = 0; i < m_spacedim; ++i)
232  {
233  for (int j = i; j < m_spacedim; ++j)
234  {
235  Vmath::Vadd(nPts, derivativesO1[i][j], 1,
236  derivativesO1[j][i], 1,
237  viscousTensor[i][j+1], 1);
238 
239  Vmath::Vmul(nPts, m_mu, 1,
240  viscousTensor[i][j+1], 1,
241  viscousTensor[i][j+1], 1);
242 
243  if (i == j)
244  {
245  // Add divergence term to diagonal
246  Vmath::Vadd(nPts, viscousTensor[i][j+1], 1,
247  divVel, 1,
248  viscousTensor[i][j+1], 1);
249  }
250  else
251  {
252  // Copy to make symmetric
253  Vmath::Vcopy(nPts, viscousTensor[i][j+1], 1,
254  viscousTensor[j][i+1], 1);
255  }
256  }
257  }
258 
259  // Terms for the energy equation
260  for (int i = 0; i < m_spacedim; ++i)
261  {
262  Vmath::Zero(nPts, viscousTensor[i][m_spacedim+1], 1);
263  // u_j * tau_ij
264  for (int j = 0; j < m_spacedim; ++j)
265  {
266  Vmath::Vvtvp(nPts, physfield[j], 1,
267  viscousTensor[i][j+1], 1,
268  viscousTensor[i][m_spacedim+1], 1,
269  viscousTensor[i][m_spacedim+1], 1);
270  }
271  // Add k*T_i
273  derivativesO1[i][m_spacedim], 1,
274  viscousTensor[i][m_spacedim+1], 1,
275  viscousTensor[i][m_spacedim+1], 1);
276  }
277  }
Array< OneD, NekDouble > m_thermalConductivity
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:445
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:216
int m_spacedim
Spatial dimension (>= expansion dim).
double NekDouble
void GetViscosityAndThermalCondFromTemp(const Array< OneD, NekDouble > &temperature, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &thermalCond)
Update viscosity todo: add artificial viscosity here.
Array< OneD, NekDouble > m_mu
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:376
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1064
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:302
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:186

◆ v_GetViscousFluxVectorDeAlias()

void Nektar::NavierStokesCFE::v_GetViscousFluxVectorDeAlias ( const Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  derivativesO1,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  viscousTensor 
)
protectedvirtual

Return the flux vector for the LDG diffusion problem.

Todo:
Complete the viscous flux vector

Reimplemented in Nektar::NavierStokesCFEAxisym.

Definition at line 283 of file NavierStokesCFE.cpp.

References GetViscosityAndThermalCondFromTemp(), Nektar::SolverUtils::EquationSystem::m_fields, m_mu, Nektar::SolverUtils::EquationSystem::m_spacedim, m_thermalConductivity, Vmath::Smul(), Vmath::Vadd(), Vmath::Vmul(), Vmath::Vvtvp(), and Vmath::Zero().

Referenced by v_InitObject().

287  {
288  // Factor to rescale 1d points in dealiasing.
289  NekDouble OneDptscale = 2;
290  // Get number of points to dealias a cubic non-linearity
291  int nScalar = physfield.num_elements();
292  int nPts = m_fields[0]->Get1DScaledTotPoints(OneDptscale);
293  int nPts_orig = physfield[0].num_elements();
294 
295  // Auxiliary variables
296  Array<OneD, NekDouble > divVel (nPts, 0.0);
297 
298  // Stokes hypothesis
299  const NekDouble lambda = -2.0/3.0;
300 
301  // Update viscosity and thermal conductivity
302  GetViscosityAndThermalCondFromTemp(physfield[nScalar-1], m_mu,
304 
305  // Interpolate inputs and initialise interpolated output
306  Array<OneD, Array<OneD, NekDouble> > vel_interp(m_spacedim);
307  Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
308  deriv_interp(m_spacedim);
309  Array<OneD, Array<OneD, Array<OneD, NekDouble> > >
310  out_interp(m_spacedim);
311  for (int i = 0; i < m_spacedim; ++i)
312  {
313  // Interpolate velocity
314  vel_interp[i] = Array<OneD, NekDouble> (nPts);
315  m_fields[0]->PhysInterp1DScaled(
316  OneDptscale, physfield[i], vel_interp[i]);
317 
318  // Interpolate derivatives
319  deriv_interp[i] = Array<OneD,Array<OneD,NekDouble> > (m_spacedim+1);
320  for (int j = 0; j < m_spacedim+1; ++j)
321  {
322  deriv_interp[i][j] = Array<OneD, NekDouble> (nPts);
323  m_fields[0]->PhysInterp1DScaled(
324  OneDptscale, derivativesO1[i][j], deriv_interp[i][j]);
325  }
326 
327  // Output (start from j=1 since flux is zero for rho)
328  out_interp[i] = Array<OneD,Array<OneD,NekDouble> > (m_spacedim+2);
329  for (int j = 1; j < m_spacedim+2; ++j)
330  {
331  out_interp[i][j] = Array<OneD, NekDouble> (nPts);
332  }
333  }
334 
335  // Velocity divergence
336  for (int j = 0; j < m_spacedim; ++j)
337  {
338  Vmath::Vadd(nPts, divVel, 1, deriv_interp[j][j], 1,
339  divVel, 1);
340  }
341 
342  // Velocity divergence scaled by lambda * mu
343  Vmath::Smul(nPts, lambda, divVel, 1, divVel, 1);
344  Vmath::Vmul(nPts, m_mu, 1, divVel, 1, divVel, 1);
345 
346  // Viscous flux vector for the rho equation = 0 (no need to dealias)
347  for (int i = 0; i < m_spacedim; ++i)
348  {
349  Vmath::Zero(nPts_orig, viscousTensor[i][0], 1);
350  }
351 
352  // Viscous stress tensor (for the momentum equations)
353  for (int i = 0; i < m_spacedim; ++i)
354  {
355  for (int j = i; j < m_spacedim; ++j)
356  {
357  Vmath::Vadd(nPts, deriv_interp[i][j], 1,
358  deriv_interp[j][i], 1,
359  out_interp[i][j+1], 1);
360 
361  Vmath::Vmul(nPts, m_mu, 1,
362  out_interp[i][j+1], 1,
363  out_interp[i][j+1], 1);
364 
365  if (i == j)
366  {
367  // Add divergence term to diagonal
368  Vmath::Vadd(nPts, out_interp[i][j+1], 1,
369  divVel, 1,
370  out_interp[i][j+1], 1);
371  }
372  else
373  {
374  // Make symmetric
375  out_interp[j][i+1] = out_interp[i][j+1];
376  }
377  }
378  }
379 
380  // Terms for the energy equation
381  for (int i = 0; i < m_spacedim; ++i)
382  {
383  Vmath::Zero(nPts, out_interp[i][m_spacedim+1], 1);
384  // u_j * tau_ij
385  for (int j = 0; j < m_spacedim; ++j)
386  {
387  Vmath::Vvtvp(nPts, vel_interp[j], 1,
388  out_interp[i][j+1], 1,
389  out_interp[i][m_spacedim+1], 1,
390  out_interp[i][m_spacedim+1], 1);
391  }
392  // Add k*T_i
394  deriv_interp[i][m_spacedim], 1,
395  out_interp[i][m_spacedim+1], 1,
396  out_interp[i][m_spacedim+1], 1);
397  }
398 
399  // Project to original space
400  for (int i = 0; i < m_spacedim; ++i)
401  {
402  for (int j = 1; j < m_spacedim+2; ++j)
403  {
404  m_fields[0]->PhysGalerkinProjection1DScaled(
405  OneDptscale,
406  out_interp[i][j],
407  viscousTensor[i][j]);
408  }
409  }
410  }
Array< OneD, NekDouble > m_thermalConductivity
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:445
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*y.
Definition: Vmath.cpp:216
int m_spacedim
Spatial dimension (>= expansion dim).
double NekDouble
void GetViscosityAndThermalCondFromTemp(const Array< OneD, NekDouble > &temperature, Array< OneD, NekDouble > &mu, Array< OneD, NekDouble > &thermalCond)
Update viscosity todo: add artificial viscosity here.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
Array< OneD, NekDouble > m_mu
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:376
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:302
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:186

◆ v_InitObject()

void Nektar::NavierStokesCFE::v_InitObject ( )
protectedvirtual

Initialization object for CompressibleFlowSystem class.

Reimplemented from Nektar::CompressibleFlowSystem.

Reimplemented in Nektar::NavierStokesCFEAxisym.

Definition at line 62 of file NavierStokesCFE.cpp.

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::SolverUtils::GetDiffusionFactory(), m_Cp, Nektar::CompressibleFlowSystem::m_diffusion, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::CompressibleFlowSystem::m_gamma, m_mu, m_muRef, m_Prandtl, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_specHP_dealiasing, m_thermalConductivity, m_thermalConductivityRef, m_ViscosityType, v_GetFluxPenalty(), v_GetViscousFluxVector(), v_GetViscousFluxVectorDeAlias(), and Nektar::CompressibleFlowSystem::v_InitObject().

Referenced by Nektar::NavierStokesCFEAxisym::v_InitObject().

63  {
65 
66  // Get gas constant from session file and compute Cp
67  NekDouble gasConstant;
68  m_session->LoadParameter ("GasConstant", gasConstant, 287.058);
69  m_Cp = m_gamma / (m_gamma - 1.0) * gasConstant;
70 
71  // Viscosity
72  int nPts = m_fields[0]->GetNpoints();
73  m_session->LoadSolverInfo("ViscosityType", m_ViscosityType, "Constant");
74  m_session->LoadParameter ("mu", m_muRef, 1.78e-05);
75  m_mu = Array<OneD, NekDouble>(nPts, m_muRef);
76 
77  // Thermal conductivity or Prandtl
78  if( m_session->DefinesParameter("thermalConductivity"))
79  {
80  ASSERTL0( !m_session->DefinesParameter("Pr"),
81  "Cannot define both Pr and thermalConductivity.");
82 
83  m_session->LoadParameter ("thermalConductivity",
86  }
87  else
88  {
89  m_session->LoadParameter ("Pr", m_Prandtl, 0.72);
91  }
93  Array<OneD, NekDouble>(nPts, m_thermalConductivityRef);
94 
95  string diffName;
96  m_session->LoadSolverInfo("DiffusionType", diffName, "LDGNS");
97 
99  .CreateInstance(diffName, diffName);
100 
102  {
103  m_diffusion->SetFluxVectorNS(
105  this);
106  }
107  else
108  {
109  m_diffusion->SetFluxVectorNS(&NavierStokesCFE::
110  v_GetViscousFluxVector, this);
111  }
112 
113  // Set up penalty term for LDGNS
114  m_diffusion->SetFluxPenaltyNS(&NavierStokesCFE::
115  v_GetFluxPenalty, this);
116 
117  // Concluding initialisation of diffusion operator
118  m_diffusion->InitObject (m_session, m_fields);
119  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:216
virtual void v_InitObject()
Initialization object for CompressibleFlowSystem class.
Array< OneD, NekDouble > m_thermalConductivity
NekDouble m_thermalConductivityRef
NavierStokesCFE(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
DiffusionFactory & GetDiffusionFactory()
Definition: Diffusion.cpp:41
bool m_specHP_dealiasing
Flag to determine if dealisising is usde for the Spectral/hp element discretisation.
virtual void v_GetFluxPenalty(const Array< OneD, Array< OneD, NekDouble > > &uFwd, const Array< OneD, Array< OneD, NekDouble > > &uBwd, Array< OneD, Array< OneD, NekDouble > > &penaltyCoeff)
Return the penalty vector for the LDGNS diffusion problem.
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
double NekDouble
SolverUtils::DiffusionSharedPtr m_diffusion
virtual void v_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.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
virtual void v_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.
Array< OneD, NekDouble > m_mu

Friends And Related Function Documentation

◆ MemoryManager< NavierStokesCFE >

friend class MemoryManager< NavierStokesCFE >
friend

Definition at line 49 of file NavierStokesCFE.h.

Member Data Documentation

◆ className

string Nektar::NavierStokesCFE::className
static
Initial value:
=
"NavierStokesCFE", NavierStokesCFE::create,
"NavierStokes equations in conservative variables.")

Definition at line 62 of file NavierStokesCFE.h.

◆ m_Cp

NekDouble Nektar::NavierStokesCFE::m_Cp
protected

Definition at line 68 of file NavierStokesCFE.h.

Referenced by GetViscosityAndThermalCondFromTemp(), and v_InitObject().

◆ m_mu

Array<OneD, NekDouble> Nektar::NavierStokesCFE::m_mu
protected

◆ m_muRef

NekDouble Nektar::NavierStokesCFE::m_muRef
protected

Definition at line 71 of file NavierStokesCFE.h.

Referenced by GetViscosityAndThermalCondFromTemp(), and v_InitObject().

◆ m_Prandtl

NekDouble Nektar::NavierStokesCFE::m_Prandtl
protected

Definition at line 69 of file NavierStokesCFE.h.

Referenced by GetViscosityAndThermalCondFromTemp(), and v_InitObject().

◆ m_thermalConductivity

Array<OneD, NekDouble> Nektar::NavierStokesCFE::m_thermalConductivity
protected

◆ m_thermalConductivityRef

NekDouble Nektar::NavierStokesCFE::m_thermalConductivityRef
protected

Definition at line 72 of file NavierStokesCFE.h.

Referenced by v_InitObject().

◆ m_ViscosityType

std::string Nektar::NavierStokesCFE::m_ViscosityType
protected

Definition at line 67 of file NavierStokesCFE.h.

Referenced by GetViscosityAndThermalCondFromTemp(), and v_InitObject().