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

Static Public Member Functions

static SolverUtils::EquationSystemSharedPtr create (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Creates an instance of this class. 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)
 
 ~MMFAdvection () override=default
 
void v_InitObject (bool DeclareFields=true) override
 Initialise the object. More...
 
void WeakDGDirectionalAdvection (const Array< OneD, Array< OneD, NekDouble > > &InField, Array< OneD, Array< OneD, NekDouble > > &OutField)
 
Array< OneD, NekDouble > & GetNormalVelocity ()
 Get the normal velocity. 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 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 EvaluateAdvectionVelocity (Array< OneD, Array< OneD, NekDouble > > &velocity)
 
NekDouble ComputeCirculatingArclength (const NekDouble zlevel, const NekDouble Rhs)
 
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)
 
void v_DoSolve () override
 Solves an unsteady problem. More...
 
void v_GenerateSummary (SolverUtils::SummaryList &s) override
 Print Summary. More...
 
void v_SetInitialConditions (const NekDouble initialtime, bool dumpInitialConditions, const int domain) override
 
void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time) override
 
- 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 void v_InitObject (bool DeclareField=true) override
 Init object for UnsteadySystem class. More...
 
SOLVER_UTILS_EXPORT void v_DoSolve () override
 Solves an unsteady problem. More...
 
virtual SOLVER_UTILS_EXPORT void v_PrintStatusInformation (const int step, const NekDouble cpuTime)
 Print Status Information. More...
 
virtual SOLVER_UTILS_EXPORT void v_PrintSummaryStatistics (const NekDouble intTime)
 Print Summary Statistics. More...
 
SOLVER_UTILS_EXPORT void v_DoInitialise (bool dumpInitialConditions=true) override
 Sets up initial conditions. More...
 
SOLVER_UTILS_EXPORT void v_GenerateSummary (SummaryList &s) override
 Print a summary of time stepping parameters. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_GetTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray)
 Return the timestep to be used for the next step in the time-marching loop. More...
 
virtual SOLVER_UTILS_EXPORT bool v_PreIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_PostIntegrate (int step)
 
virtual SOLVER_UTILS_EXPORT bool v_RequireFwdTrans ()
 
virtual SOLVER_UTILS_EXPORT void v_SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
virtual SOLVER_UTILS_EXPORT bool v_UpdateTimeStepCheck ()
 
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator ()
 Get the maximum timestep estimator for cfl control. More...
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time, int &nchk)
 
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff (const Array< OneD, Array< OneD, NekDouble > > vel, StdRegions::VarCoeffMap &varCoeffMap)
 Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity. More...
 
SOLVER_UTILS_EXPORT void DoDummyProjection (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble time)
 Perform dummy projection. More...
 
- Protected Member Functions inherited from Nektar::SolverUtils::EquationSystem
SOLVER_UTILS_EXPORT EquationSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises EquationSystem class members. More...
 
virtual SOLVER_UTILS_EXPORT void v_InitObject (bool DeclareFeld=true)
 Initialisation object for EquationSystem. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoInitialise (bool dumpInitialConditions=true)
 Virtual function for initialisation implementation. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoSolve ()
 Virtual function for solve implementation. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Virtual function for the L_inf error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the L_2 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_H1Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the H_1 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary (SummaryList &l)
 Virtual function for generating summary information. More...
 
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure (void)
 
virtual SOLVER_UTILS_EXPORT bool v_NegatedOp (void)
 Virtual function to identify if operator is negated in DoSolve. More...
 
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 
- Protected Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
SOLVER_UTILS_EXPORT void v_InitObject (bool DeclareField=true) override
 Initialisation object for EquationSystem. More...
 
SOLVER_UTILS_EXPORT bool v_PostIntegrate (int step) override
 
virtual SOLVER_UTILS_EXPORT Array< OneD, NekDoublev_GetMaxStdVelocity (const NekDouble SpeedSoundFactor=1.0)
 

Protected Attributes

SolverUtils::RiemannSolverSharedPtr m_riemannSolver
 
TestType m_TestType
 
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
 
Array< OneD, NekDoublem_traceVn
 
Array< OneD, Array< OneD, NekDouble > > m_veldotMF
 
Array< OneD, NekDoublem_vellc
 
- 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
 
- Protected Attributes inherited from Nektar::SolverUtils::UnsteadySystem
LibUtilities::TimeIntegrationSchemeSharedPtr m_intScheme
 Wrapper to the time integration scheme. More...
 
LibUtilities::TimeIntegrationSchemeOperators m_ode
 The time integration scheme operators to use. More...
 
Array< OneD, Array< OneD, NekDouble > > m_previousSolution
 Storage for previous solution for steady-state check. More...
 
std::vector< int > m_intVariables
 
NekDouble m_cflSafetyFactor
 CFL safety factor (comprise between 0 to 1). More...
 
NekDouble m_CFLGrowth
 CFL growth rate. More...
 
NekDouble m_CFLEnd
 Maximun cfl in cfl growth. More...
 
int m_abortSteps
 Number of steps between checks for abort conditions. More...
 
bool m_explicitDiffusion
 Indicates if explicit or implicit treatment of diffusion is used. More...
 
bool m_explicitAdvection
 Indicates if explicit or implicit treatment of advection is used. More...
 
bool m_explicitReaction
 Indicates if explicit or implicit treatment of reaction is used. More...
 
int m_steadyStateSteps
 Check for steady state at step interval. More...
 
NekDouble m_steadyStateTol
 Tolerance to which steady state should be evaluated at. More...
 
int m_filtersInfosteps
 Number of time steps between outputting filters information. More...
 
std::vector< std::pair< std::string, FilterSharedPtr > > m_filters
 
bool m_homoInitialFwd
 Flag to determine if simulation should start in homogeneous forward transformed state. More...
 
std::ofstream m_errFile
 
NekDouble m_epsilon
 Diffusion coefficient. More...
 
- Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
LibUtilities::CommSharedPtr m_comm
 Communicator. More...
 
bool m_verbose
 
LibUtilities::SessionReaderSharedPtr m_session
 The session reader. More...
 
std::map< std::string, SolverUtils::SessionFunctionSharedPtrm_sessionFunctions
 Map of known SessionFunctions. More...
 
LibUtilities::FieldIOSharedPtr m_fld
 Field input/output. More...
 
Array< OneD, MultiRegions::ExpListSharedPtrm_fields
 Array holding all dependent variables. More...
 
SpatialDomains::BoundaryConditionsSharedPtr m_boundaryConditions
 Pointer to boundary conditions object. More...
 
SpatialDomains::MeshGraphSharedPtr m_graph
 Pointer to graph defining mesh. More...
 
std::string m_sessionName
 Name of the session. More...
 
NekDouble m_time
 Current time of simulation. More...
 
int m_initialStep
 Number of the step where the simulation should begin. More...
 
NekDouble m_fintime
 Finish time of the simulation. More...
 
NekDouble m_timestep
 Time step size. More...
 
NekDouble m_lambda
 Lambda constant in real system if one required. More...
 
NekDouble m_checktime
 Time between checkpoints. More...
 
NekDouble m_lastCheckTime
 
NekDouble m_TimeIncrementFactor
 
int m_nchk
 Number of checkpoints written so far. More...
 
int m_steps
 Number of steps to take. More...
 
int m_checksteps
 Number of steps between checkpoints. More...
 
int m_infosteps
 Number of time steps between outputting status information. More...
 
int m_iterPIT = 0
 Number of parallel-in-time time iteration. More...
 
int m_windowPIT = 0
 Index of windows for parallel-in-time time iteration. More...
 
int m_spacedim
 Spatial dimension (>= expansion dim). More...
 
int m_expdim
 Expansion dimension. More...
 
bool m_singleMode
 Flag to determine if single homogeneous mode is used. More...
 
bool m_halfMode
 Flag to determine if half homogeneous mode is used. More...
 
bool m_multipleModes
 Flag to determine if use multiple homogenenous modes are used. More...
 
bool m_useFFT
 Flag to determine if FFT is used for homogeneous transform. More...
 
bool m_homogen_dealiasing
 Flag to determine if dealiasing is used for homogeneous simulations. More...
 
bool m_specHP_dealiasing
 Flag to determine if dealisising is usde for the Spectral/hp element discretisation. More...
 
enum MultiRegions::ProjectionType m_projectionType
 Type of projection; e.g continuous or discontinuous. More...
 
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
 Array holding trace normals for DG simulations in the forwards direction. More...
 
Array< OneD, bool > m_checkIfSystemSingular
 Flag to indicate if the fields should be checked for singularity. More...
 
LibUtilities::FieldMetaDataMap m_fieldMetaDataMap
 Map to identify relevant solver info to dump in output fields. More...
 
Array< OneD, NekDoublem_movingFrameData
 Moving reference frame status in the inertial frame X, Y, Z, Theta_x, Theta_y, Theta_z, U, V, W, Omega_x, Omega_y, Omega_z, A_x, A_y, A_z, DOmega_x, DOmega_y, DOmega_z, pivot_x, pivot_y, pivot_z. More...
 
std::vector< std::string > m_strFrameData
 variable name in m_movingFrameData More...
 
int m_NumQuadPointsError
 Number of Quadrature points used to work out the error. More...
 
enum HomogeneousType m_HomogeneousType
 
NekDouble m_LhomX
 physical length in X direction (if homogeneous) More...
 
NekDouble m_LhomY
 physical length in Y direction (if homogeneous) More...
 
NekDouble m_LhomZ
 physical length in Z direction (if homogeneous) More...
 
int m_npointsX
 number of points in X direction (if homogeneous) More...
 
int m_npointsY
 number of points in Y direction (if homogeneous) More...
 
int m_npointsZ
 number of points in Z direction (if homogeneous) More...
 
int m_HomoDirec
 number of homogenous directions More...
 
- Protected Attributes inherited from Nektar::SolverUtils::ALEHelper
Array< OneD, MultiRegions::ExpListSharedPtrm_fieldsALE
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocity
 
Array< OneD, Array< OneD, NekDouble > > m_gridVelocityTrace
 
std::vector< ALEBaseShPtrm_ALEs
 
bool m_ALESolver = false
 
bool m_ImplicitALESolver = false
 
NekDouble m_prevStageTime = 0.0
 
int m_spaceDim
 
- Protected Attributes inherited from Nektar::SolverUtils::AdvectionSystem
SolverUtils::AdvectionSharedPtr m_advObject
 Advection term. More...
 

Friends

class MemoryManager< MMFAdvection >
 

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::MMFSystem
SOLVER_UTILS_EXPORT MMFSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
SOLVER_UTILS_EXPORT ~MMFSystem () override
 
SOLVER_UTILS_EXPORT void v_GenerateSummary (SummaryList &s) override
 Virtual function for generating summary information. More...
 
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
SOLVER_UTILS_EXPORT ~UnsteadySystem () override=default
 Destructor. More...
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep (const Array< OneD, const Array< OneD, NekDouble > > &inarray)
 Calculate the larger time-step mantaining the problem stable. More...
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep ()
 
SOLVER_UTILS_EXPORT void SetTimeStep (const NekDouble timestep)
 
SOLVER_UTILS_EXPORT void SteadyStateResidual (int step, Array< OneD, NekDouble > &L2)
 
SOLVER_UTILS_EXPORT LibUtilities::TimeIntegrationSchemeSharedPtrGetTimeIntegrationScheme ()
 Returns the time integration scheme. More...
 
SOLVER_UTILS_EXPORT LibUtilities::TimeIntegrationSchemeOperatorsGetTimeIntegrationSchemeOperators ()
 Returns the time integration scheme operators. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::EquationSystem
virtual SOLVER_UTILS_EXPORT ~EquationSystem ()
 Destructor. More...
 
