Nektar++
Public Member Functions | Static 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:
[legend]

Public Member Functions

SOLVER_UTILS_EXPORT ~UnsteadySystem () override=default
 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...
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep ()
 
SOLVER_UTILS_EXPORT void SetTimeStep (const NekDouble timestep)
 
SOLVER_UTILS_EXPORT void SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
SOLVER_UTILS_EXPORT LibUtilities::TimeIntegrationSchemeSharedPtrGetTimeIntegrationScheme ()
 Returns the time integration scheme. More...
 
SOLVER_UTILS_EXPORT LibUtilities::TimeIntegrationSchemeOperatorsGetTimeIntegrationSchemeOperators ()
 Returns the time integration scheme operators. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::EquationSystem
virtual SOLVER_UTILS_EXPORT ~EquationSystem ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT void InitObject (bool DeclareField=true)
 Initialises the members of this object. More...
 
SOLVER_UTILS_EXPORT void DoInitialise (bool dumpInitialConditions=true)
 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 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 NekDouble LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Linf error computation. More...
 
SOLVER_UTILS_EXPORT NekDouble H1Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln, bool Normalised=false)
 Compute the H1 error between fields and a given exact solution. 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 GetTime ()
 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 void SetSteps (const int steps)
 
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, const Array< OneD, const NekDouble > &input)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDouble > & UpdatePhysField (const int i)
 
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 int GetInfoSteps ()
 
SOLVER_UTILS_EXPORT void SetInfoSteps (int num)
 
SOLVER_UTILS_EXPORT void SetIterationNumberPIT (int num)
 
SOLVER_UTILS_EXPORT void SetWindowNumberPIT (int num)
 
SOLVER_UTILS_EXPORT Array< OneD, const Array< OneD, NekDouble > > GetTraceNormals ()
 
SOLVER_UTILS_EXPORT void SetTime (const NekDouble time)
 
SOLVER_UTILS_EXPORT void SetTimeStep (const NekDouble timestep)
 
SOLVER_UTILS_EXPORT void SetInitialStep (const int step)
 
SOLVER_UTILS_EXPORT void SetBoundaryConditions (NekDouble time)
 Evaluates the boundary conditions at the given time. More...
 
SOLVER_UTILS_EXPORT bool NegatedOp ()
 Identify if operator is negated in DoSolve. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::ALEHelper
virtual ~ALEHelper ()=default
 
