Nektar++
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator 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)
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time)
 
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff (const Array< OneD, Array< OneD, NekDouble > > vel, StdRegions::VarCoeffMap &varCoeffMap)
 Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::EquationSystem
SOLVER_UTILS_EXPORT EquationSystem (const LibUtilities::SessionReaderSharedPtr &pSession)
 Initialises EquationSystem class members. More...
 
int nocase_cmp (const std::string &s1, const std::string &s2)
 
virtual SOLVER_UTILS_EXPORT
NekDouble 
v_LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Virtual function for the L_inf error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT
NekDouble 
v_L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the L_2 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
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 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...
 
LibUtilities::TimeIntegrationWrapperSharedPtr m_intScheme
 Wrapper to the time integration scheme. More...
 
LibUtilities::TimeIntegrationSchemeOperators m_ode
 The time integration scheme operators to use. More...
 
LibUtilities::TimeIntegrationSolutionSharedPtr m_intSoln
 
NekDouble m_epsilon
 
bool m_explicitDiffusion
 Indicates if explicit or implicit treatment of diffusion is used. More...
 
bool m_explicitAdvection
 Indicates if explicit or implicit treatment of advection is used. More...
 
bool m_explicitReaction
 Indicates if explicit or implicit treatment of reaction is used. More...
 
bool m_homoInitialFwd
 Flag to determine if simulation should start in homogeneous forward transformed state. More...
 
std::vector< int > m_intVariables
 
std::vector< FilterSharedPtrm_filters
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
LibUtilities::SessionReaderSharedPtr m_session
 The session reader. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 Field input/output. More...
 
std::map< std::string, Array
< OneD, Array< OneD, float > > > 
m_interpWeights
 Map of the interpolation weights for a specific filename. More...
 
std::map< std::string, Array
< OneD, Array< OneD, unsigned
int > > > 
m_interpInds
 Map of the interpolation indices for a specific filename. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_fields
 Array holding all dependent variables. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_base
 Base fields. More...
 
Array< OneD,
MultiRegions::ExpListSharedPtr
m_derivedfields
 Array holding all dependent variables. More...
 
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions
 Pointer to boundary conditions object. More...
 
SpatialDomains::MeshGraphSharedPtr m_graph
 Pointer to graph defining mesh. More...
 
std::string m_sessionName
 Name of the session. More...
 
NekDouble m_time
 Current time of simulation. More...
 
int m_initialStep
 Number of the step where the simulation should begin. More...
 
NekDouble m_fintime
 Finish time of the simulation. More...
 
NekDouble m_timestep
 Time step size. More...
 
NekDouble m_lambda
 Lambda constant in real system if one required. More...
 
std::set< std::string > m_loadedFields
 
NekDouble m_checktime
 Time between checkpoints. More...
 
int m_nchk
 Number of checkpoints written so far. More...
 
int m_steps
 Number of steps to take. More...
 
int m_checksteps
 Number of steps between checkpoints. More...
 
int m_spacedim
 Spatial dimension (>= expansion dim). More...
 
int m_expdim
 Expansion dimension. More...
 
bool m_singleMode
 Flag to determine if single homogeneous mode is used. More...
 
bool m_halfMode
 Flag to determine if half homogeneous mode is used. More...
 
bool m_multipleModes
 Flag to determine if use multiple homogenenous modes are used. More...
 
bool m_useFFT
 Flag to determine if FFT is used for homogeneous transform. More...
 
bool m_homogen_dealiasing
 Flag to determine if dealiasing is used for homogeneous simulations. More...
 
bool m_specHP_dealiasing
 Flag to determine if dealisising is usde for the Spectral/hp element discretisation. More...
 
enum MultiRegions::ProjectionType m_projectionType
 Type of projection; e.g continuous or discontinuous. More...
 
Array< OneD, Array< OneD,
NekDouble > > 
m_traceNormals
 Array holding trace normals for DG simulations in the forwards direction. More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_gradtan
 1 x nvariable x nq More...
 
Array< OneD, Array< OneD,
Array< OneD, NekDouble > > > 
m_tanbasis
 2 x m_spacedim x nq More...
 
Array< OneD, bool > m_checkIfSystemSingular
 Flag to indicate if the fields should be checked for singularity. More...
 
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
 Map to identify relevant solver info to dump in output fields. More...
 
int m_NumQuadPointsError
 Number of Quadrature points used to work out the error. More...
 
enum HomogeneousType m_HomogeneousType
 
NekDouble m_LhomX
 physical length in X direction (if homogeneous) More...
 