SOLVER_UTILS_EXPORT void InitObject (bool DeclareField=true)
 Initialises the members of this object. More...
 
SOLVER_UTILS_EXPORT void DoInitialise (bool dumpInitialConditions=true)
 Perform any initialisation necessary before solving the problem. More...
 
SOLVER_UTILS_EXPORT void DoSolve ()
 Solve the problem. More...
 
SOLVER_UTILS_EXPORT void TransCoeffToPhys ()
 Transform from coefficient to physical space. More...
 
SOLVER_UTILS_EXPORT void TransPhysToCoeff ()
 Transform from physical to coefficient space. More...
 
SOLVER_UTILS_EXPORT void Output ()
 Perform output operations after solve. More...
 
SOLVER_UTILS_EXPORT std::string GetSessionName ()
 Get Session name. More...
 
template<class T >
std::shared_ptr< T > as ()
 
SOLVER_UTILS_EXPORT void ResetSessionName (std::string newname)
 Reset Session name. More...
 
SOLVER_UTILS_EXPORT LibUtilities::SessionReaderSharedPtr GetSession ()
 Get Session name. More...
 
SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr GetPressure ()
 Get pressure field if available. More...
 
SOLVER_UTILS_EXPORT void ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 
SOLVER_UTILS_EXPORT void PrintSummary (std::ostream &out)
 Print a summary of parameters and solver characteristics. More...
 
SOLVER_UTILS_EXPORT void SetLambda (NekDouble lambda)
 Set parameter m_lambda. More...
 
SOLVER_UTILS_EXPORT SessionFunctionSharedPtr GetFunction (std::string name, const MultiRegions::ExpListSharedPtr &field=MultiRegions::NullExpListSharedPtr, bool cache=false)
 Get a SessionFunction by name. More...
 
SOLVER_UTILS_EXPORT void SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 Initialise the data in the dependent fields. More...
 
SOLVER_UTILS_EXPORT void EvaluateExactSolution (int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 Evaluates an exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln, bool Normalised=false)
 Compute the L2 error between fields and a given exact solution. More...
 
SOLVER_UTILS_EXPORT NekDouble L2Error (unsigned int field, bool Normalised=false)
 Compute the L2 error of the fields. More...
 
SOLVER_UTILS_EXPORT NekDouble LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Linf error computation. More...
 
SOLVER_UTILS_EXPORT NekDouble H1Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln, bool Normalised=false)
 Compute the H1 error between fields and a given exact solution. More...
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleErrorExtraPoints (unsigned int field)
 Compute error (L2 and L_inf) over an larger set of quadrature points return [L2 Linf]. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n)
 Write checkpoint file of m_fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_Output (const int n, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 Write checkpoint file of custom data fields. More...
 
SOLVER_UTILS_EXPORT void Checkpoint_BaseFlow (const int n)
 Write base flow file of m_fields. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname)
 Write field data to the given filename. More...
 
SOLVER_UTILS_EXPORT void WriteFld (const std::string &outname, MultiRegions::ExpListSharedPtr &field, std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 Write input fields to the given filename. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields)
 Input field data from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFldToMultiDomains (const std::string &infile, Array< OneD, MultiRegions::ExpListSharedPtr > &pFields, const int ndomains)
 Input field data from the given file to multiple domains. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, std::vector< std::string > &fieldStr, Array< OneD, Array< OneD, NekDouble > > &coeffs)
 Output a field. Input field data into array from the given file. More...
 
SOLVER_UTILS_EXPORT void ImportFld (const std::string &infile, MultiRegions::ExpListSharedPtr &pField, std::string &pFieldName)
 Output a field. Input field data into ExpList from the given file. More...
 
SOLVER_UTILS_EXPORT void SessionSummary (SummaryList &vSummary)
 Write out a session summary. More...
 
SOLVER_UTILS_EXPORT Array< OneD, MultiRegions::ExpListSharedPtr > & UpdateFields ()
 
SOLVER_UTILS_EXPORT LibUtilities::FieldMetaDataMapUpdateFieldMetaDataMap ()
 Get hold of FieldInfoMap so it can be updated. More...
 
SOLVER_UTILS_EXPORT NekDouble GetTime ()
 Return final time. More...
 
SOLVER_UTILS_EXPORT int GetNcoeffs ()
 
SOLVER_UTILS_EXPORT int GetNcoeffs (const int eid)
 
SOLVER_UTILS_EXPORT int GetNumExpModes ()
 
SOLVER_UTILS_EXPORT const Array< OneD, int > GetNumExpModesPerExp ()
 
SOLVER_UTILS_EXPORT int GetNvariables ()
 
SOLVER_UTILS_EXPORT const std::string GetVariable (unsigned int i)
 
SOLVER_UTILS_EXPORT int GetTraceTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTraceNpoints ()
 
SOLVER_UTILS_EXPORT int GetExpSize ()
 
SOLVER_UTILS_EXPORT int GetPhys_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetCoeff_Offset (int n)
 
SOLVER_UTILS_EXPORT int GetTotPoints ()
 
SOLVER_UTILS_EXPORT int GetTotPoints (int n)
 
SOLVER_UTILS_EXPORT int GetNpoints ()
 
SOLVER_UTILS_EXPORT int GetSteps ()
 
SOLVER_UTILS_EXPORT void SetSteps (const int steps)
 
SOLVER_UTILS_EXPORT NekDouble GetTimeStep ()
 
SOLVER_UTILS_EXPORT void CopyFromPhysField (const int i, Array< OneD, NekDouble > &output)
 
SOLVER_UTILS_EXPORT void CopyToPhysField (const int i, const Array< OneD, const NekDouble > &input)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDouble > & UpdatePhysField (const int i)
 
SOLVER_UTILS_EXPORT void ZeroPhysFields ()
 
SOLVER_UTILS_EXPORT void FwdTransFields ()
 
SOLVER_UTILS_EXPORT void SetModifiedBasis (const bool modbasis)
 