virtual SOLVER_UTILS_EXPORT void v_ALEInitObject (int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 
SOLVER_UTILS_EXPORT void InitObject (int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 
virtual SOLVER_UTILS_EXPORT void v_UpdateGridVelocity (const NekDouble &time)
 
virtual SOLVER_UTILS_EXPORT void v_ALEPreMultiplyMass (Array< OneD, Array< OneD, NekDouble > > &fields)
 
SOLVER_UTILS_EXPORT void ALEDoElmtInvMass (Array< OneD, Array< OneD, NekDouble > > &traceNormals, Array< OneD, Array< OneD, NekDouble > > &fields, NekDouble time)
 Update m_fields with u^n by multiplying by inverse mass matrix. That's then used in e.g. checkpoint output and L^2 error calculation. More...
 
SOLVER_UTILS_EXPORT void ALEDoElmtInvMassBwdTrans (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
SOLVER_UTILS_EXPORT void MoveMesh (const NekDouble &time, Array< OneD, Array< OneD, NekDouble > > &traceNormals)
 
const Array< OneD, const Array< OneD, NekDouble > > & GetGridVelocity ()
 
SOLVER_UTILS_EXPORT const Array< OneD, const Array< OneD, NekDouble > > & GetGridVelocityTrace ()
 
SOLVER_UTILS_EXPORT void ExtraFldOutputGridVelocity (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 

Static Public Attributes

static std::string cmdSetStartTime
 
static std::string cmdSetStartChkNum
 

Protected Member Functions

SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises UnsteadySystem class members. More...
 
SOLVER_UTILS_EXPORT void v_InitObject (bool DeclareField=true) override
 Init object for UnsteadySystem class. More...
 
SOLVER_UTILS_EXPORT void v_DoSolve () override
 Solves an unsteady problem. More...
 
virtual SOLVER_UTILS_EXPORT void v_PrintStatusInformation (const int step, const NekDouble cpuTime)
 Print Status Information. More...
 
virtual SOLVER_UTILS_EXPORT void v_PrintSummaryStatistics (const NekDouble intTime)
 Print Summary Statistics. More...
 
SOLVER_UTILS_EXPORT void v_DoInitialise (bool dumpInitialConditions=true) override
 Sets up initial conditions. More...
 
SOLVER_UTILS_EXPORT void v_GenerateSummary (SummaryList &s) override
 Print a summary of time stepping parameters. 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_PostIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans ()
 
virtual SOLVER_UTILS_EXPORT void v_SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
virtual SOLVER_UTILS_EXPORT bool v_UpdateTimeStepCheck ()
 
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator ()
 Get the maximum timestep estimator for cfl control. More...
 
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...
 
SOLVER_UTILS_EXPORT void DoDummyProjection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 Perform dummy projection. 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 void v_InitObject (bool DeclareFeld=true)
 Initialisation object for EquationSystem. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoInitialise (bool dumpInitialConditions=true)
 Virtual function for initialisation implementation. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoSolve ()
 Virtual function for solve implementation. 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 NekDouble v_H1Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the H_1 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_GenerateSummary (SummaryList &l)
 Virtual function for generating summary information. 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)
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure (void)
 
virtual SOLVER_UTILS_EXPORT bool v_NegatedOp (void)
 Virtual function to identify if operator is negated in DoSolve. More...
 
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 

Protected Attributes

LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
 Wrapper to the time integration scheme. More...
 
LibUtilities::TimeIntegrationSchemeOperators m_ode
 The time integration scheme operators to use. More...
 
Array< OneD, Array< OneD, NekDouble > > m_previousSolution
 Storage for previous solution for steady-state check. More...
 
std::vector< int > m_intVariables
 
NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 
NekDouble m_CFLGrowth
 CFL growth rate. More...
 
NekDouble m_CFLEnd
 Maximun cfl in cfl growth. More...
 
int m_abortSteps
 Number of steps between checks for abort conditions. More...
 
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...
 
int m_steadyStateSteps
 Check for steady state at step interval. More...
 
NekDouble m_steadyStateTol
 Tolerance to which steady state should be evaluated at. More...
 
int m_filtersInfosteps
 Number of time steps between outputting filters information. More...
 
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
 
bool m_homoInitialFwd
 Flag to determine if simulation should start in homogeneous forward transformed state. More...
 
std::ofstream m_errFile
 
NekDouble m_epsilon
 Diffusion coefficient. More...
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
bool m_verbose
 
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...
 
NekDouble m_lastCheckTime
 
NekDouble m_TimeIncrementFactor
 
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_infosteps
 Number of time steps between outputting status information. More...
 
int m_iterPIT = 0
 Number of parallel-in-time time iteration. More...
 
int m_windowPIT = 0
 Index of windows for parallel-in-time time iteration. 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...
 
Array< OneD, NekDoublem_movingFrameData
 Moving reference frame status in the inertial frame X, Y, Z, Theta_x, Theta_y, Theta_z, U, V, W, Omega_x, Omega_y, Omega_z, A_x, A_y, A_z, DOmega_x, DOmega_y, DOmega_z, pivot_x, pivot_y, pivot_z. More...
 
std::vector< std::string > m_strFrameData
 variable name in m_movingFrameData 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...
 
- Protected Attributes inherited from Nektar::SolverUtils::ALEHelper
Array< OneD, MultiRegions::ExpListSharedPtrm_fieldsALE
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocity
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
 
std::vector< ALEBaseShPtrm_ALEs
 
bool m_ALESolver = false
 
bool m_ImplicitALESolver = false
 
NekDouble m_prevStageTime = 0.0
 
int m_spaceDim
 

Private Member Functions

SOLVER_UTILS_EXPORT void AppendOutput1D (void)
 Print the solution at each solution point in a txt file. More...
 
void InitializeSteadyState ()
 
bool CheckSteadyState (int step, const NekDouble &totCPUTime=0.0)
 Calculate whether the system has reached a steady state by observing residuals to a user-defined tolerance. More...
 

Additional Inherited Members

- 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 []
 
static std::string projectionTypeLookupIds []
 

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 46 of file UnsteadySystem.h.

Constructor & Destructor Documentation

◆ ~UnsteadySystem()

SOLVER_UTILS_EXPORT Nektar::SolverUtils::UnsteadySystem::~UnsteadySystem ( )
overridedefault

Destructor.

◆ UnsteadySystem()

Nektar::SolverUtils::UnsteadySystem::UnsteadySystem ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)
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 74 of file UnsteadySystem.cpp.

77 : EquationSystem(pSession, pGraph), SolverUtils::ALEHelper()
78
79{
80}
SOLVER_UTILS_EXPORT EquationSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises EquationSystem class members.

Member Function Documentation

◆ AppendOutput1D()

void Nektar::SolverUtils::UnsteadySystem::AppendOutput1D ( void  )
private

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 668 of file UnsteadySystem.cpp.

669{
670 // Coordinates of the quadrature points in the real physical space.
671 Array<OneD, NekDouble> x(GetNpoints());
672 Array<OneD, NekDouble> y(GetNpoints());
673 Array<OneD, NekDouble> z(GetNpoints());
674 m_fields[0]->GetCoords(x, y, z);
675
676 // Print out the solution in a txt file.
677 ofstream outfile;
678 outfile.open("solution1D.txt");
679 for (int i = 0; i < GetNpoints(); i++)
680 {
681 outfile << scientific << setw(17) << setprecision(16) << x[i] << " "
682 << m_fields[m_intVariables[0]]->GetPhys()[i] << endl;
683 }
684 outfile << endl << endl;
685 outfile.close();
686}
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
SOLVER_UTILS_EXPORT int GetNpoints()
std::vector< double > z(NPUPPER)

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

Referenced by v_DoSolve().

◆ CheckForRestartTime()

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

Definition at line 691 of file UnsteadySystem.cpp.

692{
693 if (m_session->DefinesFunction("InitialConditions"))
694 {
695 for (int i = 0; i < m_fields.size(); ++i)
696 {
698
699 vType = m_session->GetFunctionType("InitialConditions",
700 m_session->GetVariable(i));
701
703 {
704 std::string filename = m_session->GetFunctionFilename(
705 "InitialConditions", m_session->GetVariable(i));
706
707 fs::path pfilename(filename);
708
709 // Redefine path for parallel file which is in directory.
710 if (fs::is_directory(pfilename))
711 {
712 fs::path metafile("Info.xml");
713 fs::path fullpath = pfilename / metafile;
714 filename = LibUtilities::PortablePath(fullpath);
715 }
718 fld->ImportFieldMetaData(filename, m_fieldMetaDataMap);
719
720 // Check to see if time defined.
722 {
723 auto iter = m_fieldMetaDataMap.find("Time");
724 if (iter != m_fieldMetaDataMap.end())
725 {
726 time = std::stod(iter->second);
727 }
728
729 iter = m_fieldMetaDataMap.find("ChkFileNum");
730 if (iter != m_fieldMetaDataMap.end())
731 {
732 nchk = std::stod(iter->second);
733 }
734 }
735
736 break;
737 }
738 }
739 }
740 if (m_session->DefinesCmdLineArgument("set-start-time"))
741 {
742 time = std::stod(
743 m_session->GetCmdLineArgument<std::string>("set-start-time")
744 .c_str());
745 }
746 if (m_session->DefinesCmdLineArgument("set-start-chknumber"))
747 {
748 nchk = std::stoi(
749 m_session->GetCmdLineArgument<std::string>("set-start-chknumber"));
750 }
751 ASSERTL0(time >= 0 && nchk >= 0,
752 "Starting time and checkpoint number should be >= 0");
753}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
static std::shared_ptr< FieldIO > CreateForFile(const LibUtilities::SessionReaderSharedPtr session, const std::string &filename)
Construct a FieldIO object for a given input filename.
Definition: FieldIO.cpp:223
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
Map to identify relevant solver info to dump in output fields.
std::shared_ptr< FieldIO > FieldIOSharedPtr
Definition: FieldIO.h:322
static std::string PortablePath(const fs::path &path)
create portable path on different platforms for std::filesystem path.
Definition: Filesystem.hpp:66
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:51

References ASSERTL0, Nektar::LibUtilities::FieldIO::CreateForFile(), Nektar::LibUtilities::eFunctionTypeFile, 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().

◆ CheckSteadyState()

bool Nektar::SolverUtils::UnsteadySystem::CheckSteadyState ( int  step,
const NekDouble totCPUTime = 0.0 
)
private

Calculate whether the system has reached a steady state by observing residuals to a user-defined tolerance.

Definition at line 886 of file UnsteadySystem.cpp.

887{
888 const int nPoints = GetTotPoints();
889 const int nFields = m_fields.size();
890
891 // Holds L2 errors.
892 Array<OneD, NekDouble> L2(nFields);
893
894 SteadyStateResidual(step, L2);
895
896 if (m_infosteps && m_comm->GetRank() == 0 &&
897 (((step + 1) % m_infosteps == 0) || ((step == m_initialStep))))
898 {
899 // Output time.
900 m_errFile << boost::format("%25.19e") % m_time;
901
902 m_errFile << " " << boost::format("%25.19e") % totCPUTime;
903
904 int stepp = step + 1;
905
906 m_errFile << " " << boost::format("%25.19e") % stepp;
907
908 // Output residuals.
909 for (int i = 0; i < nFields; ++i)
910 {
911 m_errFile << " " << boost::format("%25.19e") % L2[i];
912 }
913
914 m_errFile << endl;
915 }
916
917 // Calculate maximum L2 error.
918 NekDouble maxL2 = Vmath::Vmax(nFields, L2, 1);
919
920 if (m_infosteps && m_session->DefinesCmdLineArgument("verbose") &&
921 m_comm->GetRank() == 0 && ((step + 1) % m_infosteps == 0))
922 {
923 cout << "-- Maximum L^2 residual: " << maxL2 << endl;
924 }
925
926 if (maxL2 <= m_steadyStateTol)
927 {
928 return true;
929 }
930
931 for (int i = 0; i < m_fields.size(); ++i)
932 {
933 Vmath::Vcopy(nPoints, m_fields[i]->GetPhys(), 1, m_previousSolution[i],
934 1);
935 }
936
937 return false;
938}
LibUtilities::CommSharedPtr m_comm
Communicator.
int m_infosteps
Number of time steps between outputting status information.
NekDouble m_time
Current time of simulation.
int m_initialStep
Number of the step where the simulation should begin.
SOLVER_UTILS_EXPORT int GetTotPoints()
Array< OneD, Array< OneD, NekDouble > > m_previousSolution
Storage for previous solution for steady-state check.
NekDouble m_steadyStateTol
Tolerance to which steady state should be evaluated at.
SOLVER_UTILS_EXPORT void SteadyStateResidual(int step, Array< OneD, NekDouble > &L2)
double NekDouble
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.hpp:644
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References CellMLToNektar.pycml::format, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_comm, m_errFile, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_infosteps, Nektar::SolverUtils::EquationSystem::m_initialStep, m_previousSolution, Nektar::SolverUtils::EquationSystem::m_session, m_steadyStateTol, Nektar::SolverUtils::EquationSystem::m_time, SteadyStateResidual(), Vmath::Vcopy(), and Vmath::Vmax().

Referenced by v_DoSolve().

◆ DoDummyProjection()

void Nektar::SolverUtils::UnsteadySystem::DoDummyProjection ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble  time 
)
protected

Perform dummy projection.

Definition at line 982 of file UnsteadySystem.cpp.

986{
987
988 if (&inarray != &outarray)
989 {
990 for (int i = 0; i < inarray.size(); ++i)
991 {
992 Vmath::Vcopy(GetNpoints(), inarray[i], 1, outarray[i], 1);
993 }
994 }
995}

References Nektar::SolverUtils::EquationSystem::GetNpoints(), and Vmath::Vcopy().

Referenced by Nektar::Bidomain::v_InitObject(), Nektar::BidomainRoth::v_InitObject(), and Nektar::Monodomain::v_InitObject().

◆ GetTimeIntegrationScheme()

LibUtilities::TimeIntegrationSchemeSharedPtr & Nektar::SolverUtils::UnsteadySystem::GetTimeIntegrationScheme ( )

Returns the time integration scheme.

Definition at line 186 of file UnsteadySystem.cpp.

188{
189 return m_intScheme;
190}
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
Wrapper to the time integration scheme.

References m_intScheme.

◆ GetTimeIntegrationSchemeOperators()

LibUtilities::TimeIntegrationSchemeOperators & Nektar::SolverUtils::UnsteadySystem::GetTimeIntegrationSchemeOperators ( )

Returns the time integration scheme operators.

Definition at line 195 of file UnsteadySystem.cpp.

197{
198 return m_ode;
199}
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.

References m_ode.

◆ GetTimeStep() [1/2]

SOLVER_UTILS_EXPORT NekDouble Nektar::SolverUtils::UnsteadySystem::GetTimeStep ( )
inline

Definition at line 59 of file UnsteadySystem.h.

60 {
62 }
SOLVER_UTILS_EXPORT NekDouble GetTimeStep()

References Nektar::SolverUtils::EquationSystem::GetTimeStep().

Referenced by v_DoSolve().

◆ GetTimeStep() [2/2]

SOLVER_UTILS_EXPORT NekDouble Nektar::SolverUtils::UnsteadySystem::GetTimeStep ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray)
inline

