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

Base class for unsteady solvers. More...

#include <UnsteadySystem.h>

Inheritance diagram for Nektar::SolverUtils::UnsteadySystem:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::UnsteadySystem:
Collaboration graph
[legend]

Public Member Functions

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...
 

Public Attributes

NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 

Protected Member Functions

SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession)
 Initialises UnsteadySystem class members. More...
 
virtual SOLVER_UTILS_EXPORT void v_InitObject ()
 Init object for UnsteadySystem class. 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 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
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_PostIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_SteadyStateCheck (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)
 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...
 
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 
SOLVER_UTILS_EXPORT void EvaluateFunctionExp (std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
 
SOLVER_UTILS_EXPORT void EvaluateFunctionFld (std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
 
SOLVER_UTILS_EXPORT void EvaluateFunctionPts (std::string pFieldName, Array< OneD, NekDouble > &pArray, const std::string &pFunctionName, const NekDouble &pTime=0.0, const int domain=0)
 
SOLVER_UTILS_EXPORT void LoadPts (std::string funcFilename, std::string filename, Nektar::LibUtilities::PtsFieldSharedPtr &outPts)
 
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)
 
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 

Protected Attributes

int m_infosteps
 Number of time steps between outputting status 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...
 
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,
FieldUtils::Interpolator
m_interpolators
 Map of interpolator objects. More...
 
std::map< std::string,
std::pair< std::string,
LibUtilities::PtsFieldSharedPtr > > 
m_loadedPtsFields
 pts fields we already read from disk: {funcFilename: (filename, ptsfield)} More...
 
std::map< std::string,
std::pair< std::string,
loadedFldField > > 
m_loadedFldFields
 
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...
 
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...
 

Private Member Functions

void WeakPenaltyforScalar (const int var, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble time=0.0)
 
void WeakPenaltyforVector (const int var, const int dir, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble C11, NekDouble time=0.0)
 

Additional Inherited Members

- Protected Types inherited from Nektar::SolverUtils::EquationSystem
enum  HomogeneousType { eHomogeneous1D, eHomogeneous2D, eHomogeneous3D, eNotHomogeneous }
 Parameter for homogeneous expansions. More...
 

Detailed Description

Base class for unsteady solvers.

Provides the underlying timestepping framework for unsteady solvers including the general timestepping routines. This class is not intended to be directly instantiated, but rather is a base class on which to define unsteady solvers.

For details on implementing unsteady solvers see sectionADRSolverModuleImplementation here

Definition at line 48 of file UnsteadySystem.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::UnsteadySystem::~UnsteadySystem ( )
virtual

Destructor.

Destructor for the class UnsteadyAdvection.

Definition at line 128 of file UnsteadySystem.cpp.

129  {
130  }
Nektar::SolverUtils::UnsteadySystem::UnsteadySystem ( const LibUtilities::SessionReaderSharedPtr pSession)
protected

Initialises UnsteadySystem class members.

Processes SolverInfo parameters from the session file and sets up timestepping-specific code.

Parameters
pSessionSession object to read parameters from.

Definition at line 67 of file UnsteadySystem.cpp.

69  : EquationSystem(pSession),
70  m_infosteps(10)
71 
72  {
73  }
int m_infosteps
Number of time steps between outputting status information.
SOLVER_UTILS_EXPORT EquationSystem(const LibUtilities::SessionReaderSharedPtr &pSession)
Initialises EquationSystem class members.

Member Function Documentation

void Nektar::SolverUtils::UnsteadySystem::CheckForRestartTime ( NekDouble time,
int &  nchk 
)
protected

Definition at line 702 of file UnsteadySystem.cpp.

References Nektar::LibUtilities::FieldIO::CreateForFile(), Nektar::LibUtilities::eFunctionTypeFile, Nektar::iterator, Nektar::SolverUtils::EquationSystem::m_fieldMetaDataMap, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_session, Nektar::LibUtilities::NullFieldMetaDataMap, and Nektar::LibUtilities::PortablePath().

Referenced by v_DoInitialise().

703  {
704  if (m_session->DefinesFunction("InitialConditions"))
705  {
706  for (int i = 0; i < m_fields.num_elements(); ++i)
707  {
709 
710  vType = m_session->GetFunctionType(
711  "InitialConditions", m_session->GetVariable(i));
712 
713  if (vType == LibUtilities::eFunctionTypeFile)
714  {
715  std::string filename
716  = m_session->GetFunctionFilename(
717  "InitialConditions", m_session->GetVariable(i));
718 
719  fs::path pfilename(filename);
720 
721  // redefine path for parallel file which is in directory
722  if(fs::is_directory(pfilename))
723  {
724  fs::path metafile("Info.xml");
725  fs::path fullpath = pfilename / metafile;
726  filename = LibUtilities::PortablePath(fullpath);
727  }
730  m_session, filename);
731  fld->ImportFieldMetaData(filename, m_fieldMetaDataMap);
732 
733  // check to see if time defined
734  if (m_fieldMetaDataMap !=
736  {
738 
739  iter = m_fieldMetaDataMap.find("Time");
740  if (iter != m_fieldMetaDataMap.end())
741  {
742  time = boost::lexical_cast<NekDouble>(
743  iter->second);
744  }
745 
746  iter = m_fieldMetaDataMap.find("ChkFileNum");
747  if (iter != m_fieldMetaDataMap.end())
748  {
749  nchk = boost::lexical_cast<NekDouble>(
750  iter->second);
751  }
752  }
753 
754  break;
755  }
756  }
757  }
758  }
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
Map to identify relevant solver info to dump in output fields.
std::string PortablePath(const boost::filesystem::path &path)
create portable path on different platforms for boost::filesystem path
Definition: FileSystem.cpp:41
boost::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:309
double NekDouble
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
static boost::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:212
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:55
NekDouble Nektar::SolverUtils::UnsteadySystem::GetTimeStep ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray)

Calculate the larger time-step mantaining the problem stable.

Return the timestep to be used for the next step in the time-marching loop.

This function can be overloaded to facilitate solver which utilise a CFL (or other) parameter to determine a maximum timestep under which the problem remains stable.

Definition at line 919 of file UnsteadySystem.cpp.

References v_GetTimeStep().

921  {
922  return v_GetTimeStep(inarray);
923  }
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.
NekDouble Nektar::SolverUtils::UnsteadySystem::MaxTimeStepEstimator ( )
protected

Get the maximum timestep estimator for cfl control.

Returns the maximum time estimator for CFL control.

Definition at line 135 of file UnsteadySystem.cpp.

