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

#include <MMFAdvection.h>

Inheritance diagram for Nektar::MMFAdvection:
[legend]

Public Member Functions

virtual ~MMFAdvection ()
 Destructor. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::MMFSystem
SOLVER_UTILS_EXPORT MMFSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
virtual SOLVER_UTILS_EXPORT ~MMFSystem ()
 
SOLVER_UTILS_EXPORT void MMFInitObject (const Array< OneD, const Array< OneD, NekDouble >> &Anisotropy, const int TangentXelem=-1)
 
SOLVER_UTILS_EXPORT void CopyBoundaryTrace (const Array< OneD, const NekDouble > &Fwd, Array< OneD, NekDouble > &Bwd, const BoundaryCopyType BDCopyType, const int var=0, const std::string btype="NoUserDefined")
 
- Public Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
virtual SOLVER_UTILS_EXPORT ~UnsteadySystem ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep (const Array< OneD, const Array< OneD, NekDouble >> &inarray)
 Calculate the larger time-step mantaining the problem stable. More...
 
SOLVER_UTILS_EXPORT void SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
- 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 (bool DeclareField=true)
 Initialises the members of this object. More...
 
SOLVER_UTILS_EXPORT void DoInitialise ()
 Perform any initialisation necessary before solving the problem. More...
 
SOLVER_UTILS_EXPORT void DoSolve ()
 Solve the problem. More...
 
SOLVER_UTILS_EXPORT void TransCoeffToPhys ()
 Transform from coefficient to physical space. More...
 
SOLVER_UTILS_EXPORT void TransPhysToCoeff ()
 Transform from physical to coefficient space. More...
 
SOLVER_UTILS_EXPORT void Output ()
 Perform output operations after solve. More...
 
SOLVER_UTILS_EXPORT NekDouble LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Linf error computation. More...
 
SOLVER_UTILS_EXPORT std::string GetSessionName ()
 Get Session name. More...
 
template<class T >
std::shared_ptr< T > as ()
 
SOLVER_UTILS_EXPORT void ResetSessionName (std::string newname)
 Reset Session name. More...
 
SOLVER_UTILS_EXPORT LibUtilities::SessionReaderSharedPtr GetSession ()
 Get Session name. More...
 
SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr GetPressure ()
 Get pressure field if available. More...
 
SOLVER_UTILS_EXPORT void ExtraFldOutput (std::vector< Array< OneD, NekDouble >> &fieldcoeffs, std::vector< std::string > &variables)
 
SOLVER_UTILS_EXPORT void PrintSummary (std::ostream &out)
 Print a summary of parameters and solver characteristics. More...
 
SOLVER_UTILS_EXPORT void SetLambda (NekDouble lambda)
 Set parameter m_lambda. More...
 
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction (std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
 Get a SessionFunction by name. More...
 
SOLVER_UTILS_EXPORT void SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 Initialise the data in the dependent fields. More...
 
SOLVER_UTILS_EXPORT void EvaluateExactSolution (int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 Evaluates an exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln, bool Normalised=false)
 Compute the L2 error between fields and a given exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, bool Normalised=false)
 Compute the L2 error of the fields. More...
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleErrorExtraPoints (unsigned int field)
 Compute error (L2 and L_inf) over an larger set of quadrature points return [L2 Linf]. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n)
 Write checkpoint file of m_fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble >> &fieldcoeffs, std::vector< std::string > &variables)
 Write checkpoint file of custom data fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_BaseFlow (const int n)
 Write base flow file of m_fields. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname)
 Write field data to the given filename. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble >> &fieldcoeffs, std::vector< std::string > &variables)
 Write input fields to the given filename. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Input field data from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFldToMultiDomains (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int ndomains)
 Input field data from the given file to multiple domains. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, std::vector< std::string > &fieldStr, Array< OneD, Array< OneD, NekDouble >> &coeffs)
 Output a field. Input field data into array from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, MultiRegions::ExpListSharedPtr &pField, std::string &pFieldName)
 Output a field. Input field data into ExpList from the given file. More...
 
SOLVER_UTILS_EXPORT void SessionSummary (SummaryList &vSummary)
 Write out a session summary. More...
 
SOLVER_UTILS_EXPORT Array< OneD, MultiRegions::ExpListSharedPtr > & UpdateFields ()
 
SOLVER_UTILS_EXPORT LibUtilities::FieldMetaDataMapUpdateFieldMetaDataMap ()
 Get hold of FieldInfoMap so it can be updated. More...
 
SOLVER_UTILS_EXPORT NekDouble GetFinalTime ()
 Return final time. More...
 
SOLVER_UTILS_EXPORT int GetNcoeffs ()
 
SOLVER_UTILS_EXPORT int GetNcoeffs (const int eid)
 
SOLVER_UTILS_EXPORT int GetNumExpModes ()
 
SOLVER_UTILS_EXPORT const Array< OneD, int > GetNumExpModesPerExp ()
 
SOLVER_UTILS_EXPORT int GetNvariables ()
 
SOLVER_UTILS_EXPORT const std::string GetVariable (unsigned int i)
 
SOLVER_UTILS_EXPORT int GetTraceTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTraceNpoints ()
 
SOLVER_UTILS_EXPORT int GetExpSize ()
 
SOLVER_UTILS_EXPORT int GetPhys_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetCoeff_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTotPoints (int n)
 
SOLVER_UTILS_EXPORT int GetNpoints ()
 
SOLVER_UTILS_EXPORT int GetSteps ()
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep ()
 
SOLVER_UTILS_EXPORT void CopyFromPhysField (const int i, Array< OneD, NekDouble > &output)
 
SOLVER_UTILS_EXPORT void CopyToPhysField (const int i, Array< OneD, NekDouble > &output)
 
SOLVER_UTILS_EXPORT void SetSteps (const int steps)
 
SOLVER_UTILS_EXPORT void ZeroPhysFields ()
 
SOLVER_UTILS_EXPORT void FwdTransFields ()
 
SOLVER_UTILS_EXPORT void SetModifiedBasis (const bool modbasis)
 
SOLVER_UTILS_EXPORT int GetCheckpointNumber ()
 
SOLVER_UTILS_EXPORT void SetCheckpointNumber (int num)
 
SOLVER_UTILS_EXPORT int GetCheckpointSteps ()
 
SOLVER_UTILS_EXPORT void SetCheckpointSteps (int num)
 
SOLVER_UTILS_EXPORT Array< OneD, const Array< OneD, NekDouble > > GetTraceNormals ()
 
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 Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
SOLVER_UTILS_EXPORT AdvectionSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
virtual SOLVER_UTILS_EXPORT ~AdvectionSystem ()
 
SOLVER_UTILS_EXPORT AdvectionSharedPtr GetAdvObject ()
 Returns the advection object held by this instance. More...
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetElmtCFLVals (const bool FlagAcousticCFL=true)
 
SOLVER_UTILS_EXPORT NekDouble GetCFLEstimate (int &elmtid)
 

Static Public Member Functions

static SolverUtils::EquationSystemSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Creates an instance of this class. More...
 

Public Attributes

TestType m_TestType
 
- Public Attributes inherited from Nektar::SolverUtils::MMFSystem
NekDouble m_pi
 
int m_shapedim
 
SurfaceType m_surfaceType
 
UpwindType m_upwindType
 
TestMaxwellType m_TestMaxwellType
 
PolType m_PolType
 
IncType m_IncType
 
Array< OneD, NekDoublem_MMFfactors
 
- Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem
NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 
NekDouble m_cflNonAcoustic
 
NekDouble m_CFLGrowth
 CFL growth rate. More...
 
NekDouble m_CFLEnd
 maximun cfl in cfl growth More...
 

Static Public Attributes

static std::string className
 Name of class. More...
 
- Static Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem
static std::string cmdSetStartTime
 
static std::string cmdSetStartChkNum
 