Calculate the larger time-step mantaining the problem stable.

Definition at line 54 of file UnsteadySystem.h.

55 {
56 return v_GetTimeStep(inarray);
57 }
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.

References v_GetTimeStep().

◆ InitializeSteadyState()

void Nektar::SolverUtils::UnsteadySystem::InitializeSteadyState ( )
private

Definition at line 846 of file UnsteadySystem.cpp.

847{
848 if (m_steadyStateTol > 0.0)
849 {
850 const int nPoints = m_fields[0]->GetTotPoints();
852 Array<OneD, Array<OneD, NekDouble>>(m_fields.size());
853
854 for (int i = 0; i < m_fields.size(); ++i)
855 {
856 m_previousSolution[i] = Array<OneD, NekDouble>(nPoints);
857 Vmath::Vcopy(nPoints, m_fields[i]->GetPhys(), 1,
858 m_previousSolution[i], 1);
859 }
860
861 if (m_comm->GetRank() == 0)
862 {
863 std::string fName =
864 m_session->GetSessionName() + std::string(".resd");
865 m_errFile.open(fName.c_str());
866 m_errFile << setw(26) << left << "# Time";
867
868 m_errFile << setw(26) << left << "CPU_Time";
869
870 m_errFile << setw(26) << left << "Step";
871
872 for (int i = 0; i < m_fields.size(); ++i)
873 {
874 m_errFile << setw(26) << m_session->GetVariables()[i];
875 }
876
877 m_errFile << endl;
878 }
879 }
880}

References Nektar::SolverUtils::EquationSystem::m_comm, m_errFile, Nektar::SolverUtils::EquationSystem::m_fields, m_previousSolution, Nektar::SolverUtils::EquationSystem::m_session, m_steadyStateTol, and Vmath::Vcopy().

Referenced by v_DoInitialise().

◆ MaxTimeStepEstimator()

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 178 of file UnsteadySystem.cpp.

179{
180 return m_intScheme->GetTimeStability();
181}

References m_intScheme.

Referenced by Nektar::CompressibleFlowSystem::GetElmtTimeStep().

◆ SetTimeStep()

SOLVER_UTILS_EXPORT void Nektar::SolverUtils::UnsteadySystem::SetTimeStep ( const NekDouble  timestep)
inline

Definition at line 64 of file UnsteadySystem.h.

65 {
67 }
SOLVER_UTILS_EXPORT void SetTimeStep(const NekDouble timestep)

References Nektar::SolverUtils::EquationSystem::SetTimeStep().

◆ SteadyStateResidual()

SOLVER_UTILS_EXPORT void Nektar::SolverUtils::UnsteadySystem::SteadyStateResidual ( int  step,
Array< OneD, NekDouble > &  L2 
)
inline

Definition at line 69 of file UnsteadySystem.h.

71 {
72 v_SteadyStateResidual(step, L2);
73 }
virtual SOLVER_UTILS_EXPORT void v_SteadyStateResidual(int step, Array< OneD, NekDouble > &L2)

References v_SteadyStateResidual().

Referenced by CheckSteadyState().

◆ SVVVarDiffCoeff()

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 803 of file UnsteadySystem.cpp.

806{
807 int phystot = m_fields[0]->GetTotPoints();
808 int nvel = vel.size();
809
810 Array<OneD, NekDouble> varcoeff(phystot), tmp;
811
812 // Calculate magnitude of v.
813 Vmath::Vmul(phystot, vel[0], 1, vel[0], 1, varcoeff, 1);
814 for (int n = 1; n < nvel; ++n)
815 {
816 Vmath::Vvtvp(phystot, vel[n], 1, vel[n], 1, varcoeff, 1, varcoeff, 1);
817 }
818 Vmath::Vsqrt(phystot, varcoeff, 1, varcoeff, 1);
819
820 for (int i = 0; i < m_fields[0]->GetNumElmts(); ++i)
821 {
822 int offset = m_fields[0]->GetPhys_Offset(i);
823 int nq = m_fields[0]->GetExp(i)->GetTotPoints();
824 Array<OneD, NekDouble> unit(nq, 1.0);
825
826 int nmodes = 0;
827
828 for (int n = 0; n < m_fields[0]->GetExp(i)->GetNumBases(); ++n)
829 {
830 nmodes = max(nmodes, m_fields[0]->GetExp(i)->GetBasisNumModes(n));
831 }
832
833 NekDouble h = m_fields[0]->GetExp(i)->Integral(unit);
834 h = pow(h, (NekDouble)(1.0 / nvel)) / ((NekDouble)nmodes);
835
836 Vmath::Smul(nq, h, varcoeff + offset, 1, tmp = varcoeff + offset, 1);
837 }
838
839 // Set up map with eVarCoffLaplacian key.
840 varCoeffMap[StdRegions::eVarCoeffLaplacian] = varcoeff;
841}
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.hpp:340
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.hpp:72
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.hpp:366
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100

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