SOLVER_UTILS_EXPORT int GetCheckpointNumber ()
 
SOLVER_UTILS_EXPORT void SetCheckpointNumber (int num)
 
SOLVER_UTILS_EXPORT int GetCheckpointSteps ()
 
SOLVER_UTILS_EXPORT void SetCheckpointSteps (int num)
 
SOLVER_UTILS_EXPORT int GetInfoSteps ()
 
SOLVER_UTILS_EXPORT void SetInfoSteps (int num)
 
SOLVER_UTILS_EXPORT void SetIterationNumberPIT (int num)
 
SOLVER_UTILS_EXPORT void SetWindowNumberPIT (int num)
 
SOLVER_UTILS_EXPORT Array< OneD, const Array< OneD, NekDouble > > GetTraceNormals ()
 
SOLVER_UTILS_EXPORT void SetTime (const NekDouble time)
 
SOLVER_UTILS_EXPORT void SetTimeStep (const NekDouble timestep)
 
SOLVER_UTILS_EXPORT void SetInitialStep (const int step)
 
SOLVER_UTILS_EXPORT void SetBoundaryConditions (NekDouble time)
 Evaluates the boundary conditions at the given time. More...
 
SOLVER_UTILS_EXPORT bool NegatedOp ()
 Identify if operator is negated in DoSolve. More...
 
- Public Member Functions inherited from Nektar::SolverUtils::ALEHelper
virtual ~ALEHelper ()=default
 