NekDouble m_LhomY
 physical length in Y direction (if homogeneous) More...
 
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous) More...
 
int m_npointsX
 number of points in X direction (if homogeneous) More...
 
int m_npointsY
 number of points in Y direction (if homogeneous) More...
 
int m_npointsZ
 number of points in Z direction (if homogeneous) More...
 
int m_HomoDirec
 number of homogenous directions More...
 

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

127  {
128  }
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)
protected

Definition at line 690 of file UnsteadySystem.cpp.

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

Referenced by v_DoInitialise().

691  {
692  if (m_session->DefinesFunction("InitialConditions"))
693  {
694  for (int i = 0; i < m_fields.num_elements(); ++i)
695  {
697 
698  vType = m_session->GetFunctionType(
699  "InitialConditions", m_session->GetVariable(i));
700 
701  if (vType == LibUtilities::eFunctionTypeFile)
702  {
703  std::string filename
704  = 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  }
716  m_fld->ImportFieldMetaData(filename, m_fieldMetaDataMap);
717 
718  // check to see if time defined
719  if (m_fieldMetaDataMap !=
721  {
723 
724  iter = m_fieldMetaDataMap.find("Time");
725  if (iter != m_fieldMetaDataMap.end())
726  {
727  time = boost::lexical_cast<NekDouble>(
728  iter->second);
729  }
730  }
731 
732  break;
733  }
734  }
735  }
736  }
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
double NekDouble
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::FieldIOSharedPtr m_fld
Field input/output.
static FieldMetaDataMap NullFieldMetaDataMap
Definition: FieldIO.h:54
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 897 of file UnsteadySystem.cpp.

References v_GetTimeStep().

899  {
900  return v_GetTimeStep(inarray);
901  }
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 133 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.