Referenced by Nektar::UnsteadyViscousBurgers::DoImplicitSolve().

◆ v_DoInitialise()

void Nektar::SolverUtils::UnsteadySystem::v_DoInitialise ( bool  dumpInitialConditions = true)
overrideprotectedvirtual

Sets up initial conditions.

Sets the initial conditions.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::VCSImplicit, Nektar::SolverUtils::FileSolution, Nektar::PulseWaveSystem, Nektar::MMFSWE, Nektar::CoupledLinearNS, Nektar::VCSMapping, and Nektar::VelocityCorrectionScheme.

Definition at line 624 of file UnsteadySystem.cpp.

625{
628 SetInitialConditions(m_time, dumpInitialConditions);
629
631
633}
virtual SOLVER_UTILS_EXPORT void v_UpdateGridVelocity(const NekDouble &time)
Definition: ALEHelper.cpp:90
SOLVER_UTILS_EXPORT void SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
Initialise the data in the dependent fields.
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 CheckForRestartTime(NekDouble &time, int &nchk)

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

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

◆ v_DoSolve()

void Nektar::SolverUtils::UnsteadySystem::v_DoSolve ( void  )
overrideprotectedvirtual

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::MMFAdvection, Nektar::CFSImplicit, Nektar::MMFMaxwell, Nektar::PulseWaveSystem, Nektar::MMFSWE, and Nektar::CoupledLinearNS.

Definition at line 205 of file UnsteadySystem.cpp.