References ASSERTL0, Nektar::LibUtilities::eAdamsBashforthOrder1, Nektar::LibUtilities::eAdamsBashforthOrder2, Nektar::LibUtilities::eClassicalRungeKutta4, Nektar::LibUtilities::eForwardEuler, Nektar::LibUtilities::eMidpoint, Nektar::LibUtilities::eRungeKutta2, Nektar::LibUtilities::eRungeKutta2_ImprovedEuler, Nektar::LibUtilities::eRungeKutta2_SSP, Nektar::LibUtilities::eRungeKutta3_SSP, Nektar::LibUtilities::eRungeKutta4, and m_intScheme.

Referenced by Nektar::APE::GetCFLEstimate(), and Nektar::CompressibleFlowSystem::v_GetTimeStep().

136  {
137  NekDouble TimeStability = 0.0;
138  switch(m_intScheme->GetIntegrationMethod())
139  {
143  {
144  TimeStability = 2.784;
145  break;
146  }
153  {
154  TimeStability = 2.0;
155  break;
156  }
158  {
159  TimeStability = 1.0;
160  break;
161  }
162  default:
163  {
164  ASSERTL0(
165  false,
166  "No CFL control implementation for this time"
167  "integration scheme");
168  }
169  }
170  return TimeStability;
171  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
Adams-Bashforth Forward multi-step scheme of order 2.
Runge-Kutta multi-stage scheme 4th order explicit (old name)
Classical RungeKutta4 method (new name for eClassicalRungeKutta4)
Nonlinear SSP RungeKutta3 explicit.
Classical RungeKutta2 method (new name for eMidpoint)
Adams-Bashforth Forward multi-step scheme of order 1.
Nonlinear SSP RungeKutta2 explicit (surrogate for eRungeKutta2_ImprovedEuler)
Improved RungeKutta2 explicit (old name meaning Heun's method)
double NekDouble
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
Wrapper to the time integration scheme.
midpoint method (old name)
void Nektar::SolverUtils::UnsteadySystem::SVVVarDiffCoeff ( const Array< OneD, Array< OneD, NekDouble > >  vel,
StdRegions::VarCoeffMap varCoeffMap 
)
protected

Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity.

Definition at line 953 of file UnsteadySystem.cpp.

References Nektar::StdRegions::eVarCoeffLaplacian, Nektar::SolverUtils::EquationSystem::m_fields, Vmath::Smul(), Vmath::Vmul(), Vmath::Vsqrt(), and Vmath::Vvtvp().

Referenced by Nektar::UnsteadyViscousBurgers::DoImplicitSolve(), and Nektar::UnsteadyViscousBurgers::v_InitObject().

956  {
957  int phystot = m_fields[0]->GetTotPoints();
958  int nvel = vel.num_elements();
959 
960  Array<OneD, NekDouble> varcoeff(phystot),tmp;
961 
962  // calculate magnitude of v
963  Vmath::Vmul(phystot,vel[0],1,vel[0],1,varcoeff,1);
964  for(int n = 1; n < nvel; ++n)
965  {
966  Vmath::Vvtvp(phystot,vel[n],1,vel[n],1,varcoeff,1,varcoeff,1);
967  }
968  Vmath::Vsqrt(phystot,varcoeff,1,varcoeff,1);
969 
970  for(int i = 0; i < m_fields[0]->GetNumElmts(); ++i)
971  {
972  int offset = m_fields[0]->GetPhys_Offset(i);
973  int nq = m_fields[0]->GetExp(i)->GetTotPoints();
974  Array<OneD, NekDouble> unit(nq,1.0);
975 
976  int nmodes = 0;
977 
978  for(int n = 0; n < m_fields[0]->GetExp(i)->GetNumBases(); ++n)
979  {
980  nmodes = max(nmodes,
981  m_fields[0]->GetExp(i)->GetBasisNumModes(n));
982  }
983 
984  NekDouble h = m_fields[0]->GetExp(i)->Integral(unit);
985  h = pow(h,(NekDouble) (1.0/nvel))/((NekDouble) nmodes);
986 
987  Vmath::Smul(nq,h,varcoeff+offset,1,tmp = varcoeff+offset,1);
988  }
989 
990  // set up map with eVarCoffLaplacian key
991  varCoeffMap[StdRegions::eVarCoeffLaplacian] = varcoeff;
992  }
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.cpp:408
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:442
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:213
double NekDouble
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
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:183
void Nektar::SolverUtils::UnsteadySystem::v_AppendOutput1D ( Array< OneD, Array< OneD, NekDouble > > &  solution1D)
protectedvirtual

Print the solution at each solution point in a txt file.

Stores the solution in a file for 1D problems only. This method has been implemented to facilitate the post-processing for 1D problems.

Definition at line 506 of file UnsteadySystem.cpp.

References Nektar::SolverUtils::EquationSystem::GetNpoints(), and Nektar::SolverUtils::EquationSystem::m_fields.

Referenced by v_DoSolve().

508  {
509  // Coordinates of the quadrature points in the real physical space
510  Array<OneD,NekDouble> x(GetNpoints());
511  Array<OneD,NekDouble> y(GetNpoints());
512  Array<OneD,NekDouble> z(GetNpoints());
513  m_fields[0]->GetCoords(x, y, z);
514 
515  // Print out the solution in a txt file
516  ofstream outfile;
517  outfile.open("solution1D.txt");
518  for(int i = 0; i < GetNpoints(); i++)
519  {
520  outfile << scientific << setw (17) << setprecision(16) << x[i]
521  << " " << solution1D[0][i] << endl;
522  }
523  outfile << endl << endl;
524  outfile.close();
525  }
SOLVER_UTILS_EXPORT int GetNpoints()
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Nektar::SolverUtils::UnsteadySystem::v_DoInitialise ( void  )
protectedvirtual

Sets up initial conditions.

Sets the initial conditions.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::CoupledLinearNS, Nektar::PulseWaveSystem, Nektar::VelocityCorrectionScheme, and Nektar::VCSMapping.

Definition at line 461 of file UnsteadySystem.cpp.

References CheckForRestartTime(), Nektar::SolverUtils::EquationSystem::m_nchk, Nektar::SolverUtils::EquationSystem::m_time, Nektar::SolverUtils::EquationSystem::SetBoundaryConditions(), and Nektar::SolverUtils::EquationSystem::SetInitialConditions().

Referenced by Nektar::VCSMapping::v_DoInitialise(), and Nektar::VelocityCorrectionScheme::v_DoInitialise().

462  {
466  }
NekDouble m_time
Current time of simulation.
SOLVER_UTILS_EXPORT void CheckForRestartTime(NekDouble &time, int &nchk)
int m_nchk
Number of checkpoints written so far.
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
SOLVER_UTILS_EXPORT void SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Initialise the data in the dependent fields.
void Nektar::SolverUtils::UnsteadySystem::v_DoSolve ( void  )
protectedvirtual

Solves an unsteady problem.

Initialises the time integration scheme (as specified in the session file), and perform the time integration.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::CoupledLinearNS, and Nektar::PulseWaveSystem.

Definition at line 177 of file UnsteadySystem.cpp.

References ASSERTL0, Nektar::SolverUtils::EquationSystem::Checkpoint_Output(), Nektar::SolverUtils::EquationSystem::eNotHomogeneous, Nektar::SolverUtils::EquationSystem::GetTimeStep(), Nektar::iterator, Nektar::NekConstants::kNekZeroTol, m_cflSafetyFactor, Nektar::SolverUtils::EquationSystem::m_checksteps, Nektar::SolverUtils::EquationSystem::m_checktime, Nektar::SolverUtils::EquationSystem::m_fields, m_filters, Nektar::SolverUtils::EquationSystem::m_fintime, Nektar::SolverUtils::EquationSystem::m_HomogeneousType, m_homoInitialFwd, m_infosteps, Nektar::SolverUtils::EquationSystem::m_initialStep, m_intScheme, m_intSoln, m_intVariables, m_nanSteps, Nektar::SolverUtils::EquationSystem::m_nchk, m_ode, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_spacedim, Nektar::SolverUtils::EquationSystem::m_steps, Nektar::SolverUtils::EquationSystem::m_time, Nektar::SolverUtils::EquationSystem::m_timestep, Vmath::Nnan(), Nektar::LibUtilities::ReduceMax, Nektar::Timer::Start(), Nektar::Timer::Stop(), Nektar::Timer::TimePerTest(), v_AppendOutput1D(), v_PostIntegrate(), v_PreIntegrate(), and v_RequireFwdTrans().

Referenced by Nektar::CoupledLinearNS::v_DoSolve().

178  {
179  ASSERTL0(m_intScheme != 0, "No time integration scheme.");
180 
181  int i = 1;
182  int nvariables = 0;
183  int nfields = m_fields.num_elements();
184 
185  if (m_intVariables.empty())
186  {
187  for (i = 0; i < nfields; ++i)
188  {
189  m_intVariables.push_back(i);
190  }
191  nvariables = nfields;
192  }
193  else
194  {
195  nvariables = m_intVariables.size();
196  }
197 
198  // Integrate in wave-space if using homogeneous1D
200  {
201  for(i = 0; i < nfields; ++i)
202  {
203  m_fields[i]->HomogeneousFwdTrans(m_fields[i]->GetPhys(),
204  m_fields[i]->UpdatePhys());
205  m_fields[i]->SetWaveSpace(true);
206  m_fields[i]->SetPhysState(false);
207  }
208  }
209 
210  // Set up wrapper to fields data storage.
211  Array<OneD, Array<OneD, NekDouble> > fields(nvariables);
212  Array<OneD, Array<OneD, NekDouble> > tmp (nvariables);
213 
214  // Order storage to list time-integrated fields first.
215  for(i = 0; i < nvariables; ++i)
216  {
217  fields[i] = m_fields[m_intVariables[i]]->GetPhys();
218  m_fields[m_intVariables[i]]->SetPhysState(false);
219  }
220 
221  // Initialise time integration scheme
222  m_intSoln = m_intScheme->InitializeScheme(
223  m_timestep, fields, m_time, m_ode);
224 
225  // Initialise filters
227  for (x = m_filters.begin(); x != m_filters.end(); ++x)
228  {
229  (*x)->Initialise(m_fields, m_time);
230  }
231 
232  // Ensure that there is no conflict of parameters
233  if(m_cflSafetyFactor > 0.0)
234  {
235  // Check final condition
236  ASSERTL0(m_fintime == 0.0 || m_steps == 0,
237  "Final condition not unique: "
238  "fintime > 0.0 and Nsteps > 0");
239 
240  // Check timestep condition
241  ASSERTL0(m_timestep == 0.0,
242  "Timestep not unique: timestep > 0.0 & CFL > 0.0");
243  }
244 
245  // Check uniqueness of checkpoint output
246  ASSERTL0((m_checktime == 0.0 && m_checksteps == 0) ||
247  (m_checktime > 0.0 && m_checksteps == 0) ||
248  (m_checktime == 0.0 && m_checksteps > 0),
249  "Only one of IO_CheckTime and IO_CheckSteps "
250  "should be set!");
251 
252  Timer timer;
253  bool doCheckTime = false;
254  int step = m_initialStep;
255  int stepCounter = 0;
256  NekDouble intTime = 0.0;
257  NekDouble lastCheckTime = 0.0;
258  NekDouble cpuTime = 0.0;
259  NekDouble elapsed = 0.0;
260 
261  while (step < m_steps ||
263  {
264  if (m_cflSafetyFactor)
265  {
266  m_timestep = GetTimeStep(fields);
267 
268  // Ensure that the final timestep finishes at the final
269  // time, or at a prescribed IO_CheckTime.
270  if (m_time + m_timestep > m_fintime && m_fintime > 0.0)
271  {
273  }
274  else if (m_checktime &&
275  m_time + m_timestep - lastCheckTime >= m_checktime)
276  {
277  lastCheckTime += m_checktime;
278  m_timestep = lastCheckTime - m_time;
279  doCheckTime = true;
280  }
281  }
282 
283  // Perform any solver-specific pre-integration steps
284  timer.Start();
285  if (v_PreIntegrate(step))
286  {
287  break;
288  }
289 
290  fields = m_intScheme->TimeIntegrate(
291  stepCounter, m_timestep, m_intSoln, m_ode);
292  timer.Stop();
293 
294  m_time += m_timestep;
295  elapsed = timer.TimePerTest(1);
296  intTime += elapsed;
297  cpuTime += elapsed;
298 
299  // Write out status information
300  if (m_session->GetComm()->GetRank() == 0 &&
301  !((step+1) % m_infosteps))
302  {
303  cout << "Steps: " << setw(8) << left << step+1 << " "
304  << "Time: " << setw(12) << left << m_time;
305 
306  if (m_cflSafetyFactor)
307  {
308  cout << " Time-step: " << setw(12)
309  << left << m_timestep;
310  }
311 
312  stringstream ss;
313  ss << cpuTime << "s";
314  cout << " CPU Time: " << setw(8) << left
315  << ss.str() << endl;
316  cpuTime = 0.0;
317  }
318 
319  // Transform data into coefficient space
320  for (i = 0; i < nvariables; ++i)
321  {
322  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
323  if( v_RequireFwdTrans() )
324  {
325  m_fields[m_intVariables[i]]->FwdTrans_IterPerExp(
326  fields[i],
327  m_fields[m_intVariables[i]]->UpdateCoeffs());
328  }
329  m_fields[m_intVariables[i]]->SetPhysState(false);
330  }
331 
332  // Perform any solver-specific post-integration steps
333  if (v_PostIntegrate(step))
334  {
335  break;
336  }
337 
338  // search for NaN and quit if found
339  if (m_nanSteps && !((step+1) % m_nanSteps) )
340  {
341  int nanFound = 0;
342  for (i = 0; i < nvariables; ++i)
343  {
344  if (Vmath::Nnan(fields[i].num_elements(),
345  fields[i], 1) > 0)
346  {
347  nanFound = 1;
348  }
349  }
350  m_session->GetComm()->AllReduce(nanFound,
352  ASSERTL0 (!nanFound,
353  "NaN found during time integration.");
354  }
355  // Update filters
357  for (x = m_filters.begin(); x != m_filters.end(); ++x)
358  {
359  (*x)->Update(m_fields, m_time);
360  }
361 
362  // Write out checkpoint files
363  if ((m_checksteps && !((step + 1) % m_checksteps)) ||
364  doCheckTime)
365  {
367  {
368  vector<bool> transformed(nfields, false);
369  for(i = 0; i < nfields; i++)
370  {
371  if (m_fields[i]->GetWaveSpace())
372  {
373  m_fields[i]->SetWaveSpace(false);
374  m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
375  m_fields[i]->UpdatePhys());
376  m_fields[i]->SetPhysState(true);
377  transformed[i] = true;
378  }
379  }
381  m_nchk++;
382  for(i = 0; i < nfields; i++)
383  {
384  if (transformed[i])
385  {
386  m_fields[i]->SetWaveSpace(true);
387  m_fields[i]->HomogeneousFwdTrans(
388  m_fields[i]->GetPhys(),
389  m_fields[i]->UpdatePhys());
390  m_fields[i]->SetPhysState(false);
391  }
392  }
393  }
394  else
395  {
397  m_nchk++;
398  }
399  doCheckTime = false;
400  }
401 
402  // Step advance
403  ++step;
404  ++stepCounter;
405  }
406 
407  // Print out summary statistics
408  if (m_session->GetComm()->GetRank() == 0)
409  {
410  if (m_cflSafetyFactor > 0.0)
411  {
412  cout << "CFL safety factor : " << m_cflSafetyFactor << endl
413  << "CFL time-step : " << m_timestep << endl;
414  }
415 
416  if (m_session->GetSolverInfo("Driver") != "SteadyState")
417  {
418  cout << "Time-integration : " << intTime << "s" << endl;
419  }
420  }
421 
422  // If homogeneous, transform back into physical space if necessary.
424  {
425  for(i = 0; i < nfields; i++)
426  {
427  if (m_fields[i]->GetWaveSpace())
428  {
429  m_fields[i]->SetWaveSpace(false);
430  m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
431  m_fields[i]->UpdatePhys());
432  m_fields[i]->SetPhysState(true);
433  }
434  }
435  }
436  else
437  {
438  for(i = 0; i < nvariables; ++i)
439  {
440  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
441  m_fields[m_intVariables[i]]->SetPhysState(true);
442  }
443  }
444 
445  // Finalise filters
446  for (x = m_filters.begin(); x != m_filters.end(); ++x)
447  {
448  (*x)->Finalise(m_fields, m_time);
449  }
450 
451  // Print for 1D problems
452  if(m_spacedim == 1)
453  {
454  v_AppendOutput1D(fields);
455  }
456  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate(int step)
bool m_homoInitialFwd
Flag to determine if simulation should start in homogeneous forward transformed state.
NekDouble m_time
Current time of simulation.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
NekDouble m_timestep
Time step size.
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
virtual SOLVER_UTILS_EXPORT void v_AppendOutput1D(Array< OneD, Array< OneD, NekDouble > > &solution1D)
Print the solution at each solution point in a txt file.
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans()
NekDouble m_checktime
Time between checkpoints.
int m_nchk
Number of checkpoints written so far.
int m_checksteps
Number of steps between checkpoints.
NekDouble m_fintime
Finish time of the simulation.
int m_steps
Number of steps to take.
static const NekDouble kNekZeroTol
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate(int step)
int Nnan(int n, const T *x, const int incx)
Return number of NaN elements of x.
Definition: Vmath.cpp:892
int m_spacedim
Spatial dimension (>= expansion dim).
double NekDouble
SOLVER_UTILS_EXPORT NekDouble GetTimeStep()
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
StandardMatrixTag boost::call_traits< LhsDataType >::const_reference rhs typedef NekMatrix< LhsDataType, StandardMatrixTag >::iterator iterator
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
Wrapper to the time integration scheme.
int m_initialStep
Number of the step where the simulation should begin.
std::vector< FilterSharedPtr > m_filters
int m_infosteps
Number of time steps between outputting status information.
LibUtilities::TimeIntegrationSolutionSharedPtr m_intSoln
enum HomogeneousType m_HomogeneousType
void Nektar::SolverUtils::UnsteadySystem::v_GenerateSummary ( SummaryList s)
protectedvirtual

Print a summary of time stepping parameters.

Prints a summary with some information regards the time-stepping.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::CoupledLinearNS, Nektar::UnsteadyViscousBurgers, Nektar::UnsteadyAdvectionDiffusion, Nektar::VelocityCorrectionScheme, Nektar::CFLtester, Nektar::UnsteadyAdvection, Nektar::NonlinearPeregrine, Nektar::PulseWavePropagation, Nektar::Monodomain, Nektar::Bidomain, Nektar::BidomainRoth, Nektar::ShallowWaterSystem, Nektar::LinearSWE, Nektar::NonlinearSWE, Nektar::ImageWarpingSystem, Nektar::UnsteadyDiffusion, Nektar::VCSWeakPressure, Nektar::IsentropicVortex, and Nektar::RinglebFlow.

Definition at line 472 of file UnsteadySystem.cpp.

References Nektar::SolverUtils::AddSummaryItem(), Nektar::SolverUtils::EquationSystem::m_checksteps, m_explicitAdvection, m_explicitDiffusion, m_explicitReaction, m_intScheme, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_steps, Nektar::SolverUtils::EquationSystem::m_timestep, Nektar::LibUtilities::TimeIntegrationMethodMap, and Nektar::SolverUtils::EquationSystem::v_GenerateSummary().

Referenced by Nektar::IsentropicVortex::v_GenerateSummary(), Nektar::RinglebFlow::v_GenerateSummary(), Nektar::VCSWeakPressure::v_GenerateSummary(), Nektar::UnsteadyDiffusion::v_GenerateSummary(), Nektar::ImageWarpingSystem::v_GenerateSummary(), Nektar::ShallowWaterSystem::v_GenerateSummary(), Nektar::BidomainRoth::v_GenerateSummary(), Nektar::Bidomain::v_GenerateSummary(), Nektar::Monodomain::v_GenerateSummary(), Nektar::PulseWavePropagation::v_GenerateSummary(), Nektar::UnsteadyAdvection::v_GenerateSummary(), Nektar::CFLtester::v_GenerateSummary(), Nektar::VelocityCorrectionScheme::v_GenerateSummary(), Nektar::UnsteadyAdvectionDiffusion::v_GenerateSummary(), and Nektar::UnsteadyViscousBurgers::v_GenerateSummary().

473  {
475  AddSummaryItem(s, "Advection",
476  m_explicitAdvection ? "explicit" : "implicit");
477 
478  if(m_session->DefinesSolverInfo("AdvectionType"))
479  {
480  AddSummaryItem(s, "AdvectionType",
481  m_session->GetSolverInfo("AdvectionType"));
482  }
483 
484  AddSummaryItem(s, "Diffusion",
485  m_explicitDiffusion ? "explicit" : "implicit");
486 
487  if (m_session->GetSolverInfo("EQTYPE")
488  == "SteadyAdvectionDiffusionReaction")
489  {
490  AddSummaryItem(s, "Reaction",
491  m_explicitReaction ? "explicit" : "implicit");
492  }
493 
494  AddSummaryItem(s, "Time Step", m_timestep);
495  AddSummaryItem(s, "No. of Steps", m_steps);
496  AddSummaryItem(s, "Checkpoints (steps)", m_checksteps);
497  AddSummaryItem(s, "Integration Type",
499  m_intScheme->GetIntegrationMethod()]);
500  }
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &l)
Virtual function for generating summary information.
bool m_explicitReaction
Indicates if explicit or implicit treatment of reaction is used.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
NekDouble m_timestep
Time step size.
const char *const TimeIntegrationMethodMap[]
int m_checksteps
Number of steps between checkpoints.
int m_steps
Number of steps to take.
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
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
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
Wrapper to the time integration scheme.
NekDouble Nektar::SolverUtils::UnsteadySystem::v_GetTimeStep ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray)
protectedvirtual