134  {
135  NekDouble TimeStability = 0.0;
136  switch(m_intScheme->GetIntegrationMethod())
137  {
141  {
142  TimeStability = 2.784;
143  break;
144  }
151  {
152  TimeStability = 2.0;
153  break;
154  }
156  {
157  TimeStability = 1.0;
158  break;
159  }
160  default:
161  {
162  ASSERTL0(
163  false,
164  "No CFL control implementation for this time"
165  "integration scheme");
166  }
167  }
168  return TimeStability;
169  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 931 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().

934  {
935  int phystot = m_fields[0]->GetTotPoints();
936  int nvel = vel.num_elements();
937 
938  Array<OneD, NekDouble> varcoeff(phystot),tmp;
939 
940  // calculate magnitude of v
941  Vmath::Vmul(phystot,vel[0],1,vel[0],1,varcoeff,1);
942  for(int n = 1; n < nvel; ++n)
943  {
944  Vmath::Vvtvp(phystot,vel[n],1,vel[n],1,varcoeff,1,varcoeff,1);
945  }
946  Vmath::Vsqrt(phystot,varcoeff,1,varcoeff,1);
947 
948  for(int i = 0; i < m_fields[0]->GetNumElmts(); ++i)
949  {
950  int offset = m_fields[0]->GetPhys_Offset(i);
951  int nq = m_fields[0]->GetExp(i)->GetTotPoints();
952  Array<OneD, NekDouble> unit(nq,1.0);
953 
954  int nmodes = 0;
955 
956  for(int n = 0; n < m_fields[0]->GetExp(i)->GetNumBases(); ++n)
957  {
958  nmodes = max(nmodes,
959  m_fields[0]->GetExp(i)->GetBasisNumModes(n));
960  }
961 
962  NekDouble h = m_fields[0]->GetExp(i)->Integral(unit);
963  h = pow(h,(NekDouble) (1.0/nvel))/((NekDouble) nmodes);
964 
965  Vmath::Smul(nq,h,varcoeff+offset,1,tmp = varcoeff+offset,1);
966  }
967 
968  // set up map with eVarCoffLaplacian key
969  varCoeffMap[StdRegions::eVarCoeffLaplacian] = varcoeff;
970  }
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.cpp:394
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:428
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:199
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:169
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 494 of file UnsteadySystem.cpp.

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

Referenced by v_DoSolve().

496  {
497  // Coordinates of the quadrature points in the real physical space
498  Array<OneD,NekDouble> x(GetNpoints());
499  Array<OneD,NekDouble> y(GetNpoints());
500  Array<OneD,NekDouble> z(GetNpoints());
501  m_fields[0]->GetCoords(x, y, z);
502 
503  // Print out the solution in a txt file
504  ofstream outfile;
505  outfile.open("solution1D.txt");
506  for(int i = 0; i < GetNpoints(); i++)
507  {
508  outfile << scientific << setw (17) << setprecision(16) << x[i]
509  << " " << solution1D[0][i] << endl;
510  }
511  outfile << endl << endl;
512  outfile.close();
513  }
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 449 of file UnsteadySystem.cpp.

References CheckForRestartTime(), 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().

450  {
454  }
SOLVER_UTILS_EXPORT void CheckForRestartTime(NekDouble &time)
NekDouble m_time
Current time of simulation.
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 175 of file UnsteadySystem.cpp.

References ASSERTL0, Nektar::SolverUtils::EquationSystem::Checkpoint_Output(), Nektar::SolverUtils::EquationSystem::eHomogeneous1D, 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, 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(), and v_PreIntegrate().

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

176  {
177  ASSERTL0(m_intScheme != 0, "No time integration scheme.");
178 
179  int i = 1;
180  int nvariables = 0;
181  int nfields = m_fields.num_elements();
182 
183  if (m_intVariables.empty())
184  {
185  for (i = 0; i < nfields; ++i)
186  {
187  m_intVariables.push_back(i);
188  }
189  nvariables = nfields;
190  }
191  else
192  {
193  nvariables = m_intVariables.size();
194  }
195 
196  // Integrate in wave-space if using homogeneous1D
198  {
199  for(i = 0; i < nfields; ++i)
200  {
201  m_fields[i]->HomogeneousFwdTrans(m_fields[i]->GetPhys(),
202  m_fields[i]->UpdatePhys());
203  m_fields[i]->SetWaveSpace(true);
204  m_fields[i]->SetPhysState(false);
205  }
206  }
207 
208  // Set up wrapper to fields data storage.
209  Array<OneD, Array<OneD, NekDouble> > fields(nvariables);
210  Array<OneD, Array<OneD, NekDouble> > tmp (nvariables);
211 
212  // Order storage to list time-integrated fields first.
213  for(i = 0; i < nvariables; ++i)
214  {
215  fields[i] = m_fields[m_intVariables[i]]->GetPhys();
216  m_fields[m_intVariables[i]]->SetPhysState(false);
217  }
218 
219  // Initialise time integration scheme
220  m_intSoln = m_intScheme->InitializeScheme(
221  m_timestep, fields, m_time, m_ode);
222 
223  // Initialise filters
225  for (x = m_filters.begin(); x != m_filters.end(); ++x)
226  {
227  (*x)->Initialise(m_fields, m_time);
228  }
229 
230  // Ensure that there is no conflict of parameters
231  if(m_cflSafetyFactor > 0.0)
232  {
233  // Check final condition
234  ASSERTL0(m_fintime == 0.0 || m_steps == 0,
235  "Final condition not unique: "
236  "fintime > 0.0 and Nsteps > 0");
237 
238  // Check timestep condition
239  ASSERTL0(m_timestep == 0.0,
240  "Timestep not unique: timestep > 0.0 & CFL > 0.0");
241  }
242 
243  // Check uniqueness of checkpoint output
244  ASSERTL0((m_checktime == 0.0 && m_checksteps == 0) ||
245  (m_checktime > 0.0 && m_checksteps == 0) ||
246  (m_checktime == 0.0 && m_checksteps > 0),
247  "Only one of IO_CheckTime and IO_CheckSteps "
248  "should be set!");
249 
250  Timer timer;
251  bool doCheckTime = false;
252  int step = m_initialStep;
253  NekDouble intTime = 0.0;
254  NekDouble lastCheckTime = 0.0;
255  NekDouble cpuTime = 0.0;
256  NekDouble elapsed = 0.0;
257 
258  while (step < m_steps ||
260  {
261  if (m_cflSafetyFactor)
262  {
263  m_timestep = GetTimeStep(fields);
264 
265  // Ensure that the final timestep finishes at the final
266  // time, or at a prescribed IO_CheckTime.
267  if (m_time + m_timestep > m_fintime && m_fintime > 0.0)
268  {
270  }
271  else if (m_checktime &&
272  m_time + m_timestep - lastCheckTime >= m_checktime)
273  {
274  lastCheckTime += m_checktime;
275  m_timestep = lastCheckTime - m_time;
276  doCheckTime = true;
277  }
278  }
279 
280  // Perform any solver-specific pre-integration steps
281  timer.Start();
282  if (v_PreIntegrate(step))
283  {
284  break;
285  }
286 
287  fields = m_intScheme->TimeIntegrate(
288  step, m_timestep, m_intSoln, m_ode);
289  timer.Stop();
290 
291  m_time += m_timestep;
292  elapsed = timer.TimePerTest(1);
293  intTime += elapsed;
294  cpuTime += elapsed;
295 
296  // Write out status information
297  if (m_session->GetComm()->GetRank() == 0 &&
298  !((step+1) % m_infosteps))
299  {
300  cout << "Steps: " << setw(8) << left << step+1 << " "
301  << "Time: " << setw(12) << left << m_time;
302 
303  if (m_cflSafetyFactor)
304  {
305  cout << " Time-step: " << setw(12)
306  << left << m_timestep;
307  }
308 
309  stringstream ss;
310  ss << cpuTime << "s";
311  cout << " CPU Time: " << setw(8) << left
312  << ss.str() << endl;
313  cpuTime = 0.0;
314  }
315 
316  // Transform data into coefficient space
317  for (i = 0; i < nvariables; ++i)
318  {
319  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
320  m_fields[m_intVariables[i]]->FwdTrans_IterPerExp(
321  fields[i],
322  m_fields[m_intVariables[i]]->UpdateCoeffs());
323  m_fields[m_intVariables[i]]->SetPhysState(false);
324  }
325 
326  // Perform any solver-specific post-integration steps
327  if (v_PostIntegrate(step))
328  {
329  break;
330  }
331 
332  // search for NaN and quit if found
333  int nanFound = 0;
334  for (i = 0; i < nvariables; ++i)
335  {
336  if (Vmath::Nnan(fields[i].num_elements(), fields[i], 1) > 0)
337  {
338  nanFound = 1;
339  }
340  }
341  m_session->GetComm()->AllReduce(nanFound,
343  ASSERTL0 (!nanFound,
344  "NaN found during time integration.");
345 
346  // Update filters
348  for (x = m_filters.begin(); x != m_filters.end(); ++x)
349  {
350  (*x)->Update(m_fields, m_time);
351  }
352 
353  // Write out checkpoint files
354  if ((m_checksteps && step && !((step + 1) % m_checksteps)) ||
355  doCheckTime)
356  {
358  {
359  vector<bool> transformed(nfields, false);
360  for(i = 0; i < nfields; i++)
361  {
362  if (m_fields[i]->GetWaveSpace())
363  {
364  m_fields[i]->SetWaveSpace(false);
365  m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
366  m_fields[i]->UpdatePhys());
367  m_fields[i]->SetPhysState(true);
368  transformed[i] = true;
369  }
370  }
372  for(i = 0; i < nfields; i++)
373  {
374  if (transformed[i])
375  {
376  m_fields[i]->SetWaveSpace(true);
377  m_fields[i]->HomogeneousFwdTrans(
378  m_fields[i]->GetPhys(),
379  m_fields[i]->UpdatePhys());
380  m_fields[i]->SetPhysState(false);
381  }
382  }
383  }
384  else
385  {
387  }
388  doCheckTime = false;
389  }
390 
391  // Step advance
392  ++step;
393  }
394 
395  // Print out summary statistics
396  if (m_session->GetComm()->GetRank() == 0)
397  {
398  if (m_cflSafetyFactor > 0.0)
399  {
400  cout << "CFL safety factor : " << m_cflSafetyFactor << endl
401  << "CFL time-step : " << m_timestep << endl;
402  }
403 
404  if (m_session->GetSolverInfo("Driver") != "SteadyState")
405  {
406  cout << "Time-integration : " << intTime << "s" << endl;
407  }
408  }
409 
410  // If homogeneous, transform back into physical space if necessary.
412  {
413  for(i = 0; i < nfields; i++)
414  {
415  if (m_fields[i]->GetWaveSpace())
416  {
417  m_fields[i]->SetWaveSpace(false);
418  m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(),
419  m_fields[i]->UpdatePhys());
420  m_fields[i]->SetPhysState(true);
421  }
422  }
423  }
424  else
425  {
426  for(i = 0; i < nvariables; ++i)
427  {
428  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
429  m_fields[m_intVariables[i]]->SetPhysState(true);
430  }
431  }
432 
433  // Finalise filters
434  for (x = m_filters.begin(); x != m_filters.end(); ++x)
435  {
436  (*x)->Finalise(m_fields, m_time);
437  }
438 
439  // Print for 1D problems
440  if(m_spacedim == 1)
441  {
442  v_AppendOutput1D(fields);
443  }
444  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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.
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:878
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::CompressibleFlowSystem, Nektar::UnsteadyViscousBurgers, Nektar::VelocityCorrectionScheme, Nektar::UnsteadyAdvectionDiffusion, Nektar::CFLtester, Nektar::UnsteadyAdvection, Nektar::NonlinearPeregrine, Nektar::PulseWavePropagation, Nektar::Monodomain, Nektar::Bidomain, Nektar::BidomainRoth, Nektar::ShallowWaterSystem, Nektar::LinearSWE, Nektar::EulerADCFE, Nektar::EulerCFE, Nektar::NonlinearSWE, Nektar::NavierStokesCFE, Nektar::ImageWarpingSystem, and Nektar::UnsteadyDiffusion.

Definition at line 460 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::UnsteadyDiffusion::v_GenerateSummary(), Nektar::ImageWarpingSystem::v_GenerateSummary(), Nektar::ShallowWaterSystem::v_GenerateSummary(), Nektar::Bidomain::v_GenerateSummary(), Nektar::BidomainRoth::v_GenerateSummary(), Nektar::Monodomain::v_GenerateSummary(), Nektar::PulseWavePropagation::v_GenerateSummary(), Nektar::UnsteadyAdvection::v_GenerateSummary(), Nektar::CFLtester::v_GenerateSummary(), Nektar::UnsteadyAdvectionDiffusion::v_GenerateSummary(), Nektar::VelocityCorrectionScheme::v_GenerateSummary(), Nektar::UnsteadyViscousBurgers::v_GenerateSummary(), and Nektar::CompressibleFlowSystem::v_GenerateSummary().

461  {
463  AddSummaryItem(s, "Advection",
464  m_explicitAdvection ? "explicit" : "implicit");
465 
466  if(m_session->DefinesSolverInfo("AdvectionType"))
467  {
468  AddSummaryItem(s, "AdvectionType",
469  m_session->GetSolverInfo("AdvectionType"));
470  }
471 
472  AddSummaryItem(s, "Diffusion",
473  m_explicitDiffusion ? "explicit" : "implicit");
474 
475  if (m_session->GetSolverInfo("EQTYPE")
476  == "SteadyAdvectionDiffusionReaction")
477  {
478  AddSummaryItem(s, "Reaction",
479  m_explicitReaction ? "explicit" : "implicit");
480  }
481 
482  AddSummaryItem(s, "Time Step", m_timestep);
483  AddSummaryItem(s, "No. of Steps", m_steps);
484  AddSummaryItem(s, "Checkpoints (steps)", m_checksteps);
485  AddSummaryItem(s, "Integration Type",
487  m_intScheme->GetIntegrationMethod()]);
488  }
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 909 of file UnsteadySystem.cpp.

References ASSERTL0.

Referenced by GetTimeStep().

911  {
912  ASSERTL0(false, "Not defined for this class");
913  return 0.0;
914  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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::CompressibleFlowSystem, Nektar::PulseWaveSystem, Nektar::UnsteadyViscousBurgers, Nektar::UnsteadyAdvectionDiffusion, Nektar::IncNavierStokes, Nektar::CFLtester, Nektar::UnsteadyAdvection, Nektar::UnsteadyInviscidBurger, Nektar::ShallowWaterSystem, Nektar::NonlinearPeregrine, Nektar::APE, Nektar::EulerADCFE, Nektar::NavierStokesCFE, Nektar::EulerCFE, Nektar::PulseWavePropagation, Nektar::ImageWarpingSystem, Nektar::UnsteadyDiffusion, Nektar::Monodomain, Nektar::Bidomain, Nektar::BidomainRoth, Nektar::LinearSWE, Nektar::NonlinearSWE, Nektar::VCSMapping, Nektar::PulseWaveSystemOutput, 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, 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::APE::v_InitObject(), Nektar::ShallowWaterSystem::v_InitObject(), Nektar::CFLtester::v_InitObject(), Nektar::PulseWaveSystem::v_InitObject(), and Nektar::CompressibleFlowSystem::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  // For steady problems, we do not initialise the time integration
93  if (m_session->DefinesSolverInfo("TIMEINTEGRATIONMETHOD"))
94  {
96  CreateInstance(m_session->GetSolverInfo(
97  "TIMEINTEGRATIONMETHOD"));
98 
99  // Load generic input parameters
100  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
101  m_session->LoadParameter("CFL", m_cflSafetyFactor, 0.0);
102 
103  // Set up time to be dumped in field information
104  m_fieldMetaDataMap["Time"] =
105  boost::lexical_cast<std::string>(m_time);
106  }
107 
108  // By default attempt to forward transform initial condition.
109  m_homoInitialFwd = true;
110 
111  // Set up filters
112  LibUtilities::FilterMap::const_iterator x;
113  LibUtilities::FilterMap f = m_session->GetFilters();
114  for (x = f.begin(); x != f.end(); ++x)
115  {
116  m_filters.push_back(GetFilterFactory().CreateInstance(
117  x->first,
118  m_session,
119  x->second));
120  }
121  }
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 515 of file UnsteadySystem.cpp.

References ASSERTL0.

518  {
519  ASSERTL0(false,
520  "This function is not implemented for this equation.");
521  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 523 of file UnsteadySystem.cpp.

References ASSERTL0.

527  {
528  ASSERTL0(false,
529  "This function is not implemented for this equation.");
530  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:161
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 532 of file UnsteadySystem.cpp.

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

535  {
536  int i, j;
537  int nTraceNumPoints = GetTraceNpoints();
538  int nvariables = m_fields.num_elements();
539  int nqvar = uflux.num_elements();
540 
541  Array<OneD, NekDouble > Fwd (nTraceNumPoints);
542  Array<OneD, NekDouble > Bwd (nTraceNumPoints);
543  Array<OneD, NekDouble > Vn (nTraceNumPoints, 0.0);
544  Array<OneD, NekDouble > fluxtemp(nTraceNumPoints, 0.0);
545 
546  // Get the sign of (v \cdot n), v = an arbitrary vector
547 
548  // Evaulate upwind flux:
549  // uflux = \hat{u} \phi \cdot u = u^{(+,-)} n
550  for (j = 0; j < nqvar; ++j)
551  {
552  for (i = 0; i < nvariables ; ++i)
553  {
554  // Compute Fwd and Bwd value of ufield of i direction
555  m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
556 
557  // if Vn >= 0, flux = uFwd, i.e.,
558  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick uflux = uFwd
559  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick uflux = uFwd
560 
561  // else if Vn < 0, flux = uBwd, i.e.,
562  // edge::eForward, if V*n<0 <=> V*n_F<0, pick uflux = uBwd
563  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick uflux = uBwd
564 
565  m_fields[i]->GetTrace()->Upwind(m_traceNormals[j],
566  Fwd, Bwd, fluxtemp);
567 
568  // Imposing weak boundary condition with flux
569  // if Vn >= 0, uflux = uBwd at Neumann, i.e.,
570  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick uflux = uBwd
571  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick uflux = uBwd
572 
573  // if Vn >= 0, uflux = uFwd at Neumann, i.e.,
574  // edge::eForward, if V*n<0 <=> V*n_F<0, pick uflux = uFwd
575  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick uflux = uFwd
576 
577  if(m_fields[0]->GetBndCondExpansions().num_elements())
578  {
579  WeakPenaltyforScalar(i, ufield[i], fluxtemp);
580  }
581 
582  // if Vn >= 0, flux = uFwd*(tan_{\xi}^- \cdot \vec{n}),
583  // i.e,
584  // edge::eForward, uFwd \(\tan_{\xi}^Fwd \cdot \vec{n})
585  // edge::eBackward, uFwd \(\tan_{\xi}^Bwd \cdot \vec{n})
586 
587  // else if Vn < 0, flux = uBwd*(tan_{\xi}^- \cdot \vec{n}),
588  // i.e,
589  // edge::eForward, uBwd \(\tan_{\xi}^Fwd \cdot \vec{n})
590  // edge::eBackward, uBwd \(\tan_{\xi}^Bwd \cdot \vec{n})
591 
592  Vmath::Vmul(nTraceNumPoints,
593  m_traceNormals[j], 1,
594  fluxtemp, 1,
595  uflux[j][i], 1);
596  }
597  }
598  }
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:169
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 602 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().

606  {
607  int nTraceNumPoints = GetTraceNpoints();
608  int nvariables = m_fields.num_elements();
609  int nqvar = qfield.num_elements();
610 
611  NekDouble C11 = 1.0;
612  Array<OneD, NekDouble > Fwd(nTraceNumPoints);
613  Array<OneD, NekDouble > Bwd(nTraceNumPoints);
614  Array<OneD, NekDouble > Vn (nTraceNumPoints, 0.0);
615 
616  Array<OneD, NekDouble > qFwd (nTraceNumPoints);
617  Array<OneD, NekDouble > qBwd (nTraceNumPoints);
618  Array<OneD, NekDouble > qfluxtemp(nTraceNumPoints, 0.0);
619 
620  Array<OneD, NekDouble > uterm(nTraceNumPoints);
621 
622  // Evaulate upwind flux:
623  // qflux = \hat{q} \cdot u = q \cdot n - C_(11)*(u^+ - u^-)
624  for (int i = 0; i < nvariables; ++i)
625  {
626  qflux[i] = Array<OneD, NekDouble> (nTraceNumPoints, 0.0);
627  for (int j = 0; j < nqvar; ++j)
628  {
629  // Compute Fwd and Bwd value of ufield of jth direction
630  m_fields[i]->GetFwdBwdTracePhys(qfield[j][i],qFwd,qBwd);
631 
632  // if Vn >= 0, flux = uFwd, i.e.,
633  // edge::eForward, if V*n>=0 <=> V*n_F>=0, pick
634  // qflux = qBwd = q+
635  // edge::eBackward, if V*n>=0 <=> V*n_B<0, pick
636  // qflux = qBwd = q-
637 
638  // else if Vn < 0, flux = uBwd, i.e.,
639  // edge::eForward, if V*n<0 <=> V*n_F<0, pick
640  // qflux = qFwd = q-
641  // edge::eBackward, if V*n<0 <=> V*n_B>=0, pick
642  // qflux = qFwd = q+
643 
644  m_fields[i]->GetTrace()->Upwind(m_traceNormals[j],
645  qBwd, qFwd,
646  qfluxtemp);
647 
648  Vmath::Vmul(nTraceNumPoints,
649  m_traceNormals[j], 1,
650  qfluxtemp, 1,
651  qfluxtemp, 1);
652 
653  // Generate Stability term = - C11 ( u- - u+ )
654  m_fields[i]->GetFwdBwdTracePhys(ufield[i], Fwd, Bwd);
655 
656  Vmath::Vsub(nTraceNumPoints,
657  Fwd, 1, Bwd, 1,
658  uterm, 1);
659 
660  Vmath::Smul(nTraceNumPoints,
661  -1.0 * C11, uterm, 1,
662  uterm, 1);
663 
664  // Flux = {Fwd, Bwd} * (nx, ny, nz) + uterm * (nx, ny)
665  Vmath::Vadd(nTraceNumPoints,
666  uterm, 1,
667  qfluxtemp, 1,
668  qfluxtemp, 1);
669 
670  // Imposing weak boundary condition with flux
671  if (m_fields[0]->GetBndCondExpansions().num_elements())
672  {
673  WeakPenaltyforVector(i, j,
674  qfield[j][i],
675  qfluxtemp,
676  C11);
677  }
678 
679  // q_hat \cdot n = (q_xi \cdot n_xi) or (q_eta \cdot n_eta)
680  // n_xi = n_x * tan_xi_x + n_y * tan_xi_y + n_z * tan_xi_z
681  // n_xi = n_x * tan_eta_x + n_y * tan_eta_y + n_z*tan_eta_z
682  Vmath::Vadd(nTraceNumPoints,
683  qfluxtemp, 1,
684  qflux[i], 1,
685  qflux[i], 1);
686  }
687  }
688  }
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:199
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:329
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:285
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:169
bool Nektar::SolverUtils::UnsteadySystem::v_PostIntegrate ( int  step)
protectedvirtual

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

Definition at line 921 of file UnsteadySystem.cpp.

Referenced by v_DoSolve().

922  {
923  return false;
924  }
bool Nektar::SolverUtils::UnsteadySystem::v_PreIntegrate ( int  step)
protectedvirtual

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

Definition at line 916 of file UnsteadySystem.cpp.

Referenced by v_DoSolve().

917  {
918  return false;
919  }
bool Nektar::SolverUtils::UnsteadySystem::v_SteadyStateCheck ( int  step)
protectedvirtual

Definition at line 926 of file UnsteadySystem.cpp.

927  {
928  return false;
929  }
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 738 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().

743  {
744  int i, e, npoints, id1, id2;
745 
746  // Number of boundary regions
747  int nbnd = m_fields[var]->GetBndCondExpansions().num_elements();
748  int Nfps, numBDEdge;
749  int nTraceNumPoints = GetTraceNpoints();
750  int cnt = 0;
751 
752  Array<OneD, NekDouble > uplus(nTraceNumPoints);
753 
754  m_fields[var]->ExtractTracePhys(physfield, uplus);
755  for (i = 0; i < nbnd; ++i)
756  {
757  // Number of boundary expansion related to that region
758  numBDEdge = m_fields[var]->
759  GetBndCondExpansions()[i]->GetExpSize();
760 
761  // Evaluate boundary values g_D or g_N from input files
763  m_session->GetFunction("InitialConditions", 0);
764 
765  npoints = m_fields[var]->
766  GetBndCondExpansions()[i]->GetNpoints();
767 
768  Array<OneD,NekDouble> BDphysics(npoints);
769  Array<OneD,NekDouble> x0(npoints,0.0);
770  Array<OneD,NekDouble> x1(npoints,0.0);
771  Array<OneD,NekDouble> x2(npoints,0.0);
772 
773  m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
774  ifunc->Evaluate(x0,x1,x2,time,BDphysics);
775 
776  // Weakly impose boundary conditions by modifying flux values
777  for (e = 0; e < numBDEdge ; ++e)
778  {
779  // Number of points on the expansion
780  Nfps = m_fields[var]->
781  GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
782 
783  id1 = m_fields[var]->
784  GetBndCondExpansions()[i]->GetPhys_Offset(e);
785 
786  id2 = m_fields[0]->GetTrace()->
787  GetPhys_Offset(m_fields[0]->GetTraceMap()->
788  GetBndCondTraceToGlobalTraceMap(cnt++));
789 
790  // For Dirichlet boundary condition: uflux = g_D
791  if (m_fields[var]->GetBndConditions()[i]->
792  GetBoundaryConditionType() == SpatialDomains::eDirichlet)
793  {
794  Vmath::Vcopy(Nfps,
795  &BDphysics[id1], 1,
796  &penaltyflux[id2], 1);
797  }
798  // For Neumann boundary condition: uflux = u+
799  else if ((m_fields[var]->GetBndConditions()[i])->
800  GetBoundaryConditionType() == SpatialDomains::eNeumann)
801  {
802  Vmath::Vcopy(Nfps,
803  &uplus[id2], 1,
804  &penaltyflux[id2], 1);
805  }
806  }
807  }
808  }
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:1047
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 815 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().

822  {
823  int i, e, npoints, id1, id2;
824  int nbnd = m_fields[var]->GetBndCondExpansions().num_elements();
825  int numBDEdge, Nfps;
826  int nTraceNumPoints = GetTraceNpoints();
827  Array<OneD, NekDouble > uterm(nTraceNumPoints);
828  Array<OneD, NekDouble > qtemp(nTraceNumPoints);
829  int cnt = 0;
830 
831  m_fields[var]->ExtractTracePhys(physfield,qtemp);
832 
833  for (i = 0; i < nbnd; ++i)
834  {
835  numBDEdge = m_fields[var]->
836  GetBndCondExpansions()[i]->GetExpSize();
837 
838  // Evaluate boundary values g_D or g_N from input files
840  m_session->GetFunction("InitialConditions", 0);
841 
842  npoints = m_fields[var]->
843  GetBndCondExpansions()[i]->GetNpoints();
844 
845  Array<OneD,NekDouble> BDphysics(npoints);
846  Array<OneD,NekDouble> x0(npoints,0.0);
847  Array<OneD,NekDouble> x1(npoints,0.0);
848  Array<OneD,NekDouble> x2(npoints,0.0);
849 
850  m_fields[var]->GetBndCondExpansions()[i]->GetCoords(x0,x1,x2);
851  ifunc->Evaluate(x0,x1,x2,time,BDphysics);
852 
853  // Weakly impose boundary conditions by modifying flux values
854  for (e = 0; e < numBDEdge ; ++e)
855  {
856  Nfps = m_fields[var]->
857  GetBndCondExpansions()[i]->GetExp(e)->GetNumPoints(0);
858 
859  id1 = m_fields[var]->
860  GetBndCondExpansions()[i]->GetPhys_Offset(e);
861 
862  id2 = m_fields[0]->GetTrace()->
863  GetPhys_Offset(m_fields[0]->GetTraceMap()->
864  GetBndCondTraceToGlobalTraceMap(cnt++));
865 
866  // For Dirichlet boundary condition:
867  //qflux = q+ - C_11 (u+ - g_D) (nx, ny)
868  if(m_fields[var]->GetBndConditions()[i]->
869  GetBoundaryConditionType() == SpatialDomains::eDirichlet)
870  {
871  Vmath::Vmul(Nfps,
872  &m_traceNormals[dir][id2], 1,
873  &qtemp[id2], 1,
874  &penaltyflux[id2], 1);
875  }
876  // For Neumann boundary condition: qflux = g_N
877  else if((m_fields[var]->GetBndConditions()[i])->
878  GetBoundaryConditionType() == SpatialDomains::eNeumann)
879  {
880  Vmath::Vmul(Nfps,
881  &m_traceNormals[dir][id2], 1,
882  &BDphysics[id1], 1,
883  &penaltyflux[id2], 1);
884  }
885  }
886  }
887  }
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:169

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(), and v_InitObject().

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

Definition at line 71 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 77 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
LibUtilities::TimeIntegrationSchemeOperators Nektar::SolverUtils::UnsteadySystem::m_ode
protected