206{
207 ASSERTL0(m_intScheme != nullptr, "No time integration scheme.");
208
209 int i = 1;
210 int nvariables = 0;
211 int nfields = m_fields.size();
212 if (m_intVariables.empty())
213 {
214 for (i = 0; i < nfields; ++i)
215 {
216 m_intVariables.push_back(i);
217 }
218 nvariables = nfields;
219 }
220 else
221 {
222 nvariables = m_intVariables.size();
223 }
224
225 // Integrate in wave-space if using homogeneous1D.
227 {
228 for (i = 0; i < nfields; ++i)
229 {
230 m_fields[i]->HomogeneousFwdTrans(m_fields[i]->GetTotPoints(),
231 m_fields[i]->GetPhys(),
232 m_fields[i]->UpdatePhys());
233 m_fields[i]->SetWaveSpace(true);
234 m_fields[i]->SetPhysState(false);
235 }
236 }
237
238 // Set up wrapper to fields data storage.
239 Array<OneD, Array<OneD, NekDouble>> fields(nvariables);
240
241 // Order storage to list time-integrated fields first.
242 // @TODO: Refactor to take coeffs (FwdTrans) if boolean flag (in constructor
243 // function) says to.
244 for (i = 0; i < nvariables; ++i)
245 {
246 fields[i] = m_fields[m_intVariables[i]]->UpdatePhys();
247 m_fields[m_intVariables[i]]->SetPhysState(false);
248 }
249
250 // @TODO: Virtual function that allows to transform the field space, embed
251 // the MultiplyMassMatrix in here.
252 // @TODO: Specify what the fields variables are physical or coefficient,
253 // boolean in UnsteadySystem class...
254
255 v_ALEPreMultiplyMass(fields);
256
257 // Initialise time integration scheme.
258 m_intScheme->InitializeScheme(m_timestep, fields, m_time, m_ode);
259
260 // Initialise filters.
261 for (auto &x : m_filters)
262 {
263 x.second->Initialise(m_fields, m_time);
264 }
265
266 LibUtilities::Timer timer;
267 bool doCheckTime = false;
268 int step = m_initialStep;
269 int stepCounter = 0;
270 NekDouble intTime = 0.0;
271 NekDouble cpuTime = 0.0;
272 NekDouble cpuPrevious = 0.0;
273 NekDouble elapsed = 0.0;
274 NekDouble totFilterTime = 0.0;
275
276 m_lastCheckTime = 0.0;
277
278 Array<OneD, int> abortFlags(2, 0);
279 string abortFile = "abort";
280 if (m_session->DefinesSolverInfo("CheckAbortFile"))
281 {
282 abortFile = m_session->GetSolverInfo("CheckAbortFile");
283 }
284
285 NekDouble tmp_cflSafetyFactor = m_cflSafetyFactor;
286 while ((step < m_steps || m_time < m_fintime - NekConstants::kNekZeroTol) &&
287 abortFlags[1] == 0)
288 {
290 {
291 tmp_cflSafetyFactor =
292 min(m_CFLEnd, m_CFLGrowth * tmp_cflSafetyFactor);
293 }
294
295 // Frozen preconditioner checks.
296 if (!m_comm->IsParallelInTime())
297 {
299 {
300 m_cflSafetyFactor = tmp_cflSafetyFactor;
301
303 {
304 m_timestep = GetTimeStep(fields);
305 }
306
307 // Ensure that the final timestep finishes at the final
308 // time, or at a prescribed IO_CheckTime.
309 if (m_time + m_timestep > m_fintime && m_fintime > 0.0)
310 {
312 }
313 else if (m_checktime &&
315 {
318 doCheckTime = true;
319 }
320 }
321 }
322
323 if (m_TimeIncrementFactor > 1.0)
324 {
325 NekDouble timeincrementFactor = m_TimeIncrementFactor;
326 m_timestep *= timeincrementFactor;
327
328 if (m_time + m_timestep > m_fintime && m_fintime > 0.0)
329 {
331 }
332 }
333
334 // Perform any solver-specific pre-integration steps.
335 timer.Start();
336 if (v_PreIntegrate(
337 step)) // Could be possible to put a preintegrate step in the
338 // ALEHelper class, put in the Unsteady Advection class
339 {
340 break;
341 }
342
343 ASSERTL0(m_timestep > 0, "m_timestep < 0");
344
345 fields = m_intScheme->TimeIntegrate(stepCounter, m_timestep);
346 timer.Stop();
347
349 elapsed = timer.TimePerTest(1);
350 intTime += elapsed;
351 cpuTime += elapsed;
352
353 // Write out status information.
354 v_PrintStatusInformation(step, cpuTime);
355 if (m_infosteps &&
356 m_session->GetComm()->GetSpaceComm()->GetRank() == 0 &&
357 !((step + 1) % m_infosteps))
358 {
359 cpuPrevious = cpuTime;
360 cpuTime = 0.0;
361 }
362
363 // @TODO: Another virtual function with this in it based on if ALE or
364 // not.
365 if (m_ALESolver) // Change to advect coeffs, change flag to physical vs
366 // coefficent space
367 {
370 }
371 else
372 {
373 // Transform data into coefficient space
374 for (i = 0; i < nvariables; ++i)
375 {
376 // copy fields into ExpList::m_phys and assign the new
377 // array to fields
378 m_fields[m_intVariables[i]]->SetPhys(fields[i]);
379 fields[i] = m_fields[m_intVariables[i]]->UpdatePhys();
380 if (v_RequireFwdTrans())
381 {
382 if (m_comm->IsParallelInTime())
383 {
384 m_fields[m_intVariables[i]]->FwdTrans(
385 m_fields[m_intVariables[i]]->GetPhys(),
386 m_fields[m_intVariables[i]]->UpdateCoeffs());
387 }
388 else
389 {
390 m_fields[m_intVariables[i]]->FwdTransLocalElmt(
391 m_fields[m_intVariables[i]]->GetPhys(),
392 m_fields[m_intVariables[i]]->UpdateCoeffs());
393 }
394 }
395 m_fields[m_intVariables[i]]->SetPhysState(false);
396 }
397 }
398
399 // Perform any solver-specific post-integration steps.
400 if (v_PostIntegrate(step))
401 {
402 break;
403 }
404
405 // Check for steady-state.
407 (!((step + 1) % m_steadyStateSteps)))
408 {
409 if (CheckSteadyState(step, intTime))
410 {
411 if (m_comm->GetRank() == 0)
412 {
413 cout << "Reached Steady State to tolerance "
414 << m_steadyStateTol << endl;
415 }
416 break;
417 }
418 }
419
420 // Test for abort conditions (nan, or abort file).
421 if (m_abortSteps && !((step + 1) % m_abortSteps))
422 {
423 abortFlags[0] = 0;
424 for (i = 0; i < nvariables; ++i)
425 {
426 if (Vmath::Nnan(m_fields[m_intVariables[i]]->GetPhys().size(),
427 m_fields[m_intVariables[i]]->GetPhys(), 1) > 0)
428 {
429 abortFlags[0] = 1;
430 }
431 }
432
433 // Rank zero looks for abort file and deletes it
434 // if it exists. Communicates the abort.
435 if (m_session->GetComm()->GetSpaceComm()->GetRank() == 0)
436 {
437 if (fs::exists(abortFile))
438 {
439 fs::remove(abortFile);
440 abortFlags[1] = 1;
441 }
442 }
443
444 m_session->GetComm()->GetSpaceComm()->AllReduce(
445 abortFlags, LibUtilities::ReduceMax);
446
447 ASSERTL0(!abortFlags[0], "NaN found during time integration.");
448 }
449
450 // Update filters.
451 for (auto &x : m_filters)
452 {
453 timer.Start();
454 x.second->Update(m_fields, m_time);
455 timer.Stop();
456 elapsed = timer.TimePerTest(1);
457 totFilterTime += elapsed;
458
459 // Write out individual filter status information.
460 if (m_filtersInfosteps && m_session->GetComm()->GetRank() == 0 &&
461 !((step + 1) % m_filtersInfosteps) && !m_filters.empty() &&
462 m_session->DefinesCmdLineArgument("verbose"))
463 {
464 stringstream s0;
465 s0 << x.first << ":";
466 stringstream s1;
467 s1 << elapsed << "s";
468 stringstream s2;
469 s2 << elapsed / cpuPrevious * 100 << "%";
470 cout << "CPU time for filter " << setw(25) << left << s0.str()
471 << setw(12) << left << s1.str() << endl
472 << "\t Percentage of time integration: " << setw(10)
473 << left << s2.str() << endl;
474 }
475 }
476
477 // Write out overall filter status information.
478 if (m_filtersInfosteps && m_session->GetComm()->GetRank() == 0 &&
479 !((step + 1) % m_filtersInfosteps) && !m_filters.empty())
480 {
481 stringstream ss;
482 ss << totFilterTime << "s";
483 cout << "Total filters CPU Time:\t\t\t " << setw(10) << left
484 << ss.str() << endl;
485 }
486 totFilterTime = 0.0;
487
488 // Write out checkpoint files.
489 if ((m_checksteps && !((step + 1) % m_checksteps)) || doCheckTime)
490 {
492 {
493 // Transform to physical space for output.
494 vector<bool> transformed(nfields, false);
495 for (i = 0; i < nfields; i++)
496 {
497 if (m_fields[i]->GetWaveSpace())
498 {
499 m_fields[i]->SetWaveSpace(false);
500 m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
501 m_fields[i]->UpdatePhys());
502 transformed[i] = true;
503 }
504 }
506 m_nchk++;
507
508 // Transform back to wave-space after output.
509 for (i = 0; i < nfields; i++)
510 {
511 if (transformed[i])
512 {
513 m_fields[i]->SetWaveSpace(true);
514 m_fields[i]->HomogeneousFwdTrans(
515 m_fields[i]->GetTotPoints(), m_fields[i]->GetPhys(),
516 m_fields[i]->UpdatePhys());
517 m_fields[i]->SetPhysState(false);
518 }
519 }
520 }
521 else
522 {
524 m_nchk++;
525 }
526 doCheckTime = false;
527 }
528
529 // Step advance.
530 ++step;
531 ++stepCounter;
532 }
533
534 // Print out summary statistics.
536
537 // If homogeneous, transform back into physical space if necessary.
538 if (!m_ALESolver)
539 {
541 {
542 for (i = 0; i < nfields; i++)
543 {
544 if (m_fields[i]->GetWaveSpace())
545 {
546 m_fields[i]->SetWaveSpace(false);
547 m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
548 m_fields[i]->UpdatePhys());
549 }
550 }
551 }
552 else
553 {
554 for (i = 0; i < nvariables; ++i)
555 {
556 m_fields[m_intVariables[i]]->SetPhysState(true);
557 }
558 }
559 }
560 // Finalise filters.
561 for (auto &x : m_filters)
562 {
563 x.second->Finalise(m_fields, m_time);
564 }
565
566 // Print for 1D problems.
567 if (m_spacedim == 1)
568 {
570 }
571}
virtual SOLVER_UTILS_EXPORT void v_ALEPreMultiplyMass(Array< OneD, Array< OneD, NekDouble > > &fields)
Definition: ALEHelper.cpp:108
SOLVER_UTILS_EXPORT void ALEDoElmtInvMass(Array< OneD, Array< OneD, NekDouble > > &traceNormals, Array< OneD, Array< OneD, NekDouble > > &fields, NekDouble time)
Update m_fields with u^n by multiplying by inverse mass matrix. That's then used in e....
Definition: ALEHelper.cpp:131
int m_spacedim
Spatial dimension (>= expansion dim).
NekDouble m_timestep
Time step size.
NekDouble m_fintime
Finish time of the simulation.
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
NekDouble m_checktime
Time between checkpoints.
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
enum HomogeneousType m_HomogeneousType
int m_steps
Number of steps to take.
int m_checksteps
Number of steps between checkpoints.
NekDouble m_CFLGrowth
CFL growth rate.
virtual SOLVER_UTILS_EXPORT bool v_UpdateTimeStepCheck()
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
int m_abortSteps
Number of steps between checks for abort conditions.
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate(int step)
NekDouble m_CFLEnd
Maximun cfl in cfl growth.
virtual SOLVER_UTILS_EXPORT void v_PrintSummaryStatistics(const NekDouble intTime)
Print Summary Statistics.
virtual SOLVER_UTILS_EXPORT void v_PrintStatusInformation(const int step, const NekDouble cpuTime)
Print Status Information.
bool CheckSteadyState(int step, const NekDouble &totCPUTime=0.0)
Calculate whether the system has reached a steady state by observing residuals to a user-defined tole...
SOLVER_UTILS_EXPORT void AppendOutput1D(void)
Print the solution at each solution point in a txt file.
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate(int step)
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans()
int m_steadyStateSteps
Check for steady state at step interval.
int m_filtersInfosteps
Number of time steps between outputting filters information.
bool m_homoInitialFwd
Flag to determine if simulation should start in homogeneous forward transformed state.
SOLVER_UTILS_EXPORT NekDouble GetTimeStep()
static const NekDouble kNekZeroTol
int Nnan(int n, const T *x, const int incx)
Return number of NaN elements of x.
Definition: Vmath.hpp:743

