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 string &s1, const 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 string &s1, const 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...
 
map< std::string, Array< OneD,
Array< OneD, float > > > 
m_interpWeights
 Map of the interpolation weights for a specific filename. More...
 
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...
 
int m_NumMode
 Mode to use in case of single mode analysis. 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 124 of file UnsteadySystem.cpp.

125  {
126  }
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 65 of file UnsteadySystem.cpp.

67  : EquationSystem(pSession),
68  m_infosteps(10)
69 
70  {
71  }
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 692 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().

693  {
694  if (m_session->DefinesFunction("InitialConditions"))
695  {
696  for (int i = 0; i < m_fields.num_elements(); ++i)
697  {
699 
700  vType = m_session->GetFunctionType(
701  "InitialConditions", m_session->GetVariable(i));
702 
703  if (vType == LibUtilities::eFunctionTypeFile)
704  {
705  std::string filename
706  = m_session->GetFunctionFilename(
707  "InitialConditions", m_session->GetVariable(i));
708 
709  fs::path pfilename(filename);
710 
711  // redefine path for parallel file which is in directory
712  if(fs::is_directory(pfilename))
713  {
714  fs::path metafile("Info.xml");
715  fs::path fullpath = pfilename / metafile;
716  filename = LibUtilities::PortablePath(fullpath);
717  }
718  m_fld->ImportFieldMetaData(filename, m_fieldMetaDataMap);
719 
720  // check to see if time defined
721  if (m_fieldMetaDataMap !=
723  {
725 
726  iter = m_fieldMetaDataMap.find("Time");
727  if (iter != m_fieldMetaDataMap.end())
728  {
729  time = boost::lexical_cast<NekDouble>(
730  iter->second);
731  }
732  }
733 
734  break;
735  }
736  }
737  }
738  }
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 899 of file UnsteadySystem.cpp.

References v_GetTimeStep().

901  {
902  return v_GetTimeStep(inarray);
903  }
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 131 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.

132  {
133  NekDouble TimeStability = 0.0;
134  switch(m_intScheme->GetIntegrationMethod())
135  {
139  {
140  TimeStability = 2.784;
141  break;
142  }
149  {
150  TimeStability = 2.0;
151  break;
152  }
154  {
155  TimeStability = 1.0;
156  break;
157  }
158  default:
159  {
160  ASSERTL0(
161  false,
162  "No CFL control implementation for this time"
163  "integration scheme");
164  }
165  }
166  return TimeStability;
167  }
#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 933 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().

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

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

Referenced by v_DoSolve().

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

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

452  {
456  }
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 173 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::Timer::Start(), Nektar::Timer::Stop(), Nektar::Timer::TimePerTest(), v_AppendOutput1D(), v_PostIntegrate(), and v_PreIntegrate().

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

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

References ASSERTL0.

Referenced by GetTimeStep().

913  {
914  ASSERTL0(false, "Not defined for this class");
915  return 0.0;
916  }
#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 76 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().

77  {
79 
80  m_initialStep = 0;
81 
82  // Load SolverInfo parameters
83  m_session->MatchSolverInfo("DIFFUSIONADVANCEMENT","Explicit",
84  m_explicitDiffusion,true);
85  m_session->MatchSolverInfo("ADVECTIONADVANCEMENT","Explicit",
86  m_explicitAdvection,true);
87  m_session->MatchSolverInfo("REACTIONADVANCEMENT", "Explicit",
88  m_explicitReaction, true);
89 
90  // For steady problems, we do not initialise the time integration
91  if (m_session->DefinesSolverInfo("TIMEINTEGRATIONMETHOD"))
92  {
94  CreateInstance(m_session->GetSolverInfo(
95  "TIMEINTEGRATIONMETHOD"));
96 
97  // Load generic input parameters
98  m_session->LoadParameter("IO_InfoSteps", m_infosteps, 0);
99  m_session->LoadParameter("CFL", m_cflSafetyFactor, 0.0);
100 
101  // Set up time to be dumped in field information
102  m_fieldMetaDataMap["Time"] =
103  boost::lexical_cast<std::string>(m_time);
104  }
105 
106  // By default attempt to forward transform initial condition.
107  m_homoInitialFwd = true;
108 
109  // Set up filters
110  LibUtilities::FilterMap::const_iterator x;
111  LibUtilities::FilterMap f = m_session->GetFilters();
112  for (x = f.begin(); x != f.end(); ++x)
113  {
114  m_filters.push_back(GetFilterFactory().CreateInstance(
115  x->first,
116  m_session,
117  x->second));
118  }
119  }
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 517 of file UnsteadySystem.cpp.

References ASSERTL0.

520  {
521  ASSERTL0(false,
522  "This function is not implemented for this equation.");
523  }
#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 525 of file UnsteadySystem.cpp.

References ASSERTL0.

529  {
530  ASSERTL0(false,
531  "This function is not implemented for this equation.");
532  }
#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 534 of file UnsteadySystem.cpp.

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

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

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

Referenced by v_DoSolve().

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

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

Definition at line 918 of file UnsteadySystem.cpp.

Referenced by v_DoSolve().

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

Definition at line 928 of file UnsteadySystem.cpp.

929  {
930  return false;
931  }
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 740 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().

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

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