Return the timestep to be used for the next step in the time-marching loop.

See also
UnsteadySystem::GetTimeStep

Reimplemented in Nektar::CompressibleFlowSystem.

Definition at line 931 of file UnsteadySystem.cpp.

References ASSERTL0.

Referenced by GetTimeStep().

933  {
934  ASSERTL0(false, "Not defined for this class");
935  return 0.0;
936  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::SolverUtils::UnsteadySystem::v_InitObject ( )
protectedvirtual

Init object for UnsteadySystem class.

Initialization object for UnsteadySystem class.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::CoupledLinearNS, Nektar::IncNavierStokes, Nektar::PulseWaveSystem, Nektar::UnsteadyViscousBurgers, Nektar::UnsteadyAdvectionDiffusion, Nektar::CFLtester, Nektar::CompressibleFlowSystem, Nektar::UnsteadyAdvection, Nektar::UnsteadyInviscidBurger, Nektar::APE, Nektar::ShallowWaterSystem, Nektar::NonlinearPeregrine, Nektar::PulseWavePropagation, Nektar::ImageWarpingSystem, Nektar::UnsteadyDiffusion, Nektar::Monodomain, Nektar::Bidomain, Nektar::BidomainRoth, Nektar::VCSMapping, Nektar::LinearSWE, Nektar::NonlinearSWE, Nektar::PulseWaveSystemOutput, Nektar::NavierStokesCFE, Nektar::VelocityCorrectionScheme, and Nektar::SolverUtils::AdvectionSystem.

Definition at line 78 of file UnsteadySystem.cpp.

References Nektar::SolverUtils::GetFilterFactory(), Nektar::LibUtilities::GetTimeIntegrationWrapperFactory(), m_cflSafetyFactor, m_explicitAdvection, m_explicitDiffusion, m_explicitReaction, Nektar::SolverUtils::EquationSystem::m_fieldMetaDataMap, m_filters, m_homoInitialFwd, m_infosteps, Nektar::SolverUtils::EquationSystem::m_initialStep, m_intScheme, m_nanSteps, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_time, and Nektar::SolverUtils::EquationSystem::v_InitObject().

Referenced by Nektar::SolverUtils::AdvectionSystem::v_InitObject(), Nektar::BidomainRoth::v_InitObject(), Nektar::Bidomain::v_InitObject(), Nektar::Monodomain::v_InitObject(), Nektar::UnsteadyDiffusion::v_InitObject(), Nektar::ImageWarpingSystem::v_InitObject(), Nektar::ShallowWaterSystem::v_InitObject(), Nektar::APE::v_InitObject(), Nektar::CompressibleFlowSystem::v_InitObject(), Nektar::CFLtester::v_InitObject(), and Nektar::PulseWaveSystem::v_InitObject().

79  {
81 
82  m_initialStep = 0;
83 
84  // Load SolverInfo parameters
85  m_session->MatchSolverInfo("DIFFUSIONADVANCEMENT","Explicit",
86  m_explicitDiffusion,true);
87  m_session->MatchSolverInfo("ADVECTIONADVANCEMENT","Explicit",
88  m_explicitAdvection,true);
89  m_session->MatchSolverInfo("REACTIONADVANCEMENT", "Explicit",
90  m_explicitReaction, true);
91 
92  m_session->LoadParameter("CheckNanSteps", m_nanSteps, 1);
93 
94  // For steady problems, we do not initialise the time integration
95  if (m_session->DefinesSolverInfo("TIMEINTEGRATIONMETHOD"))
96  {
98  CreateInstance(m_session->GetSolverInfo(
99  "TIMEINTEGRATIONMETHOD"));
100 
101  // Load generic input parameters
102  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
103  m_session->LoadParameter("CFL", m_cflSafetyFactor, 0.0);
104 
105  // Set up time to be dumped in field information
106  m_fieldMetaDataMap["Time"] =
107  boost::lexical_cast<std::string>(m_time);
108  }
109 
110  // By default attempt to forward transform initial condition.
111  m_homoInitialFwd = true;
112 
113  // Set up filters
114  LibUtilities::FilterMap::const_iterator x;
115  LibUtilities::FilterMap f = m_session->GetFilters();
116  for (x = f.begin(); x != f.end(); ++x)
117  {
118  m_filters.push_back(GetFilterFactory().CreateInstance(
119  x->first,
120  m_session,
121  x->second));
122  }
123  }
bool m_explicitReaction
Indicates if explicit or implicit treatment of reaction is used.
bool m_homoInitialFwd
Flag to determine if simulation should start in homogeneous forward transformed state.
NekDouble m_time
Current time of simulation.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
virtual SOLVER_UTILS_EXPORT void v_InitObject()
Initialisation object for EquationSystem.
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
std::vector< std::pair< std::string, FilterParams > > FilterMap
Definition: SessionReader.h:66
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
Map to identify relevant solver info to dump in output fields.
TimeIntegrationWrapperFactory & GetTimeIntegrationWrapperFactory()
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
Wrapper to the time integration scheme.
int m_initialStep
Number of the step where the simulation should begin.
FilterFactory & GetFilterFactory()
Definition: Filter.cpp:42
std::vector< FilterSharedPtr > m_filters
int m_infosteps
Number of time steps between outputting status information.
void Nektar::SolverUtils::UnsteadySystem::v_NumericalFlux ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numflux 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::IncNavierStokes, Nektar::CFLtester, Nektar::PulseWavePropagation, and Nektar::ImageWarpingSystem.

Definition at line 527 of file UnsteadySystem.cpp.

References ASSERTL0.

530  {
531  ASSERTL0(false,
532  "This function is not implemented for this equation.");
533  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::SolverUtils::UnsteadySystem::v_NumericalFlux ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxX,
Array< OneD, Array< OneD, NekDouble > > &  numfluxY 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 535 of file UnsteadySystem.cpp.

References ASSERTL0.

539  {
540  ASSERTL0(false,
541  "This function is not implemented for this equation.");
542  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:198
void Nektar::SolverUtils::UnsteadySystem::v_NumFluxforScalar ( const Array< OneD, Array< OneD, NekDouble > > &  ufield,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  uflux 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 544 of file UnsteadySystem.cpp.

References Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_traceNormals, Vmath::Vmul(), and WeakPenaltyforScalar().

547  {
548  int i, j;
549  int nTraceNumPoints = GetTraceNpoints();
550  int nvariables = m_fields.num_elements();
551  int nqvar = uflux.num_elements();
552 
553  Array<OneD, NekDouble > Fwd (nTraceNumPoints);
554  Array<OneD, NekDouble > Bwd (nTraceNumPoints);
555  Array<OneD, NekDouble > Vn (nTraceNumPoints, 0.0);
556  Array<OneD, NekDouble > fluxtemp(nTraceNumPoints, 0.0);
557 
558  // Get the sign of (v \cdot n), v = an arbitrary vector
559 
560  // Evaulate upwind flux:
561  // uflux = \hat{u} \phi \cdot u = u^{(+,-)} n
562  for (j = 0; j < nqvar; ++j)
563  {
564  for (i = 0; i < nvariables ; ++i)
565  {
566  // Compute Fwd and Bwd value of ufield of i direction
567  m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
568 
569  // if Vn >= 0, flux = uFwd, i.e.,
570  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick uflux = uFwd
571  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick uflux = uFwd
572 
573  // else if Vn < 0, flux = uBwd, i.e.,
574  // edge::eForward, if V*n<0 <=> V*n_F<0, pick uflux = uBwd
575  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick uflux = uBwd
576 
577  m_fields[i]->GetTrace()->Upwind(m_traceNormals[j],
578  Fwd, Bwd, fluxtemp);
579 
580  // Imposing weak boundary condition with flux
581  // if Vn >= 0, uflux = uBwd at Neumann, i.e.,
582  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick uflux = uBwd
583  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick uflux = uBwd
584 
585  // if Vn >= 0, uflux = uFwd at Neumann, i.e.,
586  // edge::eForward, if V*n<0 <=> V*n_F<0, pick uflux = uFwd
587  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick uflux = uFwd
588 
589  if(m_fields[0]->GetBndCondExpansions().num_elements())
590  {
591  WeakPenaltyforScalar(i, ufield[i], fluxtemp);
592  }
593 
594  // if Vn >= 0, flux = uFwd*(tan_{\xi}^- \cdot \vec{n}),
595  // i.e,
596  // edge::eForward, uFwd \(\tan_{\xi}^Fwd \cdot \vec{n})
597  // edge::eBackward, uFwd \(\tan_{\xi}^Bwd \cdot \vec{n})
598 
599  // else if Vn < 0, flux = uBwd*(tan_{\xi}^- \cdot \vec{n}),
600  // i.e,
601  // edge::eForward, uBwd \(\tan_{\xi}^Fwd \cdot \vec{n})
602  // edge::eBackward, uBwd \(\tan_{\xi}^Bwd \cdot \vec{n})
603 
604  Vmath::Vmul(nTraceNumPoints,
605  m_traceNormals[j], 1,
606  fluxtemp, 1,
607  uflux[j][i], 1);
608  }
609  }
610  }
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
void WeakPenaltyforScalar(const int var, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble time=0.0)
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:183
void Nektar::SolverUtils::UnsteadySystem::v_NumFluxforVector ( const Array< OneD, Array< OneD, NekDouble > > &  ufield,
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > &  qfield,
Array< OneD, Array< OneD, NekDouble > > &  qflux 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 614 of file UnsteadySystem.cpp.

References Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_traceNormals, Vmath::Smul(), Vmath::Vadd(), Vmath::Vmul(), Vmath::Vsub(), and WeakPenaltyforVector().

618  {
619  int nTraceNumPoints = GetTraceNpoints();
620  int nvariables = m_fields.num_elements();
621  int nqvar = qfield.num_elements();
622 
623  NekDouble C11 = 1.0;
624  Array<OneD, NekDouble > Fwd(nTraceNumPoints);
625  Array<OneD, NekDouble > Bwd(nTraceNumPoints);
626  Array<OneD, NekDouble > Vn (nTraceNumPoints, 0.0);
627 
628  Array<OneD, NekDouble > qFwd (nTraceNumPoints);
629  Array<OneD, NekDouble > qBwd (nTraceNumPoints);
630  Array<OneD, NekDouble > qfluxtemp(nTraceNumPoints, 0.0);
631 
632  Array<OneD, NekDouble > uterm(nTraceNumPoints);
633 
634  // Evaulate upwind flux:
635  // qflux = \hat{q} \cdot u = q \cdot n - C_(11)*(u^+ - u^-)
636  for (int i = 0; i < nvariables; ++i)
637  {
638  qflux[i] = Array<OneD, NekDouble> (nTraceNumPoints, 0.0);
639  for (int j = 0; j < nqvar; ++j)
640  {
641  // Compute Fwd and Bwd value of ufield of jth direction
642  m_fields[i]->GetFwdBwdTracePhys(qfield[j][i],qFwd,qBwd);
643 
644  // if Vn >= 0, flux = uFwd, i.e.,
645  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick
646  // qflux = qBwd = q+
647  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick
648  // qflux = qBwd = q-
649 
650  // else if Vn < 0, flux = uBwd, i.e.,
651  // edge::eForward, if V*n<0 <=> V*n_F<0, pick
652  // qflux = qFwd = q-
653  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick
654  // qflux = qFwd = q+
655 
656  m_fields[i]->GetTrace()->Upwind(m_traceNormals[j],
657  qBwd, qFwd,
658  qfluxtemp);
659 
660  Vmath::Vmul(nTraceNumPoints,
661  m_traceNormals[j], 1,
662  qfluxtemp, 1,
663  qfluxtemp, 1);
664 
665  // Generate Stability term = - C11 ( u- - u+ )
666  m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
667 
668  Vmath::Vsub(nTraceNumPoints,
669  Fwd, 1, Bwd, 1,
670  uterm, 1);
671 
672  Vmath::Smul(nTraceNumPoints,
673  -1.0 * C11, uterm, 1,
674  uterm, 1);
675 
676  // Flux = {Fwd, Bwd} * (nx, ny, nz) + uterm * (nx, ny)
677  Vmath::Vadd(nTraceNumPoints,
678  uterm, 1,
679  qfluxtemp, 1,
680  qfluxtemp, 1);
681 
682  // Imposing weak boundary condition with flux
683  if (m_fields[0]->GetBndCondExpansions().num_elements())
684  {
685  WeakPenaltyforVector(i, j,
686  qfield[j][i],
687  qfluxtemp,
688  C11);
689  }
690 
691  // q_hat \cdot n = (q_xi \cdot n_xi) or (q_eta \cdot n_eta)
692  // n_xi = n_x * tan_xi_x + n_y * tan_xi_y + n_z * tan_xi_z
693  // n_xi = n_x * tan_eta_x + n_y * tan_eta_y + n_z*tan_eta_z
694  Vmath::Vadd(nTraceNumPoints,
695  qfluxtemp, 1,
696  qflux[i], 1,
697  qflux[i], 1);
698  }
699  }
700  }
void WeakPenaltyforVector(const int var, const int dir, const Array< OneD, const NekDouble > &physfield, Array< OneD, NekDouble > &penaltyflux, NekDouble C11, NekDouble time=0.0)
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
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:213
double NekDouble
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:343
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
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:299
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:183
bool Nektar::SolverUtils::UnsteadySystem::v_PostIntegrate ( int  step)
protectedvirtual

Reimplemented in Nektar::IncNavierStokes, Nektar::CompressibleFlowSystem, and Nektar::APE.

Definition at line 943 of file UnsteadySystem.cpp.

Referenced by v_DoSolve(), and Nektar::APE::v_PostIntegrate().

944  {
945  return false;
946  }
bool Nektar::SolverUtils::UnsteadySystem::v_PreIntegrate ( int  step)
protectedvirtual

Reimplemented in Nektar::IncNavierStokes, Nektar::UnsteadyAdvectionDiffusion, and Nektar::APE.

Definition at line 938 of file UnsteadySystem.cpp.

Referenced by v_DoSolve(), and Nektar::APE::v_PreIntegrate().

939  {
940  return false;
941  }
virtual SOLVER_UTILS_EXPORT bool Nektar::SolverUtils::UnsteadySystem::v_RequireFwdTrans ( )
inlineprotectedvirtual

Reimplemented in Nektar::VelocityCorrectionScheme.

Definition at line 140 of file UnsteadySystem.h.

Referenced by v_DoSolve().

141  {
142  return true;
143  }
bool Nektar::SolverUtils::UnsteadySystem::v_SteadyStateCheck ( int  step)
protectedvirtual

Definition at line 948 of file UnsteadySystem.cpp.

949  {
950  return false;
951  }
void Nektar::SolverUtils::UnsteadySystem::WeakPenaltyforScalar ( const int  var,
const Array< OneD, const NekDouble > &  physfield,
Array< OneD, NekDouble > &  penaltyflux,
NekDouble  time = 0.0 
)
private

Definition at line 760 of file UnsteadySystem.cpp.

References Nektar::SpatialDomains::eDirichlet, Nektar::SpatialDomains::eNeumann, Nektar::SolverUtils::EquationSystem::GetPhys_Offset(), Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_session, and Vmath::Vcopy().

Referenced by v_NumFluxforScalar().

765  {
766  int i, e, npoints, id1, id2;
767 
768  // Number of boundary regions
769  int nbnd = m_fields[var]->GetBndCondExpansions().num_elements();
770  int Nfps, numBDEdge;
771  int nTraceNumPoints = GetTraceNpoints();
772  int cnt = 0;
773 
774  Array<OneD, NekDouble > uplus(nTraceNumPoints);
775 
776  m_fields[var]->ExtractTracePhys(physfield, uplus);
777  for (i = 0; i < nbnd; ++i)
778  {
779  // Number of boundary expansion related to that region
780  numBDEdge = m_fields[var]->
781  GetBndCondExpansions()[i]->GetExpSize();
782 
783  // Evaluate boundary values g_D or g_N from input files
785  m_session->GetFunction("InitialConditions", 0);
786 
787  npoints = m_fields[var]->
788  GetBndCondExpansions()[i]->GetNpoints();
789 
790  Array<OneD,NekDouble> BDphysics(npoints);
791  Array<OneD,NekDouble> x0(npoints,0.0);
792  Array<OneD,NekDouble> x1(npoints,0.0);
793  Array<OneD,NekDouble> x2(npoints,0.0);
794 
795  m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
796  ifunc->Evaluate(x0,x1,x2,time,BDphysics);
797 
798  // Weakly impose boundary conditions by modifying flux values
799  for (e = 0; e < numBDEdge ; ++e)
800  {
801  // Number of points on the expansion
802  Nfps = m_fields[var]->
803  GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
804 
805  id1 = m_fields[var]->
806  GetBndCondExpansions()[i]->GetPhys_Offset(e);
807 
808  id2 = m_fields[0]->GetTrace()->
809  GetPhys_Offset(m_fields[0]->GetTraceMap()->
810  GetBndCondTraceToGlobalTraceMap(cnt++));
811 
812  // For Dirichlet boundary condition: uflux = g_D
813  if (m_fields[var]->GetBndConditions()[i]->
814  GetBoundaryConditionType() == SpatialDomains::eDirichlet)
815  {
816  Vmath::Vcopy(Nfps,
817  &BDphysics[id1], 1,
818  &penaltyflux[id2], 1);
819  }
820  // For Neumann boundary condition: uflux = u+
821  else if ((m_fields[var]->GetBndConditions()[i])->
822  GetBoundaryConditionType() == SpatialDomains::eNeumann)
823  {
824  Vmath::Vcopy(Nfps,
825  &uplus[id2], 1,
826  &penaltyflux[id2], 1);
827  }
828  }
829  }
830  }
boost::shared_ptr< Equation > EquationSharedPtr
SOLVER_UTILS_EXPORT int GetPhys_Offset(int n)
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1061
void Nektar::SolverUtils::UnsteadySystem::WeakPenaltyforVector ( const int  var,
const int  dir,
const Array< OneD, const NekDouble > &  physfield,
Array< OneD, NekDouble > &  penaltyflux,
NekDouble  C11,
NekDouble  time = 0.0 
)
private

Diffusion: Imposing weak boundary condition for q with flux uflux = g_D on Dirichlet boundary condition uflux = u_Fwd on Neumann boundary condition

Definition at line 837 of file UnsteadySystem.cpp.

References Nektar::SpatialDomains::eDirichlet, Nektar::SpatialDomains::eNeumann, Nektar::SolverUtils::EquationSystem::GetPhys_Offset(), Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_traceNormals, and Vmath::Vmul().

Referenced by v_NumFluxforVector().

844  {
845  int i, e, npoints, id1, id2;
846  int nbnd = m_fields[var]->GetBndCondExpansions().num_elements();
847  int numBDEdge, Nfps;
848  int nTraceNumPoints = GetTraceNpoints();
849  Array<OneD, NekDouble > uterm(nTraceNumPoints);
850  Array<OneD, NekDouble > qtemp(nTraceNumPoints);
851  int cnt = 0;
852 
853  m_fields[var]->ExtractTracePhys(physfield,qtemp);
854 
855  for (i = 0; i < nbnd; ++i)
856  {
857  numBDEdge = m_fields[var]->
858  GetBndCondExpansions()[i]->GetExpSize();
859 
860  // Evaluate boundary values g_D or g_N from input files
862  m_session->GetFunction("InitialConditions", 0);
863 
864  npoints = m_fields[var]->
865  GetBndCondExpansions()[i]->GetNpoints();
866 
867  Array<OneD,NekDouble> BDphysics(npoints);
868  Array<OneD,NekDouble> x0(npoints,0.0);
869  Array<OneD,NekDouble> x1(npoints,0.0);
870  Array<OneD,NekDouble> x2(npoints,0.0);
871 
872  m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
873  ifunc->Evaluate(x0,x1,x2,time,BDphysics);
874 
875  // Weakly impose boundary conditions by modifying flux values
876  for (e = 0; e < numBDEdge ; ++e)
877  {
878  Nfps = m_fields[var]->
879  GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
880 
881  id1 = m_fields[var]->
882  GetBndCondExpansions()[i]->GetPhys_Offset(e);
883 
884  id2 = m_fields[0]->GetTrace()->
885  GetPhys_Offset(m_fields[0]->GetTraceMap()->
886  GetBndCondTraceToGlobalTraceMap(cnt++));
887 
888  // For Dirichlet boundary condition:
889  //qflux = q+ - C_11 (u+ - g_D) (nx, ny)
890  if(m_fields[var]->GetBndConditions()[i]->
891  GetBoundaryConditionType() == SpatialDomains::eDirichlet)
892  {
893  Vmath::Vmul(Nfps,
894  &m_traceNormals[dir][id2], 1,
895  &qtemp[id2], 1,
896  &penaltyflux[id2], 1);
897  }
898  // For Neumann boundary condition: qflux = g_N
899  else if((m_fields[var]->GetBndConditions()[i])->
900  GetBoundaryConditionType() == SpatialDomains::eNeumann)
901  {
902  Vmath::Vmul(Nfps,
903  &m_traceNormals[dir][id2], 1,
904  &BDphysics[id1], 1,
905  &penaltyflux[id2], 1);
906  }
907  }
908  }
909  }
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
boost::shared_ptr< Equation > EquationSharedPtr
SOLVER_UTILS_EXPORT int GetPhys_Offset(int n)
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
SOLVER_UTILS_EXPORT int GetTraceNpoints()
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:183

Member Data Documentation

NekDouble Nektar::SolverUtils::UnsteadySystem::m_cflSafetyFactor

CFL safety factor (comprise between 0 to 1).

Definition at line 59 of file UnsteadySystem.h.

Referenced by v_DoSolve(), Nektar::CompressibleFlowSystem::v_GetTimeStep(), and v_InitObject().

NekDouble Nektar::SolverUtils::UnsteadySystem::m_epsilon
protected

Definition at line 73 of file UnsteadySystem.h.

bool Nektar::SolverUtils::UnsteadySystem::m_explicitAdvection
protected
bool Nektar::SolverUtils::UnsteadySystem::m_explicitDiffusion
protected
bool Nektar::SolverUtils::UnsteadySystem::m_explicitReaction
protected

Indicates if explicit or implicit treatment of reaction is used.

Definition at line 79 of file UnsteadySystem.h.

Referenced by v_GenerateSummary(), and v_InitObject().

std::vector<FilterSharedPtr> Nektar::SolverUtils::UnsteadySystem::m_filters
protected
bool Nektar::SolverUtils::UnsteadySystem::m_homoInitialFwd
protected
int Nektar::SolverUtils::UnsteadySystem::m_infosteps
protected
LibUtilities::TimeIntegrationWrapperSharedPtr Nektar::SolverUtils::UnsteadySystem::m_intScheme
protected
LibUtilities::TimeIntegrationSolutionSharedPtr Nektar::SolverUtils::UnsteadySystem::m_intSoln
protected
std::vector<int> Nektar::SolverUtils::UnsteadySystem::m_intVariables
protected
int Nektar::SolverUtils::UnsteadySystem::m_nanSteps
protected

Definition at line 65 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

LibUtilities::TimeIntegrationSchemeOperators Nektar::SolverUtils::UnsteadySystem::m_ode
protected