References Nektar::SolverUtils::ALEHelper::ALEDoElmtInvMass(), AppendOutput1D(), ASSERTL0, Nektar::SolverUtils::EquationSystem::Checkpoint_Output(), CheckSteadyState(), Nektar::SolverUtils::EquationSystem::eNotHomogeneous, GetTimeStep(), Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::NekConstants::kNekZeroTol, m_abortSteps, Nektar::SolverUtils::ALEHelper::m_ALESolver, m_CFLEnd, m_CFLGrowth, m_cflSafetyFactor, Nektar::SolverUtils::EquationSystem::m_checksteps, Nektar::SolverUtils::EquationSystem::m_checktime, Nektar::SolverUtils::EquationSystem::m_comm, Nektar::SolverUtils::EquationSystem::m_fields, m_filters, m_filtersInfosteps, Nektar::SolverUtils::EquationSystem::m_fintime, Nektar::SolverUtils::EquationSystem::m_HomogeneousType, m_homoInitialFwd, Nektar::SolverUtils::EquationSystem::m_infosteps, Nektar::SolverUtils::EquationSystem::m_initialStep, m_intScheme, m_intVariables, Nektar::SolverUtils::EquationSystem::m_lastCheckTime, Nektar::SolverUtils::EquationSystem::m_nchk, m_ode, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_spacedim, m_steadyStateSteps, m_steadyStateTol, Nektar::SolverUtils::EquationSystem::m_steps, Nektar::SolverUtils::EquationSystem::m_time, Nektar::SolverUtils::EquationSystem::m_TimeIncrementFactor, Nektar::SolverUtils::EquationSystem::m_timestep, Nektar::SolverUtils::EquationSystem::m_traceNormals, Vmath::Nnan(), Nektar::LibUtilities::ReduceMax, Nektar::SolverUtils::EquationSystem::SetBoundaryConditions(), Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), Nektar::LibUtilities::Timer::TimePerTest(), Nektar::SolverUtils::ALEHelper::v_ALEPreMultiplyMass(), v_PostIntegrate(), v_PreIntegrate(), v_PrintStatusInformation(), v_PrintSummaryStatistics(), v_RequireFwdTrans(), and v_UpdateTimeStepCheck().

Referenced by Nektar::CFSImplicit::v_DoSolve(), and Nektar::CoupledLinearNS::v_DoSolve().

◆ v_GenerateSummary()

void Nektar::SolverUtils::UnsteadySystem::v_GenerateSummary ( SummaryList s)
overrideprotectedvirtual

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::MMFAdvection, Nektar::UnsteadyAdvection, Nektar::UnsteadyAdvectionDiffusion, Nektar::UnsteadyInviscidBurgers, Nektar::UnsteadyViscousBurgers, Nektar::CompressibleFlowSystem, Nektar::MMFDiffusion, Nektar::ImageWarpingSystem, Nektar::CoupledLinearNS, Nektar::SmoothedProfileMethod, Nektar::VelocityCorrectionScheme, Nektar::VCSImplicit, Nektar::VCSWeakPressure, Nektar::MMFMaxwell, Nektar::PulseWavePropagation, Nektar::LinearSWE, Nektar::MMFSWE, Nektar::NonlinearPeregrine, Nektar::NonlinearSWE, Nektar::ShallowWaterSystem, Nektar::UnsteadyDiffusion, Nektar::Bidomain, Nektar::BidomainRoth, and Nektar::Monodomain.

Definition at line 639 of file UnsteadySystem.cpp.

640{
642 AddSummaryItem(s, "Advect. advancement",
643 m_explicitAdvection ? "explicit" : "implicit");
644
645 AddSummaryItem(s, "Diffuse. advancement",
646 m_explicitDiffusion ? "explicit" : "implicit");
647
648 if (m_session->GetSolverInfo("EQTYPE") ==
649 "SteadyAdvectionDiffusionReaction")
650 {
651 AddSummaryItem(s, "React. advancement",
652 m_explicitReaction ? "explicit" : "implicit");
653 }
654
655 AddSummaryItem(s, "Time Step", m_timestep);
656 AddSummaryItem(s, "No. of Steps", m_steps);
657 AddSummaryItem(s, "Checkpoints (steps)", m_checksteps);
658 if (m_intScheme)
659 {
660 AddSummaryItem(s, "Integration Type", m_intScheme->GetName());
661 }
662}
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_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion 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:47

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, and Nektar::SolverUtils::EquationSystem::v_GenerateSummary().

Referenced by Nektar::UnsteadyAdvection::v_GenerateSummary(), Nektar::UnsteadyInviscidBurgers::v_GenerateSummary(), Nektar::CompressibleFlowSystem::v_GenerateSummary(), Nektar::ImageWarpingSystem::v_GenerateSummary(), Nektar::VelocityCorrectionScheme::v_GenerateSummary(), Nektar::VCSImplicit::v_GenerateSummary(), Nektar::VCSWeakPressure::v_GenerateSummary(), Nektar::PulseWavePropagation::v_GenerateSummary(), Nektar::ShallowWaterSystem::v_GenerateSummary(), Nektar::SolverUtils::MMFSystem::v_GenerateSummary(), Nektar::UnsteadyDiffusion::v_GenerateSummary(), Nektar::Bidomain::v_GenerateSummary(), Nektar::BidomainRoth::v_GenerateSummary(), and Nektar::Monodomain::v_GenerateSummary().

◆ v_GetTimeStep()

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 761 of file UnsteadySystem.cpp.

763{
764 NEKERROR(ErrorUtil::efatal, "Not defined for this class");
765 return 0.0;
766}
#define NEKERROR(type, msg)
Assert Level 0 – Fundamental assert which is used whether in FULLDEBUG, DEBUG or OPT compilation mode...
Definition: ErrorUtil.hpp:202

References Nektar::ErrorUtil::efatal, and NEKERROR.

Referenced by GetTimeStep().

◆ v_InitObject()

void Nektar::SolverUtils::UnsteadySystem::v_InitObject ( bool  DeclareField = true)
overrideprotectedvirtual

Init object for UnsteadySystem class.

Initialization object for UnsteadySystem class.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::PulseWavePropagation, Nektar::PulseWaveSystem, Nektar::Bidomain, Nektar::BidomainRoth, Nektar::Monodomain, Nektar::NavierStokesCFE, Nektar::MMFDiffusion, Nektar::ImageWarpingSystem, Nektar::CoupledLinearNS, Nektar::IncNavierStokes, Nektar::SmoothedProfileMethod, Nektar::VCSMapping, Nektar::VelocityCorrectionScheme, Nektar::SolverUtils::FileSolution, Nektar::AcousticSystem, Nektar::APE, Nektar::LEE, Nektar::MMFAdvection, Nektar::UnsteadyAdvection, Nektar::UnsteadyAdvectionDiffusion, Nektar::UnsteadyDiffusion, Nektar::UnsteadyInviscidBurgers, Nektar::UnsteadyReactionDiffusion, Nektar::UnsteadyViscousBurgers, Nektar::CompressibleFlowSystem, Nektar::CFSImplicit, Nektar::EulerCFE, Nektar::EulerImplicitCFE, Nektar::NavierStokesCFEAxisym, Nektar::NavierStokesImplicitCFE, Nektar::Dummy, Nektar::MMFMaxwell, Nektar::LinearSWE, Nektar::MMFSWE, Nektar::NonlinearPeregrine, Nektar::NonlinearSWE, and Nektar::ShallowWaterSystem.

Definition at line 85 of file UnsteadySystem.cpp.