virtual SOLVER_UTILS_EXPORT void v_ALEInitObject (int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 
SOLVER_UTILS_EXPORT void InitObject (int spaceDim, Array< OneD, MultiRegions::ExpListSharedPtr > &fields)
 
virtual SOLVER_UTILS_EXPORT void v_UpdateGridVelocity (const NekDouble &time)
 
virtual SOLVER_UTILS_EXPORT void v_ALEPreMultiplyMass (Array< OneD, Array< OneD, NekDouble > > &fields)
 
SOLVER_UTILS_EXPORT void ALEDoElmtInvMass (Array< OneD, Array< OneD, NekDouble > > &traceNormals, Array< OneD, Array< OneD, NekDouble > > &fields, NekDouble time)
 Update m_fields with u^n by multiplying by inverse mass matrix. That's then used in e.g. checkpoint output and L^2 error calculation. More...
 
SOLVER_UTILS_EXPORT void ALEDoElmtInvMassBwdTrans (const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
 
SOLVER_UTILS_EXPORT void MoveMesh (const NekDouble &time, Array< OneD, Array< OneD, NekDouble > > &traceNormals)
 
const Array< OneD, const Array< OneD, NekDouble > > & GetGridVelocity ()
 
SOLVER_UTILS_EXPORT const Array< OneD, const Array< OneD, NekDouble > > & GetGridVelocityTrace ()
 
SOLVER_UTILS_EXPORT void ExtraFldOutputGridVelocity (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 
- Public Member Functions inherited from Nektar::SolverUtils::AdvectionSystem
SOLVER_UTILS_EXPORT AdvectionSystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 
SOLVER_UTILS_EXPORT ~AdvectionSystem () override=default
 
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)
 
- 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
 
- Protected Types inherited from Nektar::SolverUtils::EquationSystem
enum  HomogeneousType { eHomogeneous1D , eHomogeneous2D , eHomogeneous3D , eNotHomogeneous }
 Parameter for homogeneous expansions. More...
 
- Static Protected Attributes inherited from Nektar::SolverUtils::EquationSystem
static std::string equationSystemTypeLookupIds []
 
static std::string projectionTypeLookupIds []
 

Detailed Description

Definition at line 62 of file MMFAdvection.h.

Constructor & Destructor Documentation

◆ MMFAdvection()

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

Definition at line 51 of file MMFAdvection.cpp.

53 : UnsteadySystem(pSession, pGraph), MMFSystem(pSession, pGraph),
54 AdvectionSystem(pSession, pGraph)
55{
56}
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:41
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.

◆ ~MMFAdvection()

Nektar::MMFAdvection::~MMFAdvection ( )
overrideprotecteddefault

Member Function Documentation

◆ AdvectionBellPlane()

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

Definition at line 902 of file MMFAdvection.cpp.

903{
904 int nq = m_fields[0]->GetNpoints();
905
906 NekDouble dist, m_radius_limit;
907
908 Array<OneD, NekDouble> x(nq);
909 Array<OneD, NekDouble> y(nq);
910 Array<OneD, NekDouble> z(nq);
911
912 m_fields[0]->GetCoords(x, y, z);
913
914 // Sets of parameters
915 m_radius_limit = 0.5;
916
917 NekDouble x0j, x1j;
918 outfield = Array<OneD, NekDouble>(nq, 0.0);
919 for (int j = 0; j < nq; ++j)
920 {
921 x0j = x[j];
922 x1j = y[j];
923
924 dist = sqrt(x0j * x0j + x1j * x1j);
925
926 if (dist < m_radius_limit)
927 {
928 outfield[j] = 0.5 * (1.0 + cos(m_pi * dist / m_radius_limit));
929 }
930 else
931 {
932 outfield[j] = 0.0;
933 }
934 }
935}
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
std::vector< double > z(NPUPPER)
double NekDouble
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:285

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

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ AdvectionBellSphere()

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

Definition at line 937 of file MMFAdvection.cpp.

938{
939 int nq = m_fields[0]->GetNpoints();
940
941 NekDouble dist, radius, cosdiff, sin_theta, cos_theta, sin_varphi,
942 cos_varphi;
943 NekDouble m_theta_c, m_varphi_c, m_radius_limit, m_c0;
944
945 Array<OneD, NekDouble> x(nq);
946 Array<OneD, NekDouble> y(nq);
947 Array<OneD, NekDouble> z(nq);
948
949 m_fields[0]->GetCoords(x, y, z);
950
951 // Sets of parameters
952 m_theta_c = 0.0;
953 m_varphi_c = 3.0 * m_pi / 2.0;
954 m_radius_limit = 7.0 * m_pi / 64.0;
955 m_c0 = 0.0;
956
957 NekDouble x0j, x1j, x2j;
958 outfield = Array<OneD, NekDouble>(nq, 0.0);
959 for (int j = 0; j < nq; ++j)
960 {
961 x0j = x[j];
962 x1j = y[j];
963 x2j = z[j];
964
965 radius = sqrt(x0j * x0j + x1j * x1j + x2j * x2j);
966
967 sin_varphi = x1j / sqrt(x0j * x0j + x1j * x1j);
968 cos_varphi = x0j / sqrt(x0j * x0j + x1j * x1j);
969
970 sin_theta = x2j / radius;
971 cos_theta = sqrt(x0j * x0j + x1j * x1j) / radius;
972
973 cosdiff = cos_varphi * cos(m_varphi_c) + sin_varphi * sin(m_varphi_c);
974 dist = radius * acos(sin(m_theta_c) * sin_theta +
975 cos(m_theta_c) * cos_theta * cosdiff);
976
977 if (dist < m_radius_limit)
978 {
979 outfield[j] =
980 0.5 * (1.0 + cos(m_pi * dist / m_radius_limit)) + m_c0;
981 }
982 else
983 {
984 outfield[j] = m_c0;
985 }
986 }
987}

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

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ ComputeCirculatingArclength()

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

Definition at line 679 of file MMFAdvection.cpp.

681{
682
683 NekDouble Tol = 0.0001, Maxiter = 1000, N = 100;
684 NekDouble newy, F = 0.0, dF = 1.0, y0, tmp;
685
686 Array<OneD, NekDouble> xp(N + 1);
687 Array<OneD, NekDouble> yp(N + 1);
688
689 NekDouble intval = 0.0;
690 switch (m_surfaceType)
691 {
694 {
695 intval = sqrt(Rhs - zlevel * zlevel);
696 }
697 break;
698
700 {
701 intval = sqrt(0.5 * (Rhs - zlevel * zlevel * zlevel * zlevel -
702 zlevel * zlevel));
703 }
704 break;
705
707 {
708 tmp = 0.5 *
709 (Rhs - zlevel * zlevel * zlevel * zlevel - zlevel * zlevel);
710 intval = sqrt(0.5 * (1.0 + sqrt(1.0 + 4.0 * tmp)));
711 }
712 break;
713
714 default:
715 break;
716 }
717
718 switch (m_surfaceType)
719 {
720 // Find the half of all the xp and yp on zlevel ....
724 {
725 for (int j = 0; j < N + 1; ++j)
726 {
727 xp[j] = j * 2.0 * intval / N - intval;
728
729 y0 = 1.0;
730 for (int i = 0; i < Maxiter; ++i)
731 {
732 switch (m_surfaceType)
733 {
734 // Find the half of all the xp and yp on zlevel ....
737 {
738 F = xp[j] * xp[j] + y0 * y0 + zlevel * zlevel - Rhs;
739 dF = 2.0 * y0;
740 }
741 break;
742
744 {
745 F = 2.0 * xp[j] * xp[j] + y0 * y0 * y0 * y0 +
746 y0 * y0 + zlevel * zlevel * zlevel * zlevel +
747 zlevel * zlevel - Rhs;
748 dF = 4.0 * y0 * y0 * y0 + 2.0 * y0;
749 }
750 break;
751
752 default:
753 break;
754 }
755
756 newy = y0 - F / dF;
757
758 if (fabs(F / dF) < Tol)
759 {
760 yp[j] = newy;
761 break;
762 }
763
764 else
765 {
766 y0 = newy;
767 }
768
769 ASSERTL0(i < Maxiter,
770 "Advection Velocity convergence fails");
771
772 } // i-loop
773 }
774 }
775 break;
776
778 {
779 for (int j = 0; j < N + 1; ++j)
780 {
781 xp[j] = j * 2.0 * intval / N - intval;
782 tmp = 0.5 * Rhs -
783 0.5 * (zlevel * zlevel * zlevel * zlevel +
784 zlevel * zlevel) -
785 (xp[j] * xp[j] * xp[j] * xp[j] - xp[j] * xp[j]);
786 if (tmp < 0)
787 {
788 tmp = -1.0 * tmp;
789 }
790 yp[j] = sqrt(tmp);
791 } // j-loop
792 }
793 break;
794
795 default:
796 break;
797
798 } // switch-loop
799
800 NekDouble pi = 3.14159265358979323846;
801 NekDouble arclength = 0.0;
802 for (int j = 0; j < N; ++j)
803 {
804 arclength =
805 arclength + sqrt((yp[j + 1] - yp[j]) * (yp[j + 1] - yp[j]) +
806 (xp[j + 1] - xp[j]) * (xp[j + 1] - xp[j])) /
807 pi;
808 }
809
810 return arclength;
811}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208

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 1034 of file MMFAdvection.cpp.

1035{
1036 int nq = m_fields[0]->GetNpoints();
1037
1038 Array<OneD, NekDouble> velcoeff(nq, 0.0);
1039
1040 Array<OneD, NekDouble> Dtmp0(nq);
1041 Array<OneD, NekDouble> Dtmp1(nq);
1042 Array<OneD, NekDouble> Dtmp2(nq);
1043 Array<OneD, NekDouble> Drv(nq);
1044
1045 vellc = Array<OneD, NekDouble>(nq, 0.0);
1046
1047 // m_vellc = \nabla m_vel \cdot tan_i
1048 Array<OneD, NekDouble> tmp(nq);
1049 Array<OneD, NekDouble> vessel(nq);
1050
1051 for (int j = 0; j < m_shapedim; ++j)
1052 {
1053 Vmath::Zero(nq, velcoeff, 1);
1054 for (int k = 0; k < m_spacedim; ++k)
1055 {
1056 // a_j = tan_j cdot m_vel
1057 Vmath::Vvtvp(nq, &m_movingframes[j][k * nq], 1, &m_velocity[k][0],
1058 1, &velcoeff[0], 1, &velcoeff[0], 1);
1059 }
1060
1061 // d a_j / d x^k
1062 m_fields[0]->PhysDeriv(velcoeff, Dtmp0, Dtmp1, Dtmp2);
1063
1064 for (int k = 0; k < m_spacedim; ++k)
1065 {
1066 // tan_j^k ( d a_j / d x^k )
1067 switch (k)
1068 {
1069 case (0):
1070 {
1071 Vmath::Vvtvp(nq, &Dtmp0[0], 1, &m_movingframes[j][k * nq],
1072 1, &vellc[0], 1, &vellc[0], 1);
1073 }
1074 break;
1075
1076 case (1):
1077 {
1078 Vmath::Vvtvp(nq, &Dtmp1[0], 1, &m_movingframes[j][k * nq],
1079 1, &vellc[0], 1, &vellc[0], 1);
1080 }
1081 break;
1082
1083 case (2):
1084 {
1085 Vmath::Vvtvp(nq, &Dtmp2[0], 1, &m_movingframes[j][k * nq],
1086 1, &vellc[0], 1, &vellc[0], 1);
1087 }
1088 break;
1089 }
1090 }
1091 }
1092}
Array< OneD, Array< OneD, NekDouble > > m_velocity
Definition: MMFAdvection.h:89
int m_spacedim
Spatial dimension (>= expansion dim).
Array< OneD, Array< OneD, NekDouble > > m_movingframes
Definition: MMFSystem.h:183
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.hpp:366
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273

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 1094 of file MMFAdvection.cpp.

1096{
1097 int nq = m_fields[0]->GetNpoints();
1098
1099 veldotMF = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1100
1101 Array<OneD, NekDouble> magMF(nq);
1102 for (int j = 0; j < m_shapedim; ++j)
1103 {
1104 veldotMF[j] = Array<OneD, NekDouble>(nq, 0.0);
1105
1106 for (int k = 0; k < m_spacedim; ++k)
1107 {
1108 Vmath::Vvtvp(nq, &m_movingframes[j][k * nq], 1, &m_velocity[k][0],
1109 1, &veldotMF[j][0], 1, &veldotMF[j][0], 1);
1110 }
1111 }
1112}

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 69 of file MMFAdvection.h.

72 {
75 p->InitObject();
76 return p;
77 }
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 448 of file MMFAdvection.cpp.

451{
452 // Counter variable
453 int i;
454
455 // Number of fields (variables of the problem)
456 int nVariables = inarray.size();
457
458 // Set the boundary conditions
460
461 // Switch on the projection type (Discontinuous or Continuous)
462 switch (m_projectionType)
463 {
464 // Discontinuous projection
466 {
467 // Number of quadrature points
468 int nQuadraturePts = GetNpoints();
469
470 // Just copy over array
471 if (inarray != outarray)
472 {
473 for (i = 0; i < nVariables; ++i)
474 {
475 Vmath::Vcopy(nQuadraturePts, inarray[i], 1, outarray[i], 1);
476 }
477 }
478 break;
479 }
480
481 // Continuous projection
484 {
485 Array<OneD, NekDouble> coeffs(m_fields[0]->GetNcoeffs(), 0.0);
486 for (i = 0; i < nVariables; ++i)
487 {
488 m_fields[i]->FwdTrans(inarray[i], coeffs);
489 m_fields[i]->BwdTrans(coeffs, outarray[i]);
490 }
491 break;
492 }
493
494 default:
495 ASSERTL0(false, "Unknown projection scheme");
496 break;
497 }
498}
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.hpp:825

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 385 of file MMFAdvection.cpp.

389{
390 int i;
391 int nvariables = inarray.size();
392 int npoints = GetNpoints();
393
394 switch (m_projectionType)
395 {
397 {
398 int ncoeffs = inarray[0].size();
399
400 if (m_spacedim == 3)
401 {
402 Array<OneD, Array<OneD, NekDouble>> WeakAdv(nvariables);
403
404 WeakAdv[0] = Array<OneD, NekDouble>(ncoeffs * nvariables);
405 for (i = 1; i < nvariables; ++i)
406 {
407 WeakAdv[i] = WeakAdv[i - 1] + ncoeffs;
408 }
409
410 // Compute \nabla \cdot \vel u according to MMF scheme
411 WeakDGDirectionalAdvection(inarray, WeakAdv);
412
413 for (i = 0; i < nvariables; ++i)
414 {
415 m_fields[i]->MultiplyByElmtInvMass(WeakAdv[i], WeakAdv[i]);
416 m_fields[i]->BwdTrans(WeakAdv[i], outarray[i]);
417 Vmath::Neg(npoints, outarray[i], 1);
418 }
419 }
420 else
421 {
422 m_advObject->Advect(2, m_fields, m_velocity, inarray, outarray,
423 0.0);
424
425 for (i = 0; i < nvariables; ++i)
426 {
427 Vmath::Neg(npoints, outarray[i], 1);
428 }
429 }
430 }
431 break;
432
433 default:
434 {
435 ASSERTL0(false, "Unknown projection scheme");
436 }
437 break;
438 }
439}
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.hpp:292

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 597 of file MMFAdvection.cpp.

599{
600 int nq = m_fields[0]->GetNpoints();
601
602 NekDouble vel_phi, vel_theta, sin_varphi, cos_varphi, sin_theta, cos_theta;
603 NekDouble x0j, x1j, x2j;
604
605 Array<OneD, NekDouble> x0(nq);
606 Array<OneD, NekDouble> x1(nq);
607 Array<OneD, NekDouble> x2(nq);
608
609 m_fields[0]->GetCoords(x0, x1, x2);
610
611 // theta = a*sin(z/r), phi = a*tan(y/x);
612 for (int j = 0; j < nq; j++)
613 {
614 x0j = x0[j];
615 x1j = x1[j];
616 x2j = x2[j];
617
618 CartesianToSpherical(x0j, x1j, x2j, sin_varphi, cos_varphi, sin_theta,
619 cos_theta);
620
621 vel_phi = m_waveFreq * (cos_theta * cos(m_RotAngle) +
622 sin_theta * cos_varphi * sin(m_RotAngle));
623 vel_theta = -1.0 * m_waveFreq * sin_theta * sin(m_RotAngle);
624
625 velocity[0][j] =
626 -1.0 * vel_phi * sin_varphi - vel_theta * sin_theta * cos_varphi;
627 velocity[1][j] =
628 vel_phi * cos_varphi - vel_theta * sin_theta * sin_varphi;
629 velocity[2][j] = vel_theta * cos_theta;
630 }
631
632 // Project the veloicty on the tangent plane
633
634 if (m_VelProjection)
635 {
636 // Check MovingFrames \cdot SurfaceNormals = 0
637 Array<OneD, Array<OneD, NekDouble>> newvelocity(m_spacedim);
638
639 Array<OneD, Array<OneD, NekDouble>> MF1(m_spacedim);
640 Array<OneD, Array<OneD, NekDouble>> MF2(m_spacedim);
641 Array<OneD, Array<OneD, NekDouble>> SN(m_spacedim);
642
643 for (int k = 0; k < m_spacedim; ++k)
644 {
645 newvelocity[k] = Array<OneD, NekDouble>(nq);
646 MF1[k] = Array<OneD, NekDouble>(nq);
647 MF2[k] = Array<OneD, NekDouble>(nq);
648 SN[k] = Array<OneD, NekDouble>(nq);
649
650 Vmath::Vcopy(nq, &m_movingframes[0][k * nq], 1, &MF1[k][0], 1);
651 Vmath::Vcopy(nq, &m_movingframes[1][k * nq], 1, &MF2[k][0], 1);
652 }
653
654 VectorCrossProd(MF1, MF2, SN);
655 GramSchumitz(SN, m_velocity, newvelocity, true);
656
657 Array<OneD, NekDouble> tmp(nq, 0.0);
658 Array<OneD, NekDouble> Projection(nq, 0.0);
659
660 for (int k = 0; k < m_spacedim; ++k)
661 {
662 Vmath::Vsub(nq, &m_velocity[0][0], 1, &newvelocity[0][0], 1,
663 &tmp[0], 1);
664 Vmath::Vmul(nq, &tmp[0], 1, &tmp[0], 1, &tmp[0], 1);
665 Vmath::Vadd(nq, tmp, 1, Projection, 1, Projection, 1);
666 }
667
668 std::cout
669 << "Velocity vector is projected onto the tangent plane: Diff = "
670 << RootMeanSquare(Projection) << std::endl;
671
672 for (int k = 0; k < m_spacedim; ++k)
673 {
674 Vmath::Vcopy(nq, &newvelocity[k][0], 1, &m_velocity[k][0], 1);
675 }
676 }
677}
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:696
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:2399
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:792
SOLVER_UTILS_EXPORT NekDouble RootMeanSquare(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2338
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void 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.hpp:180
void Vsub(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Subtract vector z = x-y.
Definition: Vmath.hpp:220

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 507 of file MMFAdvection.cpp.

510{
511 ASSERTL1(flux[0].size() == m_velocity.size(),
512 "Dimension of flux array and velocity array do not match");
513
514 int i, j;
515 int nq = physfield[0].size();
516
517 for (i = 0; i < flux.size(); ++i)
518 {
519 for (j = 0; j < flux[0].size(); ++j)
520 {
521 Vmath::Vmul(nq, physfield[i], 1, m_velocity[j], 1, flux[i][j], 1);
522 }
523 }
524}
#define ASSERTL1(condition, msg)
Assert Level 1 – Debugging which is used whether in FULLDEBUG or DEBUG compilation mode....
Definition: ErrorUtil.hpp:242

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 355 of file MMFAdvection.cpp.

356{
357 // Number of trace (interface) points
358 int i;
359 int nTracePts = GetTraceNpoints();
360
361 // Auxiliary variable to compute the normal velocity
362 Array<OneD, NekDouble> tmp(nTracePts);
363
364 // Reset the normal velocity
365 Vmath::Zero(nTracePts, m_traceVn, 1);
366
367 for (i = 0; i < m_velocity.size(); ++i)
368 {
369 m_fields[0]->ExtractTracePhys(m_velocity[i], tmp);
370
371 Vmath::Vvtvp(nTracePts, m_traceNormals[i], 1, tmp, 1, m_traceVn, 1,
372 m_traceVn, 1);
373 }
374
375 return m_traceVn;
376}
Array< OneD, NekDouble > m_traceVn
Definition: MMFAdvection.h:90
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 989 of file MMFAdvection.cpp.

991{
992 int nq = m_fields[0]->GetNpoints();
993
994 Array<OneD, NekDouble> x0(nq);
995 Array<OneD, NekDouble> x1(nq);
996 Array<OneD, NekDouble> x2(nq);
997
998 m_fields[0]->GetCoords(x0, x1, x2);
999
1000 Array<OneD, NekDouble> u(nq);
1001
1002 for (int i = 0; i < nq; ++i)
1003 {
1004 u[i] = cos(m_pi * (x0[i] - m_advx * time)) *
1005 cos(m_pi * (x1[i] - m_advy * time));
1006 }
1007
1008 outfield = u;
1009}

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 1011 of file MMFAdvection.cpp.

1013{
1014 int nq = m_fields[0]->GetNpoints();
1015
1016 Array<OneD, NekDouble> x0(nq);
1017 Array<OneD, NekDouble> x1(nq);
1018 Array<OneD, NekDouble> x2(nq);
1019
1020 m_fields[0]->GetCoords(x0, x1, x2);
1021
1022 Array<OneD, NekDouble> u(nq);
1023
1024 for (int i = 0; i < nq; ++i)
1025 {
1026 u[i] = cos(m_pi * (x0[i] - m_advx * time)) *
1027 cos(m_pi * (x1[i] - m_advy * time)) *
1028 cos(m_pi * (x2[i] - m_advz * time));
1029 }
1030
1031 outfield = u;
1032}

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  )
overrideprotectedvirtual

Solves an unsteady problem.

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

Reimplemented from Nektar::SolverUtils::UnsteadySystem.

Definition at line 197 of file MMFAdvection.cpp.

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

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 1114 of file MMFAdvection.cpp.

1117{
1118 switch (m_TestType)
1119 {
1120 case eAdvectionBell:
1121 {
1122 AdvectionBellSphere(outfield);
1123 }
1124 break;
1125
1126 case eTestPlane:
1127 {
1128 Test2Dproblem(time, outfield);
1129 }
1130 break;
1131
1133 {
1134 AdvectionBellPlane(outfield);
1135 }
1136 break;
1137
1138 case eTestCube:
1139 {
1140 Test3Dproblem(time, outfield);
1141 }
1142 break;
1143
1144 default:
1145 break;
1146 }
1147}
@ 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)
overrideprotectedvirtual

Print Summary.

Reimplemented from Nektar::SolverUtils::MMFSystem.

Definition at line 1149 of file MMFAdvection.cpp.

1150{
1153 SolverUtils::AddSummaryItem(s, "Rotation Angle", m_RotAngle);
1154}
SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Virtual function for generating summary information.
Definition: MMFSystem.cpp:2463
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:47
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)
overrideprotectedvirtual

Initialise the object.

Initialisation object for the unsteady linear advection equation.

Reimplemented from Nektar::SolverUtils::AdvectionSystem.

Definition at line 61 of file MMFAdvection.cpp.

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

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 813 of file MMFAdvection.cpp.

816{
817 int nq = m_fields[0]->GetNpoints();
818
819 Array<OneD, NekDouble> u(nq);
820
821 switch (m_TestType)
822 {
823 case eAdvectionBell:
824 {
826 m_fields[0]->SetPhys(u);
827
828 m_Mass0 = m_fields[0]->Integral(u);
829
830 // forward transform to fill the modal coeffs
831 for (int i = 0; i < m_fields.size(); ++i)
832 {
833 m_fields[i]->SetPhysState(true);
834 m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
835 m_fields[i]->UpdateCoeffs());
836 }
837 }
838 break;
839
840 case eTestPlane:
841 {
842 Test2Dproblem(initialtime, u);
843 m_fields[0]->SetPhys(u);
844
845 m_Mass0 = m_fields[0]->Integral(u);
846
847 // forward transform to fill the modal coeffs
848 for (int i = 0; i < m_fields.size(); ++i)
849 {
850 m_fields[i]->SetPhysState(true);
851 m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
852 m_fields[i]->UpdateCoeffs());
853 }
854 }
855 break;
856
858 {
860 m_fields[0]->SetPhys(u);
861
862 m_Mass0 = m_fields[0]->Integral(u);
863 std::cout << "m_Mass0 = " << m_Mass0 << std::endl;
864
865 // forward transform to fill the modal coeffs
866 for (int i = 0; i < m_fields.size(); ++i)
867 {
868 m_fields[i]->SetPhysState(true);
869 m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
870 m_fields[i]->UpdateCoeffs());
871 }
872 }
873 break;
874
875 case eTestCube:
876 {
877 Test3Dproblem(initialtime, u);
878 m_fields[0]->SetPhys(u);
879
880 // forward transform to fill the modal coeffs
881 for (int i = 0; i < m_fields.size(); ++i)
882 {
883 m_fields[i]->SetPhysState(true);
884 m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
885 m_fields[i]->UpdateCoeffs());
886 }
887 }
888 break;
889
890 default:
891 break;
892 }
893
894 if (dumpInitialConditions)
895 {
896 // dump initial conditions to file
897 std::string outname = m_sessionName + "_initial.chk";
898 WriteFld(outname);
899 }
900}
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 526 of file MMFAdvection.cpp.

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

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.
static SolverUtils::EquationSystemSharedPtr create(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Creates an instance of this class.
Definition: MMFAdvection.h:69
EquationSystemFactory & GetEquationSystemFactory()

Name of class.

Definition at line 80 of file MMFAdvection.h.

◆ m_advx

NekDouble Nektar::MMFAdvection::m_advx
protected

Definition at line 85 of file MMFAdvection.h.

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

◆ m_advy

NekDouble Nektar::MMFAdvection::m_advy
protected

Definition at line 85 of file MMFAdvection.h.

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

◆ m_advz

NekDouble Nektar::MMFAdvection::m_advz
protected

Definition at line 85 of file MMFAdvection.h.

Referenced by Test3Dproblem(), and v_InitObject().

◆ m_Mass0

NekDouble Nektar::MMFAdvection::m_Mass0
protected

Definition at line 87 of file MMFAdvection.h.

Referenced by v_DoSolve(), and v_SetInitialConditions().

◆ m_riemannSolver

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

Definition at line 83 of file MMFAdvection.h.

Referenced by v_InitObject().

◆ m_RotAngle

NekDouble Nektar::MMFAdvection::m_RotAngle
protected

Definition at line 86 of file MMFAdvection.h.

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

◆ m_TestType

TestType Nektar::MMFAdvection::m_TestType
protected

◆ m_traceVn

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

Definition at line 90 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 91 of file MMFAdvection.h.

Referenced by v_InitObject().

◆ m_vellc

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

Definition at line 92 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 88 of file MMFAdvection.h.

Referenced by EvaluateAdvectionVelocity(), and v_InitObject().

◆ m_waveFreq

NekDouble Nektar::MMFAdvection::m_waveFreq
protected

Definition at line 86 of file MMFAdvection.h.

Referenced by EvaluateAdvectionVelocity(), and v_InitObject().