Protected Member Functions

 MMFAdvection (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Session reader. More...
 
void WeakDGDirectionalAdvection (const Array< OneD, Array< OneD, NekDouble >> &InField, Array< OneD, Array< OneD, NekDouble >> &OutField)
 
void GetFluxVector (const Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &flux)
 Evaluate the flux at each solution point. More...
 
void GetFluxVectorDeAlias (const Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &flux)
 Evaluate the flux at each solution point using dealiasing. More...
 
void DoOdeRhs (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
 Compute the RHS. More...
 
void DoOdeProjection (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
 Compute the projection. More...
 
void EvaluateAdvectionVelocity (Array< OneD, Array< OneD, NekDouble >> &velocity)
 
NekDouble ComputeCirculatingArclength (const NekDouble zlevel, const NekDouble Rhs)
 
Array< OneD, NekDouble > & GetNormalVelocity ()
 Get the normal velocity. More...
 
void ComputeNablaCdotVelocity (Array< OneD, NekDouble > &vellc)
 
void ComputeveldotMF (Array< OneD, Array< OneD, NekDouble >> &veldotMF)
 
void AdvectionBellPlane (Array< OneD, NekDouble > &outfield)
 
void AdvectionBellSphere (Array< OneD, NekDouble > &outfield)
 
void Test2Dproblem (const NekDouble time, Array< OneD, NekDouble > &outfield)
 
void Test3Dproblem (const NekDouble time, Array< OneD, NekDouble > &outfield)
 
virtual void v_InitObject (bool DeclareFields=true)
 Initialise the object. More...
 
virtual void v_DoSolve ()
 Solves an unsteady problem. More...
 
virtual void v_GenerateSummary (SolverUtils::SummaryList &s)
 Print Summary. More...
 
virtual void v_SetInitialConditions (const NekDouble initialtime, bool dumpInitialConditions, const int domain)
 
virtual void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 
- Protected Member Functions inherited from Nektar::SolverUtils::MMFSystem
void SetUpMovingFrames (const Array< OneD, const Array< OneD, NekDouble >> &Anisotropy, const int TangentXelem)
 
void CheckMovingFrames (const Array< OneD, const Array< OneD, NekDouble >> &movingframes)
 
SOLVER_UTILS_EXPORT void ComputencdotMF ()
 
SOLVER_UTILS_EXPORT void ComputeDivCurlMF ()
 
SOLVER_UTILS_EXPORT void ComputeMFtrace ()
 
SOLVER_UTILS_EXPORT void VectorDotProd (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, NekDouble > &v3)
 
SOLVER_UTILS_EXPORT void VectorCrossProd (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &v3)
 
SOLVER_UTILS_EXPORT void VectorCrossProd (const Array< OneD, NekDouble > &v1, const Array< OneD, NekDouble > &v2, Array< OneD, NekDouble > &v3)
 
SOLVER_UTILS_EXPORT void ComputeCurl (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleCartesianToMovingframes (const Array< OneD, const Array< OneD, NekDouble >> &uvec, unsigned int field)
 
SOLVER_UTILS_EXPORT void DeriveCrossProductMF (Array< OneD, Array< OneD, NekDouble >> &CrossProductMF)
 
SOLVER_UTILS_EXPORT void ComputeNtimesMF ()
 
SOLVER_UTILS_EXPORT void ComputeNtimesFz (const int dir, const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &imFwd, const Array< OneD, const NekDouble > &imBwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimesF12 (const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &im1Fwd, const Array< OneD, const NekDouble > &im1Bwd, const Array< OneD, const NekDouble > &im2Fwd, const Array< OneD, const NekDouble > &im2Bwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimestimesdFz (const int dir, const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &imFwd, const Array< OneD, const NekDouble > &imBwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimestimesdF12 (const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &im1Fwd, const Array< OneD, const NekDouble > &im1Bwd, const Array< OneD, const NekDouble > &im2Fwd, const Array< OneD, const NekDouble > &im2Bwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void CartesianToSpherical (const NekDouble x0j, const NekDouble x1j, const NekDouble x2j, NekDouble &sin_varphi, NekDouble &cos_varphi, NekDouble &sin_theta, NekDouble &cos_theta)
 
SOLVER_UTILS_EXPORT void ComputeZimYim (Array< OneD, Array< OneD, NekDouble >> &epsvec, Array< OneD, Array< OneD, NekDouble >> &muvec)
 
SOLVER_UTILS_EXPORT void AdddedtMaxwell (const Array< OneD, const Array< OneD, NekDouble >> &physarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
 
SOLVER_UTILS_EXPORT void GetMaxwellFluxVector (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void GetMaxwellFlux1D (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void GetMaxwellFlux2D (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void LaxFriedrichMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void UpwindMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void AverageMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFlux (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time=0.0)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFluxTM (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFluxTE (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetIncidentField (const int var, const NekDouble time)
 
SOLVER_UTILS_EXPORT void Computedemdxicdote ()
 
SOLVER_UTILS_EXPORT NekDouble AvgInt (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble AvgAbsInt (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble AbsIntegral (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble RootMeanSquare (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble VectorAvgMagnitude (const Array< OneD, const Array< OneD, NekDouble >> &inarray)
 
SOLVER_UTILS_EXPORT void GramSchumitz (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &outarray, bool KeepTheMagnitude=true)
 
SOLVER_UTILS_EXPORT void BubbleSort (Array< OneD, NekDouble > &refarray, Array< OneD, NekDouble > &sortarray)
 
- Protected Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises UnsteadySystem class members. More...
 
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator ()
 Get the maximum timestep estimator for cfl control. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoInitialise ()
 Sets up initial conditions. More...
 
virtual SOLVER_UTILS_EXPORT void v_AppendOutput1D (Array< OneD, Array< OneD, NekDouble >> &solution1D)
 Print the solution at each solution point in a txt file. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_GetTimeStep (const Array< OneD, const Array< OneD, NekDouble >> &inarray)
 Return the timestep to be used for the next step in the time-marching loop. More...
 
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans ()
 
virtual SOLVER_UTILS_EXPORT void v_SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time, int &nchk)
 
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff (const Array< OneD, Array< OneD, NekDouble >> vel, StdRegions::VarCoeffMap &varCoeffMap)
 Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity. More...
 
virtual SOLVER_UTILS_EXPORT bool UpdateTimeStepCheck ()
 
- Protected Member Functions inherited from Nektar::SolverUtils::EquationSystem
SOLVER_UTILS_EXPORT EquationSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises EquationSystem class members. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Virtual function for the L_inf error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the L_2 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
virtual SOLVER_UTILS_EXPORT void v_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 Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT Array< OneD, NekDoublev_GetMaxStdVelocity (const NekDouble SpeedSoundFactor=1.0)
 

Protected Attributes

SolverUtils::RiemannSolverSharedPtr m_riemannSolver
 
NekDouble m_advx
 
NekDouble m_advy
 
NekDouble m_advz
 
NekDouble m_waveFreq
 
NekDouble m_RotAngle
 
NekDouble m_Mass0
 
int m_VelProjection
 
Array< OneD, Array< OneD, NekDouble > > m_velocity
 Advection velocity. More...
 
Array< OneD, NekDoublem_traceVn
 
Array< OneD, Array< OneD, NekDouble > > m_veldotMF
 
Array< OneD, NekDoublem_vellc
 
int m_planeNumber
 
- Protected Attributes inherited from Nektar::SolverUtils::MMFSystem
NekDouble m_alpha
 
NekDouble m_Incfreq
 
int m_SmoothFactor
 
NekDouble m_SFinit
 
Array< OneD, Array< OneD, NekDouble > > m_movingframes
 
Array< OneD, Array< OneD, NekDouble > > m_surfaceNormal
 
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFFwd
 
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFBwd
 
Array< OneD, Array< OneD, NekDouble > > m_nperpcdotMFFwd
 
Array< OneD, Array< OneD, NekDouble > > m_nperpcdotMFBwd
 
Array< OneD, Array< OneD, NekDouble > > m_DivMF
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_CurlMF
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_MFtraceFwd
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_MFtraceBwd
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimesMFFwd
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimesMFBwd
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimes_ntimesMFFwd
 
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimes_ntimesMFBwd
 
Array< OneD, Array< OneD, NekDouble > > m_ZimFwd
 
Array< OneD, Array< OneD, NekDouble > > m_ZimBwd
 
Array< OneD, Array< OneD, NekDouble > > m_YimFwd
 
Array< OneD, Array< OneD, NekDouble > > m_YimBwd
 
Array< OneD, Array< OneD, NekDouble > > m_epsvec
 
Array< OneD, Array< OneD, NekDouble > > m_muvec
 
Array< OneD, Array< OneD, NekDouble > > m_negepsvecminus1
 
Array< OneD, Array< OneD, NekDouble > > m_negmuvecminus1
 
Array< OneD, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > > m_dedxi_cdot_e
 
SpatialDomains::GeomMMF m_MMFdir
 
Array< OneD, NekDoublem_MFlength
 
- Protected Attributes inherited from Nektar::SolverUtils::UnsteadySystem
int m_infosteps
 Number of time steps between outputting status information. More...
 
int m_abortSteps
 Number of steps between checks for abort conditions. More...
 
int m_filtersInfosteps
 Number of time steps between outputting filters information. More...
 
int m_nanSteps
 
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
 Wrapper to the time integration scheme. More...
 
LibUtilities::TimeIntegrationSchemeOperators m_ode
 The time integration scheme operators to use. More...
 
NekDouble m_epsilon
 
bool m_explicitDiffusion
 Indicates if explicit or implicit treatment of diffusion is used. More...
 
bool m_explicitAdvection
 Indicates if explicit or implicit treatment of advection is used. More...
 
bool m_explicitReaction
 Indicates if explicit or implicit treatment of reaction is used. More...
 
bool m_homoInitialFwd
 Flag to determine if simulation should start in homogeneous forward transformed state. More...
 
NekDouble m_steadyStateTol
 Tolerance to which steady state should be evaluated at. More...
 
int m_steadyStateSteps
 Check for steady state at step interval. More...
 
NekDouble m_steadyStateRes = 1.0
 
NekDouble m_steadyStateRes0 = 1.0
 
Array< OneD, Array< OneD, NekDouble > > m_previousSolution
 Storage for previous solution for steady-state check. More...
 
std::ofstream m_errFile
 
std::vector< int > m_intVariables
 
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
 
NekDouble m_filterTimeWarning
 Number of time steps between outputting status information. More...
 
NekDouble m_TimeIntegLambda = 0.0
 coefff of spacial derivatives(rhs or m_F in GLM) in calculating the residual of the whole equation(used in unsteady time integrations) More...
 
bool m_flagImplicitItsStatistics
 
bool m_flagImplicitSolver = false
 
Array< OneD, NekDoublem_magnitdEstimat
 estimate the magnitude of each conserved varibles More...
 
Array< OneD, NekDoublem_locTimeStep
 local time step(notice only for jfnk other see m_cflSafetyFactor) More...
 
NekDouble m_inArrayNorm = -1.0
 
int m_TotLinItePerStep = 0
 
int m_StagesPerStep = 1
 
bool m_flagUpdatePreconMat
 
int m_maxLinItePerNewton
 
int m_TotNewtonIts = 0
 
int m_TotLinIts = 0
 
int m_TotImpStages = 0
 
bool m_CalcPhysicalAV = true
 flag to update artificial viscosity More...
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
bool m_verbose
 
bool m_root
 
LibUtilities::SessionReaderSharedPtr m_session
 The session reader. More...
 
std::map< std::string, SolverUtils::SessionFunctionSharedPtrm_sessionFunctions
 Map of known SessionFunctions. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 Field input/output. More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 Array holding all dependent variables. More...
 
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions
 Pointer to boundary conditions object. More...
 
SpatialDomains::MeshGraphSharedPtr m_graph
 Pointer to graph defining mesh. More...
 
std::string m_sessionName
 Name of the session. More...
 
NekDouble m_time
 Current time of simulation. More...
 
int m_initialStep
 Number of the step where the simulation should begin. More...
 
NekDouble m_fintime
 Finish time of the simulation. More...
 
NekDouble m_timestep
 Time step size. More...
 
NekDouble m_timestepMax = -1.0
 Time step size. More...
 
NekDouble m_lambda
 Lambda constant in real system if one required. More...
 
NekDouble m_checktime
 Time between checkpoints. More...
 
NekDouble m_lastCheckTime
 
NekDouble m_TimeIncrementFactor
 
int m_nchk
 Number of checkpoints written so far. More...
 
int m_steps
 Number of steps to take. More...
 
int m_checksteps
 Number of steps between checkpoints. More...
 
int m_spacedim
 Spatial dimension (>= expansion dim). More...
 
int m_expdim
 Expansion dimension. More...
 
bool m_singleMode
 Flag to determine if single homogeneous mode is used. More...
 
bool m_halfMode
 Flag to determine if half homogeneous mode is used. More...
 
bool m_multipleModes
 Flag to determine if use multiple homogenenous modes are used. More...
 
bool m_useFFT
 Flag to determine if FFT is used for homogeneous transform. More...
 
bool m_homogen_dealiasing
 Flag to determine if dealiasing is used for homogeneous simulations. More...
 
bool m_specHP_dealiasing
 Flag to determine if dealisising is usde for the Spectral/hp element discretisation. More...
 
enum MultiRegions::ProjectionType m_projectionType
 Type of projection; e.g continuous or discontinuous. More...
 
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
 Array holding trace normals for DG simulations in the forwards direction. More...
 
Array< OneD, bool > m_checkIfSystemSingular
 Flag to indicate if the fields should be checked for singularity. More...
 
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
 Map to identify relevant solver info to dump in output fields. More...
 
Array< OneD, NekDoublem_movingFrameVelsxyz
 Moving frame of reference velocities. More...
 
Array< OneD, NekDoublem_movingFrameTheta
 Moving frame of reference angles with respect to the. More...
 
boost::numeric::ublas::matrix< NekDoublem_movingFrameProjMat
 Projection matrix for transformation between inertial and moving. More...
 
int m_NumQuadPointsError
 Number of Quadrature points used to work out the error. More...
 
enum HomogeneousType m_HomogeneousType
 
NekDouble m_LhomX
 physical length in X direction (if homogeneous) More...
 
NekDouble m_LhomY
 physical length in Y direction (if homogeneous) More...
 
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous) More...
 
int m_npointsX
 number of points in X direction (if homogeneous) More...
 
int m_npointsY
 number of points in Y direction (if homogeneous) More...
 
int m_npointsZ
 number of points in Z direction (if homogeneous) More...
 
int m_HomoDirec
 number of homogenous directions More...
 
- Protected Attributes inherited from Nektar::SolverUtils::AdvectionSystem
SolverUtils::AdvectionSharedPtr m_advObject
 Advection term. More...
 

Friends

class MemoryManager< MMFAdvection >
 

Additional Inherited Members

- Protected Types inherited from Nektar::SolverUtils::EquationSystem
enum  HomogeneousType { eHomogeneous1D , eHomogeneous2D , eHomogeneous3D , eNotHomogeneous }
 Parameter for homogeneous expansions. More...
 
- Static Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
static std::string equationSystemTypeLookupIds []
 

Detailed Description

Definition at line 61 of file MMFAdvection.h.

Constructor & Destructor Documentation

◆ ~MMFAdvection()

Nektar::MMFAdvection::~MMFAdvection ( )
virtual

Destructor.

Unsteady linear advection equation destructor.

Definition at line 202 of file MMFAdvection.cpp.

203 {
204 }

◆ MMFAdvection()

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

Session reader.

Definition at line 51 of file MMFAdvection.cpp.

53  : UnsteadySystem(pSession, pGraph), MMFSystem(pSession, pGraph),
54  AdvectionSystem(pSession, pGraph)
55 {
56  m_planeNumber = 0;
57 }
SOLVER_UTILS_EXPORT AdvectionSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
SOLVER_UTILS_EXPORT MMFSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Definition: MMFSystem.cpp:43
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.

References m_planeNumber.

Member Function Documentation

◆ AdvectionBellPlane()

void Nektar::MMFAdvection::AdvectionBellPlane ( Array< OneD, NekDouble > &  outfield)
protected

Definition at line 1089 of file MMFAdvection.cpp.

1090 {
1091  int nq = m_fields[0]->GetNpoints();
1092 
1093  NekDouble dist, m_radius_limit;
1094 
1095  Array<OneD, NekDouble> x(nq);
1096  Array<OneD, NekDouble> y(nq);
1097  Array<OneD, NekDouble> z(nq);
1098 
1099  m_fields[0]->GetCoords(x, y, z);
1100 
1101  // Sets of parameters
1102  m_radius_limit = 0.5;
1103 
1104  NekDouble x0j, x1j;
1105  outfield = Array<OneD, NekDouble>(nq, 0.0);
1106  for (int j = 0; j < nq; ++j)
1107  {
1108  x0j = x[j];
1109  x1j = y[j];
1110 
1111  dist = sqrt(x0j * x0j + x1j * x1j);
1112 
1113  if (dist < m_radius_limit)
1114  {
1115  outfield[j] = 0.5 * (1.0 + cos(m_pi * dist / m_radius_limit));
1116  }
1117  else
1118  {
1119  outfield[j] = 0.0;
1120  }
1121  }
1122 }
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:291

References Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, and tinysimd::sqrt().

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ AdvectionBellSphere()

void Nektar::MMFAdvection::AdvectionBellSphere ( Array< OneD, NekDouble > &  outfield)
protected

Definition at line 1124 of file MMFAdvection.cpp.

1125 {
1126  int nq = m_fields[0]->GetNpoints();
1127 
1128  NekDouble dist, radius, cosdiff, sin_theta, cos_theta, sin_varphi,
1129  cos_varphi;
1130  NekDouble m_theta_c, m_varphi_c, m_radius_limit, m_c0;
1131 
1132  Array<OneD, NekDouble> x(nq);
1133  Array<OneD, NekDouble> y(nq);
1134  Array<OneD, NekDouble> z(nq);
1135 
1136  m_fields[0]->GetCoords(x, y, z);
1137 
1138  // Sets of parameters
1139  m_theta_c = 0.0;
1140  m_varphi_c = 3.0 * m_pi / 2.0;
1141  m_radius_limit = 7.0 * m_pi / 64.0;
1142  m_c0 = 0.0;
1143 
1144  NekDouble x0j, x1j, x2j;
1145  outfield = Array<OneD, NekDouble>(nq, 0.0);
1146  for (int j = 0; j < nq; ++j)
1147  {
1148  x0j = x[j];
1149  x1j = y[j];
1150  x2j = z[j];
1151 
1152  radius = sqrt(x0j * x0j + x1j * x1j + x2j * x2j);
1153 
1154  sin_varphi = x1j / sqrt(x0j * x0j + x1j * x1j);
1155  cos_varphi = x0j / sqrt(x0j * x0j + x1j * x1j);
1156 
1157  sin_theta = x2j / radius;
1158  cos_theta = sqrt(x0j * x0j + x1j * x1j) / radius;
1159 
1160  cosdiff = cos_varphi * cos(m_varphi_c) + sin_varphi * sin(m_varphi_c);
1161  dist = radius * acos(sin(m_theta_c) * sin_theta +
1162  cos(m_theta_c) * cos_theta * cosdiff);
1163 
1164  if (dist < m_radius_limit)
1165  {
1166  outfield[j] =
1167  0.5 * (1.0 + cos(m_pi * dist / m_radius_limit)) + m_c0;
1168  }
1169  else
1170  {
1171  outfield[j] = m_c0;
1172  }
1173  }
1174 }

References Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, and tinysimd::sqrt().

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ ComputeCirculatingArclength()

NekDouble Nektar::MMFAdvection::ComputeCirculatingArclength ( const NekDouble  zlevel,
const NekDouble  Rhs 
)
protected

Definition at line 864 of file MMFAdvection.cpp.

866 {
867 
868  NekDouble Tol = 0.0001, Maxiter = 1000, N = 100;
869  NekDouble newy, F = 0.0, dF = 1.0, y0, tmp;
870 
871  Array<OneD, NekDouble> xp(N + 1);
872  Array<OneD, NekDouble> yp(N + 1);
873 
874  NekDouble intval = 0.0;
875  switch (m_surfaceType)
876  {
879  {
880  intval = sqrt(Rhs - zlevel * zlevel);
881  }
882  break;
883 
885  {
886  intval = sqrt(0.5 * (Rhs - zlevel * zlevel * zlevel * zlevel -
887  zlevel * zlevel));
888  }
889  break;
890 
892  {
893  tmp = 0.5 *
894  (Rhs - zlevel * zlevel * zlevel * zlevel - zlevel * zlevel);
895  intval = sqrt(0.5 * (1.0 + sqrt(1.0 + 4.0 * tmp)));
896  }
897  break;
898 
899  default:
900  break;
901  }
902 
903  switch (m_surfaceType)
904  {
905  // Find the half of all the xp and yp on zlevel ....
909  {
910  for (int j = 0; j < N + 1; ++j)
911  {
912  xp[j] = j * 2.0 * intval / N - intval;
913 
914  y0 = 1.0;
915  for (int i = 0; i < Maxiter; ++i)
916  {
917  switch (m_surfaceType)
918  {
919  // Find the half of all the xp and yp on zlevel ....
922  {
923  F = xp[j] * xp[j] + y0 * y0 + zlevel * zlevel - Rhs;
924  dF = 2.0 * y0;
925  }
926  break;
927 
929  {
930  F = 2.0 * xp[j] * xp[j] + y0 * y0 * y0 * y0 +
931  y0 * y0 + zlevel * zlevel * zlevel * zlevel +
932  zlevel * zlevel - Rhs;
933  dF = 4.0 * y0 * y0 * y0 + 2.0 * y0;
934  }
935  break;
936 
937  default:
938  break;
939  }
940 
941  newy = y0 - F / dF;
942 
943  if (fabs(F / dF) < Tol)
944  {
945  yp[j] = newy;
946  break;
947  }
948 
949  else
950  {
951  y0 = newy;
952  }
953 
954  ASSERTL0(i < Maxiter,
955  "Advection Velocity convergence fails");
956 
957  } // i-loop
958  }
959  }
960  break;
961 
963  {
964  for (int j = 0; j < N + 1; ++j)
965  {
966  xp[j] = j * 2.0 * intval / N - intval;
967  tmp = 0.5 * Rhs -
968  0.5 * (zlevel * zlevel * zlevel * zlevel +
969  zlevel * zlevel) -
970  (xp[j] * xp[j] * xp[j] * xp[j] - xp[j] * xp[j]);
971  if (tmp < 0)
972  {
973  tmp = -1.0 * tmp;
974  }
975  yp[j] = sqrt(tmp);
976  } // j-loop
977  }
978  break;
979 
980  default:
981  break;
982 
983  } // switch-loop
984 
985  NekDouble pi = 3.14159265358979323846;
986  NekDouble arclength = 0.0;
987  for (int j = 0; j < N; ++j)
988  {
989  arclength =
990  arclength + sqrt((yp[j + 1] - yp[j]) * (yp[j + 1] - yp[j]) +
991  (xp[j + 1] - xp[j]) * (xp[j + 1] - xp[j])) /
992  pi;
993  }
994 
995  return arclength;
996 }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:215

References ASSERTL0, Nektar::SolverUtils::eIrregular, Nektar::SolverUtils::eNonconvex, Nektar::SolverUtils::eSphere, Nektar::SolverUtils::eTRSphere, Nektar::SolverUtils::MMFSystem::m_surfaceType, and tinysimd::sqrt().

◆ ComputeNablaCdotVelocity()

void Nektar::MMFAdvection::ComputeNablaCdotVelocity ( Array< OneD, NekDouble > &  vellc)
protected

Definition at line 1221 of file MMFAdvection.cpp.

1222 {
1223  int nq = m_fields[0]->GetNpoints();
1224 
1225  Array<OneD, NekDouble> velcoeff(nq, 0.0);
1226 
1227  Array<OneD, NekDouble> Dtmp0(nq);
1228  Array<OneD, NekDouble> Dtmp1(nq);
1229  Array<OneD, NekDouble> Dtmp2(nq);
1230  Array<OneD, NekDouble> Drv(nq);
1231 
1232  vellc = Array<OneD, NekDouble>(nq, 0.0);
1233 
1234  // m_vellc = \nabla m_vel \cdot tan_i
1235  Array<OneD, NekDouble> tmp(nq);
1236  Array<OneD, NekDouble> vessel(nq);
1237 
1238  for (int j = 0; j < m_shapedim; ++j)
1239  {
1240  Vmath::Zero(nq, velcoeff, 1);
1241  for (int k = 0; k < m_spacedim; ++k)
1242  {
1243  // a_j = tan_j cdot m_vel
1244  Vmath::Vvtvp(nq, &m_movingframes[j][k * nq], 1, &m_velocity[k][0],
1245  1, &velcoeff[0], 1, &velcoeff[0], 1);
1246  }
1247 
1248  // d a_j / d x^k
1249  m_fields[0]->PhysDeriv(velcoeff, Dtmp0, Dtmp1, Dtmp2);
1250 
1251  for (int k = 0; k < m_spacedim; ++k)
1252  {
1253  // tan_j^k ( d a_j / d x^k )
1254  switch (k)
1255  {
1256  case (0):
1257  {
1258  Vmath::Vvtvp(nq, &Dtmp0[0], 1, &m_movingframes[j][k * nq],
1259  1, &vellc[0], 1, &vellc[0], 1);
1260  }
1261  break;
1262 
1263  case (1):
1264  {
1265  Vmath::Vvtvp(nq, &Dtmp1[0], 1, &m_movingframes[j][k * nq],
1266  1, &vellc[0], 1, &vellc[0], 1);
1267  }
1268  break;
1269 
1270  case (2):
1271  {
1272  Vmath::Vvtvp(nq, &Dtmp2[0], 1, &m_movingframes[j][k * nq],
1273  1, &vellc[0], 1, &vellc[0], 1);
1274  }
1275  break;
1276  }
1277  }
1278  }
1279 }
Array< OneD, Array< OneD, NekDouble > > m_velocity
Advection velocity.
Definition: MMFAdvection.h:95
int m_spacedim
Spatial dimension (>= expansion dim).
Array< OneD, Array< OneD, NekDouble > > m_movingframes
Definition: MMFSystem.h:186
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:574
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492

References Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::MMFSystem::m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, m_velocity, Vmath::Vvtvp(), and Vmath::Zero().

Referenced by v_InitObject().

◆ ComputeveldotMF()

void Nektar::MMFAdvection::ComputeveldotMF ( Array< OneD, Array< OneD, NekDouble >> &  veldotMF)
protected

Definition at line 1281 of file MMFAdvection.cpp.

1283 {
1284  int nq = m_fields[0]->GetNpoints();
1285 
1286  veldotMF = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1287 
1288  Array<OneD, NekDouble> magMF(nq);
1289  for (int j = 0; j < m_shapedim; ++j)
1290  {
1291  veldotMF[j] = Array<OneD, NekDouble>(nq, 0.0);
1292 
1293  for (int k = 0; k < m_spacedim; ++k)
1294  {
1295  Vmath::Vvtvp(nq, &m_movingframes[j][k * nq], 1, &m_velocity[k][0],
1296  1, &veldotMF[j][0], 1, &veldotMF[j][0], 1);
1297  }
1298  }
1299 }

References Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::MMFSystem::m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, m_velocity, and Vmath::Vvtvp().

Referenced by v_InitObject().

◆ create()

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

Creates an instance of this class.

Definition at line 68 of file MMFAdvection.h.

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

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

◆ DoOdeProjection()

void Nektar::MMFAdvection::DoOdeProjection ( const Array< OneD, const Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const NekDouble  time 
)
protected

Compute the projection.

Compute the projection for the linear advection equation.

Parameters
inarrayGiven fields.
outarrayCalculated solution.
timeTime.

Definition at line 461 of file MMFAdvection.cpp.

464 {
465  // Counter variable
466  int i;
467 
468  // Number of fields (variables of the problem)
469  int nVariables = inarray.size();
470 
471  // Set the boundary conditions
472  SetBoundaryConditions(time);
473 
474  // Switch on the projection type (Discontinuous or Continuous)
475  switch (m_projectionType)
476  {
477  // Discontinuous projection
479  {
480  // Number of quadrature points
481  int nQuadraturePts = GetNpoints();
482 
483  // Just copy over array
484  for (i = 0; i < nVariables; ++i)
485  {
486  Vmath::Vcopy(nQuadraturePts, inarray[i], 1, outarray[i], 1);
487  }
488  break;
489  }
490 
491  // Continuous projection
494  {
495  Array<OneD, NekDouble> coeffs(m_fields[0]->GetNcoeffs(), 0.0);
496  for (i = 0; i < nVariables; ++i)
497  {
498  m_fields[i]->FwdTrans(inarray[i], coeffs);
499  m_fields[i]->BwdTrans(coeffs, outarray[i]);
500  }
501  break;
502  }
503 
504  default:
505  ASSERTL0(false, "Unknown projection scheme");
506  break;
507  }
508 }
SOLVER_UTILS_EXPORT int GetNpoints()
SOLVER_UTILS_EXPORT int GetNcoeffs()
enum MultiRegions::ProjectionType m_projectionType
Type of projection; e.g continuous or discontinuous.
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255

References ASSERTL0, Nektar::MultiRegions::eDiscontinuous, Nektar::MultiRegions::eGalerkin, Nektar::MultiRegions::eMixed_CG_Discontinuous, Nektar::SolverUtils::EquationSystem::GetNcoeffs(), Nektar::SolverUtils::EquationSystem::GetNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_projectionType, Nektar::SolverUtils::EquationSystem::SetBoundaryConditions(), and Vmath::Vcopy().

Referenced by v_InitObject().

◆ DoOdeRhs()

void Nektar::MMFAdvection::DoOdeRhs ( const Array< OneD, const Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
const NekDouble  time 
)
protected

Compute the RHS.

Compute the right-hand side for the linear advection equation.

Parameters
inarrayGiven fields.
outarrayCalculated solution.
timeTime.

Definition at line 393 of file MMFAdvection.cpp.

396 {
397  boost::ignore_unused(time);
398 
399  int i;
400  int nvariables = inarray.size();
401  int npoints = GetNpoints();
402 
403  switch (m_projectionType)
404  {
406  {
407  int ncoeffs = inarray[0].size();
408 
409  if (m_spacedim == 3)
410  {
411  Array<OneD, Array<OneD, NekDouble>> WeakAdv(nvariables);
412 
413  WeakAdv[0] = Array<OneD, NekDouble>(ncoeffs * nvariables);
414  for (i = 1; i < nvariables; ++i)
415  {
416  WeakAdv[i] = WeakAdv[i - 1] + ncoeffs;
417  }
418 
419  // Compute \nabla \cdot \vel u according to MMF scheme
420  WeakDGDirectionalAdvection(inarray, WeakAdv);
421 
422  for (i = 0; i < nvariables; ++i)
423  {
424  m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
425  m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
426 
427  // Add m_vellc * inarray[i] = \nabla v^m \cdot e^m to
428  // outarray[i]
429  // Vmath::Vvtvp(npoints, &m_vellc[0], 1, &inarray[i][0], 1,
430  // &outarray[i][0], 1, &outarray[i][0], 1);
431  Vmath::Neg(npoints, outarray[i], 1);
432  }
433  }
434  else
435  {
436  m_advObject->Advect(2, m_fields, m_velocity, inarray, outarray,
437  0.0);
438 
439  for (i = 0; i < nvariables; ++i)
440  {
441  Vmath::Neg(npoints, outarray[i], 1);
442  }
443  }
444  }
445  break;
446 
447  default:
448  {
449  ASSERTL0(false, "Unknown projection scheme");
450  }
451  break;
452  }
453 }
void WeakDGDirectionalAdvection(const Array< OneD, Array< OneD, NekDouble >> &InField, Array< OneD, Array< OneD, NekDouble >> &OutField)
SolverUtils::AdvectionSharedPtr m_advObject
Advection term.
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:518

References ASSERTL0, Nektar::MultiRegions::eDiscontinuous, Nektar::SolverUtils::EquationSystem::GetNpoints(), Nektar::SolverUtils::AdvectionSystem::m_advObject, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_projectionType, Nektar::SolverUtils::EquationSystem::m_spacedim, m_velocity, Vmath::Neg(), and WeakDGDirectionalAdvection().

Referenced by v_InitObject().

◆ EvaluateAdvectionVelocity()

void Nektar::MMFAdvection::EvaluateAdvectionVelocity ( Array< OneD, Array< OneD, NekDouble >> &  velocity)
protected

Definition at line 607 of file MMFAdvection.cpp.

609 {
610  int nq = m_fields[0]->GetNpoints();
611 
612  NekDouble vel_phi, vel_theta, sin_varphi, cos_varphi, sin_theta, cos_theta;
613  NekDouble x0j, x1j, x2j;
614 
615  Array<OneD, NekDouble> x0(nq);
616  Array<OneD, NekDouble> x1(nq);
617  Array<OneD, NekDouble> x2(nq);
618 
619  m_fields[0]->GetCoords(x0, x1, x2);
620 
621  // theta = a*sin(z/r), phi = a*tan(y/x);
622  for (int j = 0; j < nq; j++)
623  {
624  x0j = x0[j];
625  x1j = x1[j];
626  x2j = x2[j];
627 
628  CartesianToSpherical(x0j, x1j, x2j, sin_varphi, cos_varphi, sin_theta,
629  cos_theta);
630 
631  vel_phi = m_waveFreq * (cos_theta * cos(m_RotAngle) +
632  sin_theta * cos_varphi * sin(m_RotAngle));
633  vel_theta = -1.0 * m_waveFreq * sin_theta * sin(m_RotAngle);
634 
635  velocity[0][j] =
636  -1.0 * vel_phi * sin_varphi - vel_theta * sin_theta * cos_varphi;
637  velocity[1][j] =
638  vel_phi * cos_varphi - vel_theta * sin_theta * sin_varphi;
639  velocity[2][j] = vel_theta * cos_theta;
640  }
641 
642  // Project the veloicty on the tangent plane
643 
644  if (m_VelProjection)
645  {
646  // Check MovingFrames \cdot SurfaceNormals = 0
647  Array<OneD, Array<OneD, NekDouble>> newvelocity(m_spacedim);
648 
649  Array<OneD, Array<OneD, NekDouble>> MF1(m_spacedim);
650  Array<OneD, Array<OneD, NekDouble>> MF2(m_spacedim);
651  Array<OneD, Array<OneD, NekDouble>> SN(m_spacedim);
652 
653  for (int k = 0; k < m_spacedim; ++k)
654  {
655  newvelocity[k] = Array<OneD, NekDouble>(nq);
656  MF1[k] = Array<OneD, NekDouble>(nq);
657  MF2[k] = Array<OneD, NekDouble>(nq);
658  SN[k] = Array<OneD, NekDouble>(nq);
659 
660  Vmath::Vcopy(nq, &m_movingframes[0][k * nq], 1, &MF1[k][0], 1);
661  Vmath::Vcopy(nq, &m_movingframes[1][k * nq], 1, &MF2[k][0], 1);
662  }
663 
664  VectorCrossProd(MF1, MF2, SN);
665  GramSchumitz(SN, m_velocity, newvelocity, true);
666 
667  Array<OneD, NekDouble> tmp(nq, 0.0);
668  Array<OneD, NekDouble> Projection(nq, 0.0);
669 
670  for (int k = 0; k < m_spacedim; ++k)
671  {
672  Vmath::Vsub(nq, &m_velocity[0][0], 1, &newvelocity[0][0], 1,
673  &tmp[0], 1);
674  Vmath::Vmul(nq, &tmp[0], 1, &tmp[0], 1, &tmp[0], 1);
675  Vmath::Vadd(nq, tmp, 1, Projection, 1, Projection, 1);
676  }
677 
678  std::cout
679  << "Velocity vector is projected onto the tangent plane: Diff = "
680  << RootMeanSquare(Projection) << std::endl;
681 
682  for (int k = 0; k < m_spacedim; ++k)
683  {
684  Vmath::Vcopy(nq, &newvelocity[k][0], 1, &m_velocity[k][0], 1);
685  }
686  }
687 }
SOLVER_UTILS_EXPORT void GramSchumitz(const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &outarray, bool KeepTheMagnitude=true)
Definition: MMFSystem.cpp:2405
SOLVER_UTILS_EXPORT void VectorCrossProd(const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &v3)
Definition: MMFSystem.cpp:699
SOLVER_UTILS_EXPORT void CartesianToSpherical(const NekDouble x0j, const NekDouble x1j, const NekDouble x2j, NekDouble &sin_varphi, NekDouble &cos_varphi, NekDouble &sin_theta, NekDouble &cos_theta)
Definition: MMFSystem.cpp:795
SOLVER_UTILS_EXPORT NekDouble RootMeanSquare(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2344
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:209
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:359
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:419

References Nektar::SolverUtils::MMFSystem::CartesianToSpherical(), Nektar::SolverUtils::MMFSystem::GramSchumitz(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_movingframes, m_RotAngle, Nektar::SolverUtils::EquationSystem::m_spacedim, m_velocity, m_VelProjection, m_waveFreq, Nektar::SolverUtils::MMFSystem::RootMeanSquare(), Vmath::Vadd(), Vmath::Vcopy(), Nektar::SolverUtils::MMFSystem::VectorCrossProd(), Vmath::Vmul(), and Vmath::Vsub().

Referenced by v_InitObject().

◆ GetFluxVector()

void Nektar::MMFAdvection::GetFluxVector ( const Array< OneD, Array< OneD, NekDouble >> &  physfield,
Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &  flux 
)
protected

Evaluate the flux at each solution point.

Return the flux vector for the linear advection equation.

Parameters
iComponent of the flux vector to calculate.
physfieldFields.
fluxResulting flux.

Definition at line 517 of file MMFAdvection.cpp.

520 {
521  ASSERTL1(flux[0].size() == m_velocity.size(),
522  "Dimension of flux array and velocity array do not match");
523 
524  int i, j;
525  int nq = physfield[0].size();
526 
527  for (i = 0; i < flux.size(); ++i)
528  {
529  for (j = 0; j < flux[0].size(); ++j)
530  {
531  Vmath::Vmul(nq, physfield[i], 1, m_velocity[j], 1, flux[i][j], 1);
532  }
533  }
534 }
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:249

References ASSERTL1, m_velocity, and Vmath::Vmul().

Referenced by v_InitObject().

◆ GetFluxVectorDeAlias()

void Nektar::MMFAdvection::GetFluxVectorDeAlias ( const Array< OneD, Array< OneD, NekDouble >> &  physfield,
Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &  flux 
)
protected

Evaluate the flux at each solution point using dealiasing.

◆ GetNormalVelocity()

Array< OneD, NekDouble > & Nektar::MMFAdvection::GetNormalVelocity ( )
protected

Get the normal velocity.

Get the normal velocity for the linear advection equation.

Definition at line 363 of file MMFAdvection.cpp.

364 {
365  // Number of trace (interface) points
366  int i;
367  int nTracePts = GetTraceNpoints();
368 
369  // Auxiliary variable to compute the normal velocity
370  Array<OneD, NekDouble> tmp(nTracePts);
371 
372  // Reset the normal velocity
373  Vmath::Zero(nTracePts, m_traceVn, 1);
374 
375  for (i = 0; i < m_velocity.size(); ++i)
376  {
377  m_fields[0]->ExtractTracePhys(m_velocity[i], tmp);
378 
379  Vmath::Vvtvp(nTracePts, m_traceNormals[i], 1, tmp, 1, m_traceVn, 1,
380  m_traceVn, 1);
381  }
382 
383  return m_traceVn;
384 }
Array< OneD, NekDouble > m_traceVn
Definition: MMFAdvection.h:96
SOLVER_UTILS_EXPORT int GetTraceNpoints()
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.

References Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_traceNormals, m_traceVn, m_velocity, Vmath::Vvtvp(), and Vmath::Zero().

Referenced by v_InitObject().

◆ Test2Dproblem()

void Nektar::MMFAdvection::Test2Dproblem ( const NekDouble  time,
Array< OneD, NekDouble > &  outfield 
)
protected

Definition at line 1176 of file MMFAdvection.cpp.

1178 {
1179  int nq = m_fields[0]->GetNpoints();
1180 
1181  Array<OneD, NekDouble> x0(nq);
1182  Array<OneD, NekDouble> x1(nq);
1183  Array<OneD, NekDouble> x2(nq);
1184 
1185  m_fields[0]->GetCoords(x0, x1, x2);
1186 
1187  Array<OneD, NekDouble> u(nq);
1188 
1189  for (int i = 0; i < nq; ++i)
1190  {
1191  u[i] = cos(m_pi * (x0[i] - m_advx * time)) *
1192  cos(m_pi * (x1[i] - m_advy * time));
1193  }
1194 
1195  outfield = u;
1196 }

References m_advx, m_advy, Nektar::SolverUtils::EquationSystem::m_fields, and Nektar::SolverUtils::MMFSystem::m_pi.

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ Test3Dproblem()

void Nektar::MMFAdvection::Test3Dproblem ( const NekDouble  time,
Array< OneD, NekDouble > &  outfield 
)
protected

Definition at line 1198 of file MMFAdvection.cpp.

1200 {
1201  int nq = m_fields[0]->GetNpoints();
1202 
1203  Array<OneD, NekDouble> x0(nq);
1204  Array<OneD, NekDouble> x1(nq);
1205  Array<OneD, NekDouble> x2(nq);
1206 
1207  m_fields[0]->GetCoords(x0, x1, x2);
1208 
1209  Array<OneD, NekDouble> u(nq);
1210 
1211  for (int i = 0; i < nq; ++i)
1212  {
1213  u[i] = cos(m_pi * (x0[i] - m_advx * time)) *
1214  cos(m_pi * (x1[i] - m_advy * time)) *
1215  cos(m_pi * (x2[i] - m_advz * time));
1216  }
1217 
1218  outfield = u;
1219 }

References m_advx, m_advy, m_advz, Nektar::SolverUtils::EquationSystem::m_fields, and Nektar::SolverUtils::MMFSystem::m_pi.

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ v_DoSolve()

void Nektar::MMFAdvection::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::UnsteadySystem.

Definition at line 206 of file MMFAdvection.cpp.

207 {
208  ASSERTL0(m_intScheme != 0, "No time integration scheme.");
209 
210  int i, nchk = 1;
211  int nvariables = 0;
212  int nfields = m_fields.size();
213  int nq = m_fields[0]->GetNpoints();
214 
215  if (m_intVariables.empty())
216  {
217  for (i = 0; i < nfields; ++i)
218  {
219  m_intVariables.push_back(i);
220  }
221  nvariables = nfields;
222  }
223  else
224  {
225  nvariables = m_intVariables.size();
226  }
227 
228  // Set up wrapper to fields data storage.
229  Array<OneD, Array<OneD, NekDouble>> fields(nvariables);
230  Array<OneD, Array<OneD, NekDouble>> tmp(nvariables);
231 
232  // Order storage to list time-integrated fields first.
233  for (i = 0; i < nvariables; ++i)
234  {
235  fields[i] = m_fields[m_intVariables[i]]->GetPhys();
236  m_fields[m_intVariables[i]]->SetPhysState(false);
237  }
238 
239  // Initialise time integration scheme
240  m_intScheme->InitializeScheme(m_timestep, fields, m_time, m_ode);
241 
242  // Check uniqueness of checkpoint output
243  ASSERTL0((m_checktime == 0.0 && m_checksteps == 0) ||
244  (m_checktime > 0.0 && m_checksteps == 0) ||
245  (m_checktime == 0.0 && m_checksteps > 0),
246  "Only one of IO_CheckTime and IO_CheckSteps "
247  "should be set!");
248 
249  LibUtilities::Timer timer;
250  bool doCheckTime = false;
251  int step = 0;
252  NekDouble intTime = 0.0;
253  NekDouble cpuTime = 0.0;
254  NekDouble elapsed = 0.0;
255 
256  int Ntot, indx;
257  // Perform integration in time.
258  Ntot = m_steps / m_checksteps + 1;
259 
260  Array<OneD, NekDouble> dMass(Ntot);
261 
262  Array<OneD, NekDouble> zeta(nq);
263  Array<OneD, Array<OneD, NekDouble>> fieldsprimitive(nvariables);
264  for (int i = 0; i < nvariables; ++i)
265  {
266  fieldsprimitive[i] = Array<OneD, NekDouble>(nq);
267  }
268 
269  while (step < m_steps || m_time < m_fintime - NekConstants::kNekZeroTol)
270  {
271  timer.Start();
272  fields = m_intScheme->TimeIntegrate(step, m_timestep, m_ode);
273  timer.Stop();
274 
275  m_time += m_timestep;
276  elapsed = timer.TimePerTest(1);
277  intTime += elapsed;
278  cpuTime += elapsed;
279 
280  // Write out status information
281  if (m_session->GetComm()->GetRank() == 0 && !((step + 1) % m_infosteps))
282  {
283  std::cout << "Steps: " << std::setw(8) << std::left << step + 1
284  << " "
285  << "Time: " << std::setw(12) << std::left << m_time;
286 
287  std::stringstream ss;
288  ss << cpuTime << "s";
289  std::cout << " CPU Time: " << std::setw(8) << std::left << ss.str()
290  << std::endl;
291 
292  // Masss = h^*
293  indx = (step + 1) / m_checksteps;
294  dMass[indx] =
295  (m_fields[0]->Integral(fields[0]) - m_Mass0) / m_Mass0;
296 
297  std::cout << "dMass = " << std::setw(8) << std::left << dMass[indx]
298  << std::endl;
299 
300  cpuTime = 0.0;
301  }
302 
303  // Transform data into coefficient space
304  for (i = 0; i < nvariables; ++i)
305  {
306  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
307  m_fields[m_intVariables[i]]->FwdTransLocalElmt(
308  fields[i], m_fields[m_intVariables[i]]->UpdateCoeffs());
309  m_fields[m_intVariables[i]]->SetPhysState(false);
310  }
311 
312  // Write out checkpoint files
313  if ((m_checksteps && step && !((step + 1) % m_checksteps)) ||
314  doCheckTime)
315  {
316  Checkpoint_Output(nchk++);
317  doCheckTime = false;
318  }
319 
320  // Step advance
321  ++step;
322  }
323 
324  std::cout << "dMass = ";
325  for (i = 0; i < Ntot; ++i)
326  {
327  std::cout << dMass[i] << " , ";
328  }
329  std::cout << std::endl << std::endl;
330 
331  // Print out summary statistics
332  if (m_session->GetComm()->GetRank() == 0)
333  {
334  if (m_cflSafetyFactor > 0.0)
335  {
336  std::cout << "CFL safety factor : " << m_cflSafetyFactor
337  << std::endl
338  << "CFL time-step : " << m_timestep << std::endl;
339  }
340 
341  if (m_session->GetSolverInfo("Driver") != "SteadyState")
342  {
343  std::cout << "Time-integration : " << intTime << "s" << std::endl;
344  }
345  }
346 
347  for (i = 0; i < nvariables; ++i)
348  {
349  m_fields[m_intVariables[i]]->SetPhys(fields[i]);
350  m_fields[m_intVariables[i]]->SetPhysState(true);
351  }
352 
353  for (i = 0; i < nvariables; ++i)
354  {
355  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
356  m_fields[i]->UpdateCoeffs());
357  }
358 }
NekDouble m_timestep
Time step size.
NekDouble m_time
Current time of simulation.
NekDouble m_fintime
Finish time of the simulation.
SOLVER_UTILS_EXPORT void Checkpoint_Output(const int n)
Write checkpoint file of m_fields.
NekDouble m_checktime
Time between checkpoints.
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
int m_steps
Number of steps to take.
int m_checksteps
Number of steps between checkpoints.
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
NekDouble m_cflSafetyFactor
CFL safety factor (comprise between 0 to 1).
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
Wrapper to the time integration scheme.
int m_infosteps
Number of time steps between outputting status information.
static const NekDouble kNekZeroTol

References ASSERTL0, Nektar::SolverUtils::EquationSystem::Checkpoint_Output(), Nektar::NekConstants::kNekZeroTol, Nektar::SolverUtils::UnsteadySystem::m_cflSafetyFactor, Nektar::SolverUtils::EquationSystem::m_checksteps, Nektar::SolverUtils::EquationSystem::m_checktime, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_fintime, Nektar::SolverUtils::UnsteadySystem::m_infosteps, Nektar::SolverUtils::UnsteadySystem::m_intScheme, Nektar::SolverUtils::UnsteadySystem::m_intVariables, m_Mass0, Nektar::SolverUtils::UnsteadySystem::m_ode, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_steps, Nektar::SolverUtils::EquationSystem::m_time, Nektar::SolverUtils::EquationSystem::m_timestep, Nektar::LibUtilities::Timer::Start(), Nektar::LibUtilities::Timer::Stop(), and Nektar::LibUtilities::Timer::TimePerTest().

◆ v_EvaluateExactSolution()

void Nektar::MMFAdvection::v_EvaluateExactSolution ( unsigned int  field,
Array< OneD, NekDouble > &  outfield,
const NekDouble  time 
)
protectedvirtual

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 1301 of file MMFAdvection.cpp.

1304 {
1305  boost::ignore_unused(field);
1306 
1307  switch (m_TestType)
1308  {
1309  case eAdvectionBell:
1310  {
1311  AdvectionBellSphere(outfield);
1312  }
1313  break;
1314 
1315  case eTestPlane:
1316  {
1317  Test2Dproblem(time, outfield);
1318  }
1319  break;
1320 
1321  case eTestPlaneMassConsv:
1322  {
1323  AdvectionBellPlane(outfield);
1324  }
1325  break;
1326 
1327  case eTestCube:
1328  {
1329  Test3Dproblem(time, outfield);
1330  }
1331  break;
1332 
1333  default:
1334  break;
1335  }
1336 }
@ eTestPlaneMassConsv
Definition: MMFAdvection.h:46
@ eAdvectionBell
Definition: MMFAdvection.h:48
void Test2Dproblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
void Test3Dproblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
void AdvectionBellSphere(Array< OneD, NekDouble > &outfield)
void AdvectionBellPlane(Array< OneD, NekDouble > &outfield)
@ eTestPlane
Definition: MMFDiffusion.h:48
@ eTestCube
Definition: MMFDiffusion.h:49

References AdvectionBellPlane(), AdvectionBellSphere(), eAdvectionBell, Nektar::eTestCube, Nektar::eTestPlane, eTestPlaneMassConsv, m_TestType, Test2Dproblem(), and Test3Dproblem().

◆ v_GenerateSummary()

void Nektar::MMFAdvection::v_GenerateSummary ( SolverUtils::SummaryList s)
protectedvirtual

Print Summary.

Reimplemented from Nektar::SolverUtils::MMFSystem.

Definition at line 1338 of file MMFAdvection.cpp.

1339 {
1342  SolverUtils::AddSummaryItem(s, "Rotation Angle", m_RotAngle);
1343 }
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s)
Print a summary of time stepping parameters.
Definition: MMFSystem.cpp:2469
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:49
const char *const TestTypeMap[]
Definition: MMFDiffusion.h:58

References Nektar::SolverUtils::AddSummaryItem(), m_RotAngle, m_TestType, Nektar::TestTypeMap, and Nektar::SolverUtils::MMFSystem::v_GenerateSummary().

◆ v_InitObject()

void Nektar::MMFAdvection::v_InitObject ( bool  DeclareFields = true)
protectedvirtual

Initialise the object.

Initialisation object for the unsteady linear advection equation.

Reimplemented from Nektar::SolverUtils::AdvectionSystem.

Definition at line 62 of file MMFAdvection.cpp.

63 {
64  // Call to the initialisation object
65  UnsteadySystem::v_InitObject(DeclareFields);
66 
67  int nq = m_fields[0]->GetNpoints();
68  int shapedim = m_fields[0]->GetShapeDimension();
69  Array<OneD, Array<OneD, NekDouble>> Anisotropy(shapedim);
70  for (int j = 0; j < shapedim; ++j)
71  {
72  Anisotropy[j] = Array<OneD, NekDouble>(nq, 1.0);
73  }
74 
75  MMFSystem::MMFInitObject(Anisotropy);
76 
77  // Define TestType
78  ASSERTL0(m_session->DefinesSolverInfo("TESTTYPE"),
79  "No TESTTYPE defined in session.");
80  std::string TestTypeStr = m_session->GetSolverInfo("TESTTYPE");
81  for (int i = 0; i < (int)SIZE_TestType; ++i)
82  {
83  if (boost::iequals(TestTypeMap[i], TestTypeStr))
84  {
85  m_TestType = (TestType)i;
86  break;
87  }
88  }
89 
90  m_session->LoadParameter("Angular Frequency", m_waveFreq, m_pi);
91  m_session->LoadParameter("Rotational Angle", m_RotAngle, 0.0);
92  m_session->LoadParameter("Velocity Projection", m_VelProjection, 0);
93 
94  // Read the advection velocities from session file
95  m_session->LoadParameter("advx", m_advx, 1.0);
96  m_session->LoadParameter("advy", m_advy, 1.0);
97  m_session->LoadParameter("advz", m_advz, 1.0);
98 
99  std::vector<std::string> vel;
100  vel.push_back("Vx");
101  vel.push_back("Vy");
102  vel.push_back("Vz");
103 
104  // Resize the advection velocities vector to dimension of the problem
105  vel.resize(m_spacedim);
106 
107  // Store in the global variable m_velocity the advection velocities
108 
109  m_velocity = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
110  for (int k = 0; k < m_spacedim; ++k)
111  {
112  m_velocity[k] = Array<OneD, NekDouble>(nq);
113  }
114 
115  switch (m_surfaceType)
116  {
121  {
122  // true = project velocity onto the tangent plane
124  }
125  break;
126 
127  case SolverUtils::ePlane:
128  case SolverUtils::eCube:
129  {
130  GetFunction("AdvectionVelocity")->Evaluate(vel, m_velocity);
131  }
132  break;
133 
134  default:
135  break;
136  }
137 
138  std::cout << "|Velocity vector| = ( " << RootMeanSquare(m_velocity[0])
139  << " , " << RootMeanSquare(m_velocity[1]) << " , "
140  << RootMeanSquare(m_velocity[2]) << " ) " << std::endl;
141 
142  // Define the normal velocity fields
143  if (m_fields[0]->GetTrace())
144  {
145  m_traceVn = Array<OneD, NekDouble>(GetTraceNpoints());
146  }
147 
148  std::string advName;
149  std::string riemName;
150  m_session->LoadSolverInfo("AdvectionType", advName, "WeakDG");
151  m_advObject =
153  m_advObject->SetFluxVector(&MMFAdvection::GetFluxVector, this);
154  m_session->LoadSolverInfo("UpwindType", riemName, "Upwind");
156  riemName, m_session);
157  m_riemannSolver->SetScalar("Vn", &MMFAdvection::GetNormalVelocity, this);
158 
159  m_advObject->SetRiemannSolver(m_riemannSolver);
160  m_advObject->InitObject(m_session, m_fields);
161 
162  // Compute m_traceVn = n \cdot v
164 
165  // Compute m_vellc = nabal a^j \cdot m_vel
167  std::cout << "m_vellc is generated with mag = " << AvgInt(m_vellc)
168  << std::endl;
169 
170  // Compute vel \cdot MF
172 
173  // Modify e^i as v^i e^i
174  for (int j = 0; j < m_shapedim; j++)
175  {
176  for (int k = 0; k < m_spacedim; k++)
177  {
178  Vmath::Vmul(nq, &m_veldotMF[j][0], 1, &m_movingframes[j][k * nq], 1,
179  &m_movingframes[j][k * nq], 1);
180  }
181  }
182 
183  // Reflect it into m_ncdotMFFwd and Bwd
184  ComputencdotMF();
185 
186  // If explicit it computes RHS and PROJECTION for the time integration
188  {
191  }
192  // Otherwise it gives an error (no implicit integration)
193  else
194  {
195  ASSERTL0(false, "Implicit unsteady Advection not set up.");
196  }
197 }
tBaseSharedPtr CreateInstance(tKey idKey, tParam... args)
Create an instance of the class referred to by idKey.
Definition: NekFactory.hpp:144
void DefineProjection(FuncPointerT func, ObjectPointerT obj)
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
Array< OneD, Array< OneD, NekDouble > > m_veldotMF
Definition: MMFAdvection.h:98
void ComputeNablaCdotVelocity(Array< OneD, NekDouble > &vellc)
void EvaluateAdvectionVelocity(Array< OneD, Array< OneD, NekDouble >> &velocity)
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
Compute the RHS.
void DoOdeProjection(const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
Compute the projection.
void GetFluxVector(const Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, Array< OneD, NekDouble >>> &flux)
Evaluate the flux at each solution point.
void ComputeveldotMF(Array< OneD, Array< OneD, NekDouble >> &veldotMF)
SolverUtils::RiemannSolverSharedPtr m_riemannSolver
Definition: MMFAdvection.h:86
Array< OneD, NekDouble > & GetNormalVelocity()
Get the normal velocity.
Array< OneD, NekDouble > m_vellc
Definition: MMFAdvection.h:99
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction(std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
Get a SessionFunction by name.
SOLVER_UTILS_EXPORT NekDouble AvgInt(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2298
SOLVER_UTILS_EXPORT void MMFInitObject(const Array< OneD, const Array< OneD, NekDouble >> &Anisotropy, const int TangentXelem=-1)
Definition: MMFSystem.cpp:53
SOLVER_UTILS_EXPORT void ComputencdotMF()
Definition: MMFSystem.cpp:320
bool m_explicitAdvection
Indicates if explicit or implicit treatment of advection is used.
virtual SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true)
Init object for UnsteadySystem class.
AdvectionFactory & GetAdvectionFactory()
Gets the factory for initialising advection objects.
Definition: Advection.cpp:47
RiemannSolverFactory & GetRiemannSolverFactory()
@ SIZE_TestType
Length of enum list.
Definition: MMFDiffusion.h:55

References ASSERTL0, Nektar::SolverUtils::MMFSystem::AvgInt(), ComputeNablaCdotVelocity(), Nektar::SolverUtils::MMFSystem::ComputencdotMF(), ComputeveldotMF(), Nektar::LibUtilities::NekFactory< tKey, tBase, tParam >::CreateInstance(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineProjection(), DoOdeProjection(), DoOdeRhs(), Nektar::SolverUtils::eCube, Nektar::SolverUtils::eIrregular, Nektar::SolverUtils::eNonconvex, Nektar::SolverUtils::ePlane, Nektar::SolverUtils::eSphere, Nektar::SolverUtils::eTRSphere, EvaluateAdvectionVelocity(), Nektar::SolverUtils::GetAdvectionFactory(), GetFluxVector(), Nektar::SolverUtils::EquationSystem::GetFunction(), GetNormalVelocity(), Nektar::SolverUtils::GetRiemannSolverFactory(), Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::AdvectionSystem::m_advObject, m_advx, m_advy, m_advz, Nektar::SolverUtils::UnsteadySystem::m_explicitAdvection, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::UnsteadySystem::m_ode, Nektar::SolverUtils::MMFSystem::m_pi, m_riemannSolver, m_RotAngle, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::MMFSystem::m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, Nektar::SolverUtils::MMFSystem::m_surfaceType, m_TestType, m_traceVn, m_veldotMF, m_vellc, m_velocity, m_VelProjection, m_waveFreq, Nektar::SolverUtils::MMFSystem::MMFInitObject(), Nektar::SolverUtils::MMFSystem::RootMeanSquare(), Nektar::SIZE_TestType, Nektar::TestTypeMap, Nektar::SolverUtils::UnsteadySystem::v_InitObject(), and Vmath::Vmul().

◆ v_SetInitialConditions()

void Nektar::MMFAdvection::v_SetInitialConditions ( const NekDouble  initialtime,
bool  dumpInitialConditions,
const int  domain 
)
protectedvirtual

Set the physical fields based on a restart file, or a function describing the initial condition given in the session.

Parameters
initialtimeTime at which to evaluate the function.
dumpInitialConditionsWrite the initial condition to file?

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 998 of file MMFAdvection.cpp.

1001 {
1002  boost::ignore_unused(domain);
1003 
1004  int nq = m_fields[0]->GetNpoints();
1005 
1006  Array<OneD, NekDouble> u(nq);
1007 
1008  switch (m_TestType)
1009  {
1010  case eAdvectionBell:
1011  {
1013  m_fields[0]->SetPhys(u);
1014 
1015  m_Mass0 = m_fields[0]->Integral(u);
1016 
1017  // forward transform to fill the modal coeffs
1018  for (int i = 0; i < m_fields.size(); ++i)
1019  {
1020  m_fields[i]->SetPhysState(true);
1021  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
1022  m_fields[i]->UpdateCoeffs());
1023  }
1024  }
1025  break;
1026 
1027  case eTestPlane:
1028  {
1029  Test2Dproblem(initialtime, u);
1030  m_fields[0]->SetPhys(u);
1031 
1032  m_Mass0 = m_fields[0]->Integral(u);
1033 
1034  // forward transform to fill the modal coeffs
1035  for (int i = 0; i < m_fields.size(); ++i)
1036  {
1037  m_fields[i]->SetPhysState(true);
1038  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
1039  m_fields[i]->UpdateCoeffs());
1040  }
1041  }
1042  break;
1043 
1044  case eTestPlaneMassConsv:
1045  {
1046  AdvectionBellPlane(u);
1047  m_fields[0]->SetPhys(u);
1048 
1049  m_Mass0 = m_fields[0]->Integral(u);
1050  std::cout << "m_Mass0 = " << m_Mass0 << std::endl;
1051 
1052  // forward transform to fill the modal coeffs
1053  for (int i = 0; i < m_fields.size(); ++i)
1054  {
1055  m_fields[i]->SetPhysState(true);
1056  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
1057  m_fields[i]->UpdateCoeffs());
1058  }
1059  }
1060  break;
1061 
1062  case eTestCube:
1063  {
1064  Test3Dproblem(initialtime, u);
1065  m_fields[0]->SetPhys(u);
1066 
1067  // forward transform to fill the modal coeffs
1068  for (int i = 0; i < m_fields.size(); ++i)
1069  {
1070  m_fields[i]->SetPhysState(true);
1071  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
1072  m_fields[i]->UpdateCoeffs());
1073  }
1074  }
1075  break;
1076 
1077  default:
1078  break;
1079  }
1080 
1081  if (dumpInitialConditions)
1082  {
1083  // dump initial conditions to file
1084  std::string outname = m_sessionName + "_initial.chk";
1085  WriteFld(outname);
1086  }
1087 }
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
std::string m_sessionName
Name of the session.

References AdvectionBellPlane(), AdvectionBellSphere(), eAdvectionBell, Nektar::eTestCube, Nektar::eTestPlane, eTestPlaneMassConsv, Nektar::SolverUtils::EquationSystem::m_fields, m_Mass0, Nektar::SolverUtils::EquationSystem::m_sessionName, m_TestType, Test2Dproblem(), Test3Dproblem(), and Nektar::SolverUtils::EquationSystem::WriteFld().

◆ WeakDGDirectionalAdvection()

void Nektar::MMFAdvection::WeakDGDirectionalAdvection ( const Array< OneD, Array< OneD, NekDouble >> &  InField,
Array< OneD, Array< OneD, NekDouble >> &  OutField 
)
protected

Definition at line 536 of file MMFAdvection.cpp.

539 {
540  int i, j;
541  int ncoeffs = GetNcoeffs();
542  int nTracePointsTot = GetTraceNpoints();
543  int nvariables = m_fields.size();
544 
545  Array<OneD, Array<OneD, NekDouble>> physfield(nvariables);
546 
547  // Get the variables in physical space
548  // already in physical space
549  for (i = 0; i < nvariables; ++i)
550  {
551  physfield[i] = InField[i];
552  }
553 
554  Array<OneD, Array<OneD, NekDouble>> WeakDeriv(m_shapedim);
555  for (i = 0; i < nvariables; ++i)
556  {
557  for (j = 0; j < m_shapedim; ++j)
558  {
559  WeakDeriv[j] = Array<OneD, NekDouble>(ncoeffs, 0.0);
560 
561  // Directional derivation with respect to the j'th moving frame
562  // tmp[j] = \nabla \physfield[i] \cdot \mathbf{e}^j
563  // Implemented at TriExp::v_IProductWRTDirectionalDerivBase_SumFa
564  m_fields[0]->IProductWRTDirectionalDerivBase(
565  m_movingframes[j], physfield[i], WeakDeriv[j]);
566  }
567 
568  // Get the numerical flux and add to the modal coeffs
569  // if the NumericalFluxs function already includes the
570  // normal in the output
571 
572  Array<OneD, NekDouble> Fwd(nTracePointsTot);
573  Array<OneD, NekDouble> Bwd(nTracePointsTot);
574 
575  Array<OneD, NekDouble> flux(nTracePointsTot, 0.0);
576  Array<OneD, NekDouble> fluxFwd(nTracePointsTot);
577  Array<OneD, NekDouble> fluxBwd(nTracePointsTot);
578 
579  // Evaluate numerical flux in physical space which may in
580  // general couple all component of vectors
581  m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd, Bwd);
582 
583  // evaulate upwinded m_fields[i]
584  m_fields[i]->GetTrace()->Upwind(m_traceVn, Fwd, Bwd, flux);
585 
586  OutField[i] = Array<OneD, NekDouble>(ncoeffs, 0.0);
587  for (j = 0; j < m_shapedim; ++j)
588  {
589  // calculate numflux = (n \cdot MF)*flux
590  Vmath::Vmul(nTracePointsTot, &flux[0], 1, &m_ncdotMFFwd[j][0], 1,
591  &fluxFwd[0], 1);
592  Vmath::Vmul(nTracePointsTot, &flux[0], 1, &m_ncdotMFBwd[j][0], 1,
593  &fluxBwd[0], 1);
594  Vmath::Neg(ncoeffs, WeakDeriv[j], 1);
595 
596  // FwdBwdtegral because generallize (N \cdot MF)_{FWD} \neq -(N
597  // \cdot MF)_{BWD}
598  m_fields[i]->AddFwdBwdTraceIntegral(fluxFwd, fluxBwd, WeakDeriv[j]);
599  m_fields[i]->SetPhysState(false);
600 
601  Vmath::Vadd(ncoeffs, &WeakDeriv[j][0], 1, &OutField[i][0], 1,
602  &OutField[i][0], 1);
603  }
604  }
605 }
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFBwd
Definition: MMFSystem.h:190
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFFwd
Definition: MMFSystem.h:189

References Nektar::SolverUtils::EquationSystem::GetNcoeffs(), Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::MMFSystem::m_ncdotMFBwd, Nektar::SolverUtils::MMFSystem::m_ncdotMFFwd, Nektar::SolverUtils::MMFSystem::m_shapedim, m_traceVn, Vmath::Neg(), Vmath::Vadd(), and Vmath::Vmul().

Referenced by DoOdeRhs().

Friends And Related Function Documentation

◆ MemoryManager< MMFAdvection >

friend class MemoryManager< MMFAdvection >
friend

Definition at line 52 of file MMFAdvection.h.

Member Data Documentation

◆ className

std::string Nektar::MMFAdvection::className
static
Initial value:
=
"MMFAdvection", MMFAdvection::create, "MMFAdvection equation.")
tKey RegisterCreatorFunction(tKey idKey, CreatorFunction classCreator, std::string pDesc="")
Register a class with the factory.
Definition: NekFactory.hpp:198
static SolverUtils::EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
Definition: MMFAdvection.h:68
EquationSystemFactory & GetEquationSystemFactory()

Name of class.

Definition at line 78 of file MMFAdvection.h.

◆ m_advx

NekDouble Nektar::MMFAdvection::m_advx
protected

Definition at line 88 of file MMFAdvection.h.

Referenced by Test2Dproblem(), Test3Dproblem(), and v_InitObject().

◆ m_advy

NekDouble Nektar::MMFAdvection::m_advy
protected

Definition at line 88 of file MMFAdvection.h.

Referenced by Test2Dproblem(), Test3Dproblem(), and v_InitObject().

◆ m_advz

NekDouble Nektar::MMFAdvection::m_advz
protected

Definition at line 88 of file MMFAdvection.h.

Referenced by Test3Dproblem(), and v_InitObject().

◆ m_Mass0

NekDouble Nektar::MMFAdvection::m_Mass0
protected

Definition at line 91 of file MMFAdvection.h.

Referenced by v_DoSolve(), and v_SetInitialConditions().

◆ m_planeNumber

int Nektar::MMFAdvection::m_planeNumber
protected

Definition at line 103 of file MMFAdvection.h.

Referenced by MMFAdvection().

◆ m_riemannSolver

SolverUtils::RiemannSolverSharedPtr Nektar::MMFAdvection::m_riemannSolver
protected

Definition at line 86 of file MMFAdvection.h.

Referenced by v_InitObject().

◆ m_RotAngle

NekDouble Nektar::MMFAdvection::m_RotAngle
protected

Definition at line 89 of file MMFAdvection.h.

Referenced by EvaluateAdvectionVelocity(), v_GenerateSummary(), and v_InitObject().

◆ m_TestType

TestType Nektar::MMFAdvection::m_TestType

◆ m_traceVn

Array<OneD, NekDouble> Nektar::MMFAdvection::m_traceVn
protected

Definition at line 96 of file MMFAdvection.h.

Referenced by GetNormalVelocity(), v_InitObject(), and WeakDGDirectionalAdvection().

◆ m_veldotMF

Array<OneD, Array<OneD, NekDouble> > Nektar::MMFAdvection::m_veldotMF
protected

Definition at line 98 of file MMFAdvection.h.

Referenced by v_InitObject().

◆ m_vellc

Array<OneD, NekDouble> Nektar::MMFAdvection::m_vellc
protected

Definition at line 99 of file MMFAdvection.h.

Referenced by v_InitObject().

◆ m_velocity

Array<OneD, Array<OneD, NekDouble> > Nektar::MMFAdvection::m_velocity
protected

◆ m_VelProjection

int Nektar::MMFAdvection::m_VelProjection
protected

Definition at line 92 of file MMFAdvection.h.

Referenced by EvaluateAdvectionVelocity(), and v_InitObject().

◆ m_waveFreq

NekDouble Nektar::MMFAdvection::m_waveFreq
protected

Definition at line 89 of file MMFAdvection.h.

Referenced by EvaluateAdvectionVelocity(), and v_InitObject().