86{
87 EquationSystem::v_InitObject(DeclareField);
89 m_initialStep = 0;
90
91 // Load SolverInfo parameters.
92 m_session->MatchSolverInfo("DIFFUSIONADVANCEMENT", "Explicit",
94 m_session->MatchSolverInfo("ADVECTIONADVANCEMENT", "Explicit",
96 m_session->MatchSolverInfo("REACTIONADVANCEMENT", "Explicit",
97 m_explicitReaction, true);
98 m_session->LoadParameter("CheckAbortSteps", m_abortSteps, 1);
99
100 // Steady state tolerance.
101 m_session->LoadParameter("SteadyStateTol", m_steadyStateTol, 0.0);
102
103 // Frequency for checking steady state.
104 m_session->LoadParameter("SteadyStateSteps", m_steadyStateSteps, 1);
105
106 // For steady problems, we do not initialise the time integration.
107 if (m_session->DefinesSolverInfo("TimeIntegrationMethod") ||
108 m_session->DefinesTimeIntScheme())
109 {
110 LibUtilities::TimeIntScheme timeInt;
111 if (m_session->DefinesTimeIntScheme())
112 {
113 timeInt = m_session->GetTimeIntScheme();
114 }
115 else
116 {
117 std::cout << "TimeIntegrationMethod is deprecated, please use "
118 "TIMEINTEGRATIONSCHEME";
119 timeInt.method = m_session->GetSolverInfo("TimeIntegrationMethod");
120 }
121
124 timeInt.method, timeInt.variant, timeInt.order,
125 timeInt.freeParams);
126
127 // Load generic input parameters.
128 m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
129 m_session->LoadParameter("IO_FiltersInfoSteps", m_filtersInfosteps,
130 10 * m_infosteps);
131 m_session->LoadParameter("CFL", m_cflSafetyFactor, 0.0);
132 m_session->LoadParameter("CFLEnd", m_CFLEnd, 0.0);
133 m_session->LoadParameter("CFLGrowth", m_CFLGrowth, 1.0);
134
135 // Ensure that there is no conflict of parameters.
136 if (m_cflSafetyFactor > 0.0)
137 {
138 // Check final condition.
139 ASSERTL0(m_fintime == 0.0 || m_steps == 0,
140 "Final condition not unique: "
141 "fintime > 0.0 and Nsteps > 0");
142 // Check timestep condition.
143 ASSERTL0(m_timestep == 0.0,
144 "Timestep not unique: timestep > 0.0 & CFL > 0.0");
145 }
146 else
147 {
148 ASSERTL0(m_timestep != 0.0, "Need to set either TimeStep or CFL");
149 }
150
151 // Ensure that there is no conflict of parameters.
152 if (m_CFLGrowth > 1.0)
153 {
154 // Check final condition.
156 "m_CFLEnd >= m_cflSafetyFactor required");
157 }
158
159 // Set up time to be dumped in field information.
160 m_fieldMetaDataMap["Time"] = boost::lexical_cast<std::string>(m_time);
161 }
162
163 // By default attempt to forward transform initial condition.
164 m_homoInitialFwd = true;
165
166 // Set up filters.
167 for (auto &x : m_session->GetFilters())
168 {
169 m_filters.push_back(make_pair(
170 x.name, GetFilterFactory().CreateInstance(
171 x.name, m_session, shared_from_this(), x.params)));
172 }
173}
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
virtual SOLVER_UTILS_EXPORT void v_ALEInitObject(int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
Definition: ALEHelper.cpp:42
virtual SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareFeld=true)
Initialisation object for EquationSystem.
TimeIntegrationSchemeFactory & GetTimeIntegrationSchemeFactory()
FilterFactory & GetFilterFactory()

References ASSERTL0, Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::LibUtilities::TimeIntScheme::freeParams, Nektar::SolverUtils::GetFilterFactory(), Nektar::LibUtilities::GetTimeIntegrationSchemeFactory(), m_abortSteps, m_CFLEnd, m_CFLGrowth, m_cflSafetyFactor, m_explicitAdvection, m_explicitDiffusion, m_explicitReaction, Nektar::SolverUtils::EquationSystem::m_fieldMetaDataMap, Nektar::SolverUtils::EquationSystem::m_fields, m_filters, m_filtersInfosteps, Nektar::SolverUtils::EquationSystem::m_fintime, m_homoInitialFwd, Nektar::SolverUtils::EquationSystem::m_infosteps, Nektar::SolverUtils::EquationSystem::m_initialStep, m_intScheme, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_spacedim, m_steadyStateSteps, m_steadyStateTol, Nektar::SolverUtils::EquationSystem::m_steps, Nektar::SolverUtils::EquationSystem::m_time, Nektar::SolverUtils::EquationSystem::m_timestep, Nektar::LibUtilities::TimeIntScheme::method, Nektar::LibUtilities::TimeIntScheme::order, Nektar::SolverUtils::ALEHelper::v_ALEInitObject(), Nektar::SolverUtils::EquationSystem::v_InitObject(), and Nektar::LibUtilities::TimeIntScheme::variant.

Referenced by Nektar::PulseWaveSystem::v_InitObject(), Nektar::SolverUtils::AdvectionSystem::v_InitObject(), Nektar::Bidomain::v_InitObject(), Nektar::BidomainRoth::v_InitObject(), Nektar::Monodomain::v_InitObject(), Nektar::MMFDiffusion::v_InitObject(), Nektar::MMFAdvection::v_InitObject(), Nektar::UnsteadyDiffusion::v_InitObject(), Nektar::Dummy::v_InitObject(), Nektar::MMFMaxwell::v_InitObject(), Nektar::MMFSWE::v_InitObject(), and Nektar::ShallowWaterSystem::v_InitObject().

◆ v_PostIntegrate()

bool Nektar::SolverUtils::UnsteadySystem::v_PostIntegrate ( int  step)
protectedvirtual

◆ v_PreIntegrate()

bool Nektar::SolverUtils::UnsteadySystem::v_PreIntegrate ( int  step)
protectedvirtual

◆ v_PrintStatusInformation()

void Nektar::SolverUtils::UnsteadySystem::v_PrintStatusInformation ( const int  step,
const NekDouble  cpuTime 
)
protectedvirtual

Print Status Information.

Reimplemented in Nektar::CFSImplicit.

Definition at line 573 of file UnsteadySystem.cpp.

575{
576 if (m_infosteps && m_session->GetComm()->GetSpaceComm()->GetRank() == 0 &&
577 !((step + 1) % m_infosteps))
578 {
579 if (m_comm->IsParallelInTime())
580 {
581 cout << "RANK " << m_session->GetComm()->GetTimeComm()->GetRank()
582 << " Steps: " << setw(8) << left << step + 1 << " "
583 << "Time: " << setw(12) << left << m_time;
584 }
585 else
586 {
587 cout << "Steps: " << setw(8) << left << step + 1 << " "
588 << "Time: " << setw(12) << left << m_time;
589 }
590
592 {
593 cout << " Time-step: " << setw(12) << left << m_timestep;
594 }
595
596 stringstream ss;
597 ss << cpuTime << "s";
598 cout << " CPU Time: " << setw(8) << left << ss.str() << endl;
599 }
600}

References m_cflSafetyFactor, Nektar::SolverUtils::EquationSystem::m_comm, Nektar::SolverUtils::EquationSystem::m_infosteps, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_time, and Nektar::SolverUtils::EquationSystem::m_timestep.

Referenced by v_DoSolve(), and Nektar::CFSImplicit::v_PrintStatusInformation().

◆ v_PrintSummaryStatistics()

void Nektar::SolverUtils::UnsteadySystem::v_PrintSummaryStatistics ( const NekDouble  intTime)
protectedvirtual

Print Summary Statistics.

Reimplemented in Nektar::CFSImplicit.

Definition at line 602 of file UnsteadySystem.cpp.

603{
604 if (m_session->GetComm()->GetRank() == 0)
605 {
606 if (m_cflSafetyFactor > 0.0)
607 {
608 cout << "CFL safety factor : " << m_cflSafetyFactor << endl
609 << "CFL time-step : " << m_timestep << endl;
610 }
611
612 if (m_session->GetSolverInfo("Driver") != "SteadyState" &&
613 m_session->GetSolverInfo("Driver") != "Parareal" &&
614 m_session->GetSolverInfo("Driver") != "PFASST")
615 {
616 cout << "Time-integration : " << intTime << "s" << endl;
617 }
618 }
619}

References m_cflSafetyFactor, Nektar::SolverUtils::EquationSystem::m_session, and Nektar::SolverUtils::EquationSystem::m_timestep.

Referenced by v_DoSolve(), and Nektar::CFSImplicit::v_PrintSummaryStatistics().

◆ v_RequireFwdTrans()

bool Nektar::SolverUtils::UnsteadySystem::v_RequireFwdTrans ( void  )
protectedvirtual

Reimplemented in Nektar::Dummy, Nektar::VelocityCorrectionScheme, and Nektar::SolverUtils::FileSolution.

Definition at line 787 of file UnsteadySystem.cpp.

788{
789 return true;
790}

Referenced by v_DoSolve().

◆ v_SteadyStateResidual()

void Nektar::SolverUtils::UnsteadySystem::v_SteadyStateResidual ( int  step,
Array< OneD, NekDouble > &  L2 
)
protectedvirtual

Reimplemented in Nektar::CompressibleFlowSystem.

Definition at line 943 of file UnsteadySystem.cpp.

945{
946 const int nPoints = GetTotPoints();
947 const int nFields = m_fields.size();
948
949 // Holds L2 errors.
950 Array<OneD, NekDouble> RHSL2(nFields);
951 Array<OneD, NekDouble> residual(nFields);
952 Array<OneD, NekDouble> reference(nFields);
953
954 for (int i = 0; i < nFields; ++i)
955 {
956 Array<OneD, NekDouble> tmp(nPoints);
957
958 Vmath::Vsub(nPoints, m_fields[i]->GetPhys(), 1, m_previousSolution[i],
959 1, tmp, 1);
960 Vmath::Vmul(nPoints, tmp, 1, tmp, 1, tmp, 1);
961 residual[i] = Vmath::Vsum(nPoints, tmp, 1);
962
964 tmp, 1);
965 reference[i] = Vmath::Vsum(nPoints, tmp, 1);
966 }
967
968 m_comm->GetSpaceComm()->AllReduce(residual, LibUtilities::ReduceSum);
969 m_comm->GetSpaceComm()->AllReduce(reference, LibUtilities::ReduceSum);
970
971 // L2 error.
972 for (int i = 0; i < nFields; ++i)
973 {
974 reference[i] = (reference[i] == 0) ? 1 : reference[i];
975 L2[i] = sqrt(residual[i] / reference[i]);
976 }
977}
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.hpp:608
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.hpp:220
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

References Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_comm, Nektar::SolverUtils::EquationSystem::m_fields, m_previousSolution, Nektar::LibUtilities::ReduceSum, tinysimd::sqrt(), Vmath::Vmul(), Vmath::Vsub(), and Vmath::Vsum().

Referenced by SteadyStateResidual().

◆ v_UpdateTimeStepCheck()

bool Nektar::SolverUtils::UnsteadySystem::v_UpdateTimeStepCheck ( void  )
protectedvirtual

Reimplemented in Nektar::CFSImplicit.

Definition at line 795 of file UnsteadySystem.cpp.

796{
797 return true;
798}

Referenced by v_DoSolve().

Member Data Documentation

◆ cmdSetStartChkNum

std::string Nektar::SolverUtils::UnsteadySystem::cmdSetStartChkNum
static
Initial value:
=
"set-start-chknumber", "",
"Set the starting number of the checkpoint file.")
static std::string RegisterCmdLineArgument(const std::string &pName, const std::string &pShortName, const std::string &pDescription)
Registers a command-line argument with the session reader.

Definition at line 82 of file UnsteadySystem.h.

◆ cmdSetStartTime

std::string Nektar::SolverUtils::UnsteadySystem::cmdSetStartTime
static
Initial value:
=
"set-start-time", "", "Set the starting time of the simulation.")

Definition at line 81 of file UnsteadySystem.h.

◆ m_abortSteps

int Nektar::SolverUtils::UnsteadySystem::m_abortSteps
protected

Number of steps between checks for abort conditions.

Definition at line 104 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

◆ m_CFLEnd

NekDouble Nektar::SolverUtils::UnsteadySystem::m_CFLEnd
protected

Maximun cfl in cfl growth.

Definition at line 101 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

◆ m_CFLGrowth

NekDouble Nektar::SolverUtils::UnsteadySystem::m_CFLGrowth
protected

CFL growth rate.

Definition at line 99 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

◆ m_cflSafetyFactor

NekDouble Nektar::SolverUtils::UnsteadySystem::m_cflSafetyFactor
protected

◆ m_epsilon

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

Diffusion coefficient.

Definition at line 130 of file UnsteadySystem.h.

◆ m_errFile

std::ofstream Nektar::SolverUtils::UnsteadySystem::m_errFile
protected

Definition at line 127 of file UnsteadySystem.h.

Referenced by CheckSteadyState(), and InitializeSteadyState().

◆ m_explicitAdvection

bool Nektar::SolverUtils::UnsteadySystem::m_explicitAdvection
protected

◆ m_explicitDiffusion

bool Nektar::SolverUtils::UnsteadySystem::m_explicitDiffusion
protected

◆ m_explicitReaction

bool Nektar::SolverUtils::UnsteadySystem::m_explicitReaction
protected

Indicates if explicit or implicit treatment of reaction is used.

Definition at line 111 of file UnsteadySystem.h.

Referenced by v_GenerateSummary(), and v_InitObject().

◆ m_filters

std::vector<std::pair<std::string, FilterSharedPtr> > Nektar::SolverUtils::UnsteadySystem::m_filters
protected

◆ m_filtersInfosteps

int Nektar::SolverUtils::UnsteadySystem::m_filtersInfosteps
protected

Number of time steps between outputting filters information.

Definition at line 119 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

◆ m_homoInitialFwd

bool Nektar::SolverUtils::UnsteadySystem::m_homoInitialFwd
protected

◆ m_intScheme

LibUtilities::TimeIntegrationSchemeSharedPtr Nektar::SolverUtils::UnsteadySystem::m_intScheme
protected

◆ m_intVariables

std::vector<int> Nektar::SolverUtils::UnsteadySystem::m_intVariables
protected

◆ m_ode

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

◆ m_previousSolution

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::UnsteadySystem::m_previousSolution
protected

Storage for previous solution for steady-state check.

Definition at line 92 of file UnsteadySystem.h.

Referenced by CheckSteadyState(), InitializeSteadyState(), and v_SteadyStateResidual().

◆ m_steadyStateSteps

int Nektar::SolverUtils::UnsteadySystem::m_steadyStateSteps
protected

Check for steady state at step interval.

Definition at line 114 of file UnsteadySystem.h.

Referenced by v_DoSolve(), and v_InitObject().

◆ m_steadyStateTol

NekDouble Nektar::SolverUtils::UnsteadySystem::m_steadyStateTol
protected

Tolerance to which steady state should be evaluated at.

Definition at line 116 of file UnsteadySystem.h.

Referenced by CheckSteadyState(), InitializeSteadyState(), v_DoSolve(), and v_InitObject().