Nektar++
Public Member Functions | Public Attributes | Protected Member Functions | Protected Attributes | List of all members
Nektar::SolverUtils::MMFSystem Class Reference

A base class for PDEs which include an advection component. More...

#include <MMFSystem.h>

Inheritance diagram for Nektar::SolverUtils::MMFSystem:
[legend]

Public Member Functions

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
 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 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 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 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 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 Attributes

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 Member Functions

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 void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary (SummaryList &l)
 Virtual function for generating summary information. More...
 
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions (NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
 
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure (void)
 
virtual SOLVER_UTILS_EXPORT bool v_NegatedOp (void)
 Virtual function to identify if operator is negated in DoSolve. More...
 
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble > > &fieldcoeffs, std::vector< std::string > &variables)
 

Protected Attributes

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_movingFrameVelsxyz
 Moving frame of reference velocities (u, v, w, omega_x, omega_y, omega_z, a_x, a_y, a_z, domega_x, domega_y, domega_z) More...
 
Array< OneD, NekDoublem_movingFrameData
 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...
 

Additional Inherited Members

- Static Public Attributes inherited from Nektar::SolverUtils::UnsteadySystem
static std::string cmdSetStartTime
 
static std::string cmdSetStartChkNum
 
- 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

A base class for PDEs which include an advection component.

Definition at line 143 of file MMFSystem.h.

Constructor & Destructor Documentation

◆ MMFSystem()

Nektar::SolverUtils::MMFSystem::MMFSystem ( const LibUtilities::SessionReaderSharedPtr pSession,
const SpatialDomains::MeshGraphSharedPtr pGraph 
)

Definition at line 41 of file MMFSystem.cpp.

43 : UnsteadySystem(pSession, pGraph)
44{
45}
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.

◆ ~MMFSystem()

Nektar::SolverUtils::MMFSystem::~MMFSystem ( )
override

Definition at line 47 of file MMFSystem.cpp.

48{
49}

Member Function Documentation

◆ AbsIntegral()

NekDouble Nektar::SolverUtils::MMFSystem::AbsIntegral ( const Array< OneD, const NekDouble > &  inarray)
protected

Definition at line 2324 of file MMFSystem.cpp.

2325{
2326 int nq = m_fields[0]->GetNpoints();
2327 Array<OneD, NekDouble> tmp(nq);
2328
2329 if (inarray.size() != nq)
2330 {
2331 ASSERTL0(false, "AbsIntegral Error: Vector size is not correct");
2332 }
2333
2334 Vmath::Vabs(nq, inarray, 1, tmp, 1);
2335 return m_fields[0]->Integral(tmp);
2336}
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:208
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
void Vabs(int n, const T *x, const int incx, T *y, const int incy)
vabs: y = |x|
Definition: Vmath.hpp:352

References ASSERTL0, Nektar::SolverUtils::EquationSystem::m_fields, and Vmath::Vabs().

◆ AdddedtMaxwell()

void Nektar::SolverUtils::MMFSystem::AdddedtMaxwell ( const Array< OneD, const Array< OneD, NekDouble > > &  physarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protected

Definition at line 1364 of file MMFSystem.cpp.

1367{
1368 int nq = GetTotPoints();
1369
1370 // m_dedxi_cdot_e[m][j][n][] = de^m / d \xi^j \cdot e^n
1371 Array<OneD, NekDouble> dedtej(nq);
1372 Array<OneD, NekDouble> de1dtcdotej(nq);
1373 Array<OneD, NekDouble> de2dtcdotej(nq);
1374
1375 Array<OneD, NekDouble> normH(nq);
1376 Array<OneD, NekDouble> NormedH1(nq, 0.0);
1377 Array<OneD, NekDouble> NormednegH2(nq, 0.0);
1378
1379 Vmath::Vmul(nq, &physarray[0][0], 1, &physarray[0][0], 1, &normH[0], 1);
1380 Vmath::Vvtvp(nq, &physarray[1][0], 1, &physarray[1][0], 1, &normH[0], 1,
1381 &normH[0], 1);
1382 Vmath::Vsqrt(nq, normH, 1, normH, 1);
1383
1384 NekDouble Tol = 0.001;
1385 for (int i = 0; i < nq; ++i)
1386 {
1387 if (normH[i] > Tol)
1388 {
1389 NormedH1[i] = physarray[0][i] / normH[i];
1390 NormednegH2[i] = -1.0 * physarray[1][i] / normH[i];
1391 }
1392 }
1393
1394 for (int j = 0; j < m_shapedim; ++j)
1395 {
1396 // Compute de1 / dt \cdot ej = (-H2 de^1/d\xi1 \cdot e^j + H1 de^1/d\xi2
1397 // \cdot e^j) / sqrt{ H1^2 + H2^2 }
1398 Vmath::Vmul(nq, &NormednegH2[0], 1, &m_dedxi_cdot_e[0][0][j][0], 1,
1399 &de1dtcdotej[0], 1);
1400 Vmath::Vvtvp(nq, &NormedH1[0], 1, &m_dedxi_cdot_e[0][1][j][0], 1,
1401 &de1dtcdotej[0], 1, &de1dtcdotej[0], 1);
1402
1403 // Compute de2 / dt \cdot ej = (-H2 de2/d\xi1 \cdot e^j + H1 de2/d\xi2
1404 // \cdot e^j) / sqrt{ H1^2 + H2^2 }
1405 Vmath::Vmul(nq, &NormednegH2[0], 1, &m_dedxi_cdot_e[1][0][j][0], 1,
1406 &de2dtcdotej[0], 1);
1407 Vmath::Vvtvp(nq, &NormedH1[0], 1, &m_dedxi_cdot_e[1][1][j][0], 1,
1408 &de2dtcdotej[0], 1, &de2dtcdotej[0], 1);
1409
1410 // Add dedt component: (H1 (de1/dt) + H2 (de2/dt) ) \cdot ej
1411 Vmath::Vmul(nq, &physarray[0][0], 1, &de1dtcdotej[0], 1, &dedtej[0], 1);
1412 Vmath::Vvtvp(nq, &physarray[1][0], 1, &de2dtcdotej[0], 1, &dedtej[0], 1,
1413 &dedtej[0], 1);
1414
1415 Vmath::Neg(nq, dedtej, 1);
1416
1417 switch (m_PolType)
1418 {
1420 {
1421 if (j == 0)
1422 {
1423 Vmath::Vmul(nq, m_muvec[0], 1, dedtej, 1, dedtej, 1);
1424 }
1425
1426 else if (j == 1)
1427 {
1428 Vmath::Vmul(nq, m_muvec[1], 1, dedtej, 1, dedtej, 1);
1429 }
1430 }
1431 break;
1432
1434 {
1435 if (j == 0)
1436 {
1437 Vmath::Vmul(nq, m_epsvec[0], 1, dedtej, 1, dedtej, 1);
1438 }
1439
1440 else if (j == 1)
1441 {
1442 Vmath::Vmul(nq, m_epsvec[1], 1, dedtej, 1, dedtej, 1);
1443 }
1444 }
1445 break;
1446
1447 default:
1448 break;
1449 }
1450
1451 Vmath::Vadd(nq, &dedtej[0], 1, &outarray[j][0], 1, &outarray[j][0], 1);
1452 }
1453}
SOLVER_UTILS_EXPORT int GetTotPoints()
Array< OneD, Array< OneD, Array< OneD, Array< OneD, NekDouble > > > > m_dedxi_cdot_e
Definition: MMFSystem.h:217
Array< OneD, Array< OneD, NekDouble > > m_muvec
Definition: MMFSystem.h:210
Array< OneD, Array< OneD, NekDouble > > m_epsvec
Definition: MMFSystem.h:209
double NekDouble
void Vsqrt(int n, const T *x, const int incx, T *y, const int incy)
sqrt y = sqrt(x)
Definition: Vmath.hpp:340
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.hpp:72
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.hpp:292
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 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

References Nektar::SolverUtils::eTransElectric, Nektar::SolverUtils::eTransMagnetic, Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_dedxi_cdot_e, m_epsvec, m_muvec, m_PolType, m_shapedim, Vmath::Neg(), Vmath::Vadd(), Vmath::Vmul(), Vmath::Vsqrt(), and Vmath::Vvtvp().

Referenced by Nektar::MMFMaxwell::DoOdeRhs().

◆ AverageMaxwellFlux1D()

void Nektar::SolverUtils::MMFSystem::AverageMaxwellFlux1D ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd 
)
protected

Definition at line 1562 of file MMFSystem.cpp.

1566{
1567 int i;
1568 int nTraceNumPoints = GetTraceTotPoints();
1569 int nvar = 2;
1570
1571 // get temporary arrays
1572 Array<OneD, Array<OneD, NekDouble>> Fwd(nvar);
1573 Array<OneD, Array<OneD, NekDouble>> Bwd(nvar);
1574
1575 for (i = 0; i < nvar; ++i)
1576 {
1577 Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1578 Bwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1579 }
1580
1581 // get the physical values at the trace from the dependent variables
1582 for (i = 0; i < nvar; ++i)
1583 {
1584 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd[i], Bwd[i]);
1585 }
1586
1587 // E = 0 at the boundaries
1588 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0, "PEC");
1589
1590 // d H / d n = 0 at the boundaries
1591 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PEC");
1592
1593 for (i = 0; i < nTraceNumPoints; ++i)
1594 {
1595 numfluxFwd[0][i] = 0.5 * m_traceNormals[0][i] * (Fwd[1][i] + Bwd[1][i]);
1596 numfluxFwd[1][i] = 0.5 * m_traceNormals[0][i] * (Fwd[0][i] + Bwd[0][i]);
1597
1598 numfluxBwd[0][i] = 0.5 * m_traceNormals[0][i] * (Fwd[1][i] + Bwd[1][i]);
1599 numfluxBwd[1][i] = 0.5 * m_traceNormals[0][i] * (Fwd[0][i] + Bwd[0][i]);
1600 }
1601}
SOLVER_UTILS_EXPORT int GetTraceTotPoints()
Array< OneD, Array< OneD, NekDouble > > m_traceNormals
Array holding trace normals for DG simulations in the forwards direction.
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")
Definition: MMFSystem.cpp:835

References CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, and Nektar::SolverUtils::EquationSystem::m_traceNormals.

Referenced by NumericalMaxwellFlux().

◆ AvgAbsInt()

NekDouble Nektar::SolverUtils::MMFSystem::AvgAbsInt ( const Array< OneD, const NekDouble > &  inarray)
protected

Definition at line 2307 of file MMFSystem.cpp.

2308{
2309 int nq = m_fields[0]->GetNpoints();
2310 Array<OneD, NekDouble> Ones(nq, 1.0);
2311 Array<OneD, NekDouble> tmp(nq);
2312
2313 if (inarray.size() != nq)
2314 {
2315 ASSERTL0(false, "AvgAbsInt Error: Vector size is not correct");
2316 }
2317
2318 NekDouble jac = m_fields[0]->Integral(Ones);
2319
2320 Vmath::Vabs(nq, inarray, 1, tmp, 1);
2321 return (m_fields[0]->Integral(tmp)) / jac;
2322}

References ASSERTL0, Nektar::SolverUtils::EquationSystem::m_fields, and Vmath::Vabs().

Referenced by CheckMovingFrames(), Nektar::MMFSWE::EvaluateWaterDepth(), and Nektar::MMFSWE::TestVorticityComputation().

◆ AvgInt()

NekDouble Nektar::SolverUtils::MMFSystem::AvgInt ( const Array< OneD, const NekDouble > &  inarray)
protected

Definition at line 2292 of file MMFSystem.cpp.

2293{
2294 int nq = m_fields[0]->GetNpoints();
2295 Array<OneD, NekDouble> Ones(nq, 1.0);
2296
2297 if (inarray.size() != nq)
2298 {
2299 ASSERTL0(false, "AvgInt Error: Vector size is not correct");
2300 }
2301
2302 NekDouble jac = m_fields[0]->Integral(Ones);
2303
2304 return (m_fields[0]->Integral(inarray)) / jac;
2305}

References ASSERTL0, and Nektar::SolverUtils::EquationSystem::m_fields.

Referenced by ComputeDivCurlMF(), and Nektar::MMFAdvection::v_InitObject().

◆ BubbleSort()

void Nektar::SolverUtils::MMFSystem::BubbleSort ( Array< OneD, NekDouble > &  refarray,
Array< OneD, NekDouble > &  sortarray 
)
protected

Definition at line 2368 of file MMFSystem.cpp.

2370{
2371 int nq = refarray.size();
2372
2373 bool swapped = true;
2374 int j = 0;
2375 NekDouble tmp;
2376
2377 while (swapped)
2378 {
2379 swapped = false;
2380 j++;
2381 for (int i = 0; i < nq - j; i++)
2382 {
2383 if (refarray[i] > refarray[i + 1])
2384 {
2385 tmp = refarray[i];
2386 refarray[i] = refarray[i + 1];
2387 refarray[i + 1] = tmp;
2388
2389 tmp = sortarray[i];
2390 sortarray[i] = sortarray[i + 1];
2391 sortarray[i + 1] = tmp;
2392
2393 swapped = true;
2394 }
2395 }
2396 }
2397}

◆ CartesianToMovingframes()

Array< OneD, NekDouble > Nektar::SolverUtils::MMFSystem::CartesianToMovingframes ( const Array< OneD, const Array< OneD, NekDouble > > &  uvec,
unsigned int  field 
)
protected

Definition at line 771 of file MMFSystem.cpp.

773{
774 int nq = m_fields[0]->GetNpoints();
775
776 Array<OneD, NekDouble> outarray(nq, 0.0);
777
778 // new u0 = ( [u v] \cdot e^1 )/ |e^1|^2
779 Vmath::Vmul(nq, &m_movingframes[field][0], 1, &uvec[0][0], 1, &outarray[0],
780 1);
781 Vmath::Vvtvp(nq, &m_movingframes[field][nq], 1, &uvec[1][0], 1,
782 &outarray[0], 1, &outarray[0], 1);
783 Vmath::Vvtvp(nq, &m_movingframes[field][2 * nq], 1, &uvec[2][0], 1,
784 &outarray[0], 1, &outarray[0], 1);
785
786 return outarray;
787}
Array< OneD, Array< OneD, NekDouble > > m_movingframes
Definition: MMFSystem.h:183

References Nektar::SolverUtils::EquationSystem::m_fields, m_movingframes, Vmath::Vmul(), and Vmath::Vvtvp().

Referenced by Nektar::MMFSWE::IsolatedMountainFlow(), Nektar::MMFSWE::RossbyWave(), Nektar::MMFSWE::SteadyZonalFlow(), Nektar::MMFMaxwell::TestMaxwellSphere(), Nektar::MMFSWE::TestSWE2Dproblem(), Nektar::MMFSWE::TestVorticityComputation(), Nektar::MMFSWE::UnstableJetFlow(), and Nektar::MMFSWE::UnsteadyZonalFlow().

◆ CartesianToSpherical()

void Nektar::SolverUtils::MMFSystem::CartesianToSpherical ( const NekDouble  x0j,
const NekDouble  x1j,
const NekDouble  x2j,
NekDouble sin_varphi,
NekDouble cos_varphi,
NekDouble sin_theta,
NekDouble cos_theta 
)
protected

Definition at line 792 of file MMFSystem.cpp.

796{
797 NekDouble radius;
798 NekDouble radxy;
799 NekDouble Tol = 0.0000001;
800
801 NekDouble m_Xscale = 1.0;
802 NekDouble m_Yscale = 1.0;
803 NekDouble m_Zscale = 1.0;
804
805 radius = sqrt(x0j * x0j / (m_Xscale * m_Xscale) +
806 x1j * x1j / (m_Yscale * m_Yscale) +
807 x2j * x2j / (m_Zscale * m_Zscale));
808 radxy = sqrt(x0j * x0j / (m_Xscale * m_Xscale) +
809 x1j * x1j / (m_Yscale * m_Yscale));
810
811 if (radxy > Tol)
812 {
813 sin_varphi = x1j / (radxy * m_Yscale);
814 cos_varphi = x0j / (radxy * m_Xscale);
815 }
816
817 else
818 {
819 sin_varphi = 0.0;
820 if (x2j > 0)
821 {
822 cos_varphi = 1.0;
823 }
824
825 else
826 {
827 cos_varphi = -1.0;
828 }
829 }
830
831 sin_theta = x2j / (radius * m_Zscale);
832 cos_theta = radxy / radius;
833}
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:294

References tinysimd::sqrt().

Referenced by Nektar::MMFAdvection::EvaluateAdvectionVelocity(), Nektar::MMFMaxwell::EvaluateCoriolis(), Nektar::MMFSWE::EvaluateCoriolisForZonalFlow(), Nektar::MMFSWE::EvaluateStandardCoriolis(), Nektar::MMFSWE::EvaluateWaterDepth(), Nektar::MMFSWE::IsolatedMountainFlow(), Nektar::MMFSWE::RossbyWave(), Nektar::MMFSWE::SteadyZonalFlow(), Nektar::MMFMaxwell::TestMaxwellSphere(), Nektar::MMFSWE::TestVorticityComputation(), Nektar::MMFSWE::UnstableJetFlow(), and Nektar::MMFSWE::UnsteadyZonalFlow().

◆ CheckMovingFrames()

void Nektar::SolverUtils::MMFSystem::CheckMovingFrames ( const Array< OneD, const Array< OneD, NekDouble > > &  movingframes)
protected

Definition at line 230 of file MMFSystem.cpp.

232{
233 NekDouble t1x, t1y, t1z, t2x, t2y, t2z, t3x, t3y, t3z;
234 NekDouble dot12 = 0.0, dot23 = 0.0, dot31 = 0.0;
235 NekDouble Tol = 0.0001;
236
237 int nq = m_fields[0]->GetNpoints();
238
239 for (int i = 0; i < nq; ++i)
240 {
241 t1x = movingframes[0][i];
242 t1y = movingframes[0][i + nq];
243 t1z = movingframes[0][i + 2 * nq];
244
245 t2x = movingframes[1][i];
246 t2y = movingframes[1][i + nq];
247 t2z = movingframes[1][i + 2 * nq];
248
249 t3x = movingframes[2][i];
250 t3y = movingframes[2][i + nq];
251 t3z = movingframes[2][i + 2 * nq];
252
253 dot12 = t1x * t2x + t1y * t2y + t1z * t2z;
254 dot23 = t2x * t3x + t2y * t3y + t2z * t3z;
255 dot31 = t3x * t1x + t3y * t1y + t3z * t1z;
256 }
257
258 std::cout << "======================================================"
259 << std::endl;
260 std::cout << "======================================================"
261 << std::endl;
262 std::cout << "*** The first moving frame is alinged along"
263 << SpatialDomains::GeomMMFMap[m_MMFdir] << std::endl;
264
265 Array<OneD, NekDouble> tmpx(nq), tmpy(nq), tmpz(nq);
266
267 Vmath::Vcopy(nq, &movingframes[0][0], 1, &tmpx[0], 1);
268 Vmath::Vcopy(nq, &movingframes[0][nq], 1, &tmpy[0], 1);
269 Vmath::Vcopy(nq, &movingframes[0][2 * nq], 1, &tmpz[0], 1);
270 std::cout << nq << " , "
271 << "*** Avg MF1 = ( " << AvgAbsInt(tmpx) << " , "
272 << AvgAbsInt(tmpy) << " , " << AvgAbsInt(tmpz) << " ) "
273 << std::endl;
274
275 Vmath::Vcopy(nq, &movingframes[1][0], 1, &tmpx[0], 1);
276 Vmath::Vcopy(nq, &movingframes[1][nq], 1, &tmpy[0], 1);
277 Vmath::Vcopy(nq, &movingframes[1][2 * nq], 1, &tmpz[0], 1);
278 std::cout << "*** Avg MF2 = ( " << AvgAbsInt(tmpx) << " , "
279 << AvgAbsInt(tmpy) << " , " << AvgAbsInt(tmpz) << " ) "
280 << std::endl;
281
282 if (m_shapedim == 3)
283 {
284 Vmath::Vcopy(nq, &movingframes[2][0], 1, &tmpx[0], 1);
285 Vmath::Vcopy(nq, &movingframes[2][nq], 1, &tmpy[0], 1);
286 Vmath::Vcopy(nq, &movingframes[2][2 * nq], 1, &tmpz[0], 1);
287 std::cout << "*** Avg MF3 = ( " << AvgAbsInt(tmpx) << " , "
288 << AvgAbsInt(tmpy) << " , " << AvgAbsInt(tmpz) << " ) "
289 << std::endl;
290 }
291
292 if ((fabs(dot12) + fabs(dot23) + fabs(dot31)) < Tol)
293 {
294 std::cout << "*** Moving frames are Orthogonal" << std::endl;
295 }
296
297 else
298 {
299 std::cout << "*** Moving frames are NOT Orthogonal" << std::endl;
300 }
301
302 Array<OneD, NekDouble> tmp;
303 for (int j = 0; j < m_shapedim; ++j)
304 {
305 tmp = Array<OneD, NekDouble>(nq, 0.0);
306 for (int k = 0; k < m_spacedim; ++k)
307 {
308 Vmath::Vvtvp(nq, &movingframes[j][k * nq], 1,
309 &movingframes[j][k * nq], 1, &tmp[0], 1, &tmp[0], 1);
310 }
311 Vmath::Vsqrt(nq, tmp, 1, tmp, 1);
312 std::cout << "*** Avg. Magnitude of MF" << j << " = " << AvgAbsInt(tmp)
313 << std::endl;
314 }
315}
int m_spacedim
Spatial dimension (>= expansion dim).
SpatialDomains::GeomMMF m_MMFdir
Definition: MMFSystem.h:219
SOLVER_UTILS_EXPORT NekDouble AvgAbsInt(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2307
const char *const GeomMMFMap[]
Session file names associated with tangent principle directions.
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.hpp:825

References AvgAbsInt(), Nektar::SpatialDomains::GeomMMFMap, Nektar::SolverUtils::EquationSystem::m_fields, m_MMFdir, m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vcopy(), Vmath::Vsqrt(), and Vmath::Vvtvp().

Referenced by SetUpMovingFrames().

◆ ComputeCurl()

void Nektar::SolverUtils::MMFSystem::ComputeCurl ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray 
)
protected

Definition at line 731 of file MMFSystem.cpp.

735{
736 int nq = inarray[0].size();
737
738 Array<OneD, NekDouble> tmpx, tmpy, tmpz;
739 Array<OneD, NekDouble> Dtmpzdx, Dtmpydx, Dtmpxdy, Dtmpzdy, Dtmpxdz, Dtmpydz;
740
741 tmpx = Array<OneD, NekDouble>(nq);
742 tmpy = Array<OneD, NekDouble>(nq);
743 tmpz = Array<OneD, NekDouble>(nq);
744
745 Dtmpzdx = Array<OneD, NekDouble>(nq);
746 Dtmpydx = Array<OneD, NekDouble>(nq);
747 Dtmpxdy = Array<OneD, NekDouble>(nq);
748 Dtmpzdy = Array<OneD, NekDouble>(nq);
749 Dtmpxdz = Array<OneD, NekDouble>(nq);
750 Dtmpydz = Array<OneD, NekDouble>(nq);
751
752 for (int k = 0; k < m_spacedim; ++k)
753 {
754 Vmath::Vcopy(nq, &inarray[0][0], 1, &tmpx[0], 1);
755 Vmath::Vcopy(nq, &inarray[1][0], 1, &tmpy[0], 1);
756 Vmath::Vcopy(nq, &inarray[2][0], 1, &tmpz[0], 1);
757
758 m_fields[0]->PhysDeriv(0, tmpz, Dtmpzdx);
759 m_fields[0]->PhysDeriv(0, tmpy, Dtmpydx);
760 m_fields[0]->PhysDeriv(1, tmpx, Dtmpxdy);
761 m_fields[0]->PhysDeriv(1, tmpz, Dtmpzdy);
762 m_fields[0]->PhysDeriv(2, tmpx, Dtmpxdz);
763 m_fields[0]->PhysDeriv(2, tmpy, Dtmpydz);
764
765 Vmath::Vsub(nq, &Dtmpzdy[0], 1, &Dtmpydz[0], 1, &outarray[0][0], 1);
766 Vmath::Vsub(nq, &Dtmpxdz[0], 1, &Dtmpzdx[0], 1, &outarray[1][0], 1);
767 Vmath::Vsub(nq, &Dtmpydx[0], 1, &Dtmpxdy[0], 1, &outarray[2][0], 1);
768 }
769}
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::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vcopy(), and Vmath::Vsub().

Referenced by ComputeDivCurlMF().

◆ Computedemdxicdote()

void Nektar::SolverUtils::MMFSystem::Computedemdxicdote ( )
protected

Definition at line 1274 of file MMFSystem.cpp.

1275{
1276 int MFdim = 3;
1277 int nq = GetTotPoints();
1278
1279 // Initialization
1281 Array<OneD, Array<OneD, Array<OneD, Array<OneD, NekDouble>>>>(MFdim);
1282 for (int indm = 0; indm < MFdim; ++indm)
1283 {
1284 m_dedxi_cdot_e[indm] =
1285 Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
1286 for (int indj = 0; indj < MFdim; ++indj)
1287 {
1288 m_dedxi_cdot_e[indm][indj] =
1289 Array<OneD, Array<OneD, NekDouble>>(MFdim);
1290 for (int indn = 0; indn < MFdim; ++indn)
1291 {
1292 m_dedxi_cdot_e[indm][indj][indn] =
1293 Array<OneD, NekDouble>(nq, 0.0);
1294 }
1295 }
1296 }
1297
1298 Array<OneD, NekDouble> tmp(nq);
1299 Array<OneD, NekDouble> tmpderiv(nq);
1300 Array<OneD, NekDouble> dedt;
1301 for (int indm = 0; indm < MFdim; ++indm)
1302 {
1303 for (int indj = 0; indj < MFdim; ++indj)
1304 {
1305 for (int indn = 0; indn < MFdim; ++indn)
1306 {
1307 dedt = Array<OneD, NekDouble>(nq, 0.0);
1308 for (int k = 0; k < m_spacedim; ++k)
1309 {
1310 // Compute d e^m / d \xi_j cdot e^n
1311 Vmath::Vcopy(nq, &m_movingframes[indm][k * nq], 1, &tmp[0],
1312 1);
1313 m_fields[0]->PhysDirectionalDeriv(m_movingframes[indj], tmp,
1314 tmpderiv);
1315
1316 Vmath::Vvtvp(nq, &tmpderiv[0], 1,
1317 &m_movingframes[indn][k * nq], 1, &dedt[0], 1,
1318 &dedt[0], 1);
1319 }
1320
1321 Vmath::Vcopy(nq, &dedt[0], 1,
1322 &m_dedxi_cdot_e[indm][indj][indn][0], 1);
1323 }
1324 }
1325 }
1326
1327 int indx = 0;
1328 std::cout << "*** m_dedxi_cdot_e[0]/dxi1 = ( "
1329 << RootMeanSquare(m_dedxi_cdot_e[indx][0][0]) << " , "
1330 << RootMeanSquare(m_dedxi_cdot_e[indx][0][1]) << " , "
1331 << RootMeanSquare(m_dedxi_cdot_e[indx][0][2]) << " )_1, "
1332 << std::endl;
1333 std::cout << "*** m_dedxi_cdot_e[0]/dxi2 = ( "
1334 << RootMeanSquare(m_dedxi_cdot_e[indx][1][0]) << " , "
1335 << RootMeanSquare(m_dedxi_cdot_e[indx][1][1]) << " , "
1336 << RootMeanSquare(m_dedxi_cdot_e[indx][1][2]) << " )_2 "
1337 << std::endl;
1338
1339 indx = 1;
1340 std::cout << "*** m_dedxi_cdot_e[1]/dxi1 = ( "
1341 << RootMeanSquare(m_dedxi_cdot_e[indx][0][0]) << " , "
1342 << RootMeanSquare(m_dedxi_cdot_e[indx][0][1]) << " , "
1343 << RootMeanSquare(m_dedxi_cdot_e[indx][0][2]) << " )_1, "
1344 << std::endl;
1345 std::cout << "*** m_dedxi_cdot_e[1]/dxi2 = ( "
1346 << RootMeanSquare(m_dedxi_cdot_e[indx][1][0]) << " , "
1347 << RootMeanSquare(m_dedxi_cdot_e[indx][1][1]) << " , "
1348 << RootMeanSquare(m_dedxi_cdot_e[indx][1][2]) << " )_2 "
1349 << std::endl;
1350
1351 indx = 2;
1352 std::cout << "*** m_dedxi_cdot_e[2]/dxi1 = ( "
1353 << RootMeanSquare(m_dedxi_cdot_e[indx][0][0]) << " , "
1354 << RootMeanSquare(m_dedxi_cdot_e[indx][0][1]) << " , "
1355 << RootMeanSquare(m_dedxi_cdot_e[indx][0][2]) << " )_1, "
1356 << std::endl;
1357 std::cout << "*** m_dedxi_cdot_e[2]/dxi2 = ( "
1358 << RootMeanSquare(m_dedxi_cdot_e[indx][1][0]) << " , "
1359 << RootMeanSquare(m_dedxi_cdot_e[indx][1][1]) << " , "
1360 << RootMeanSquare(m_dedxi_cdot_e[indx][1][2]) << " )_2 "
1361 << std::endl;
1362}
SOLVER_UTILS_EXPORT NekDouble RootMeanSquare(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2338

References Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_dedxi_cdot_e, Nektar::SolverUtils::EquationSystem::m_fields, m_movingframes, Nektar::SolverUtils::EquationSystem::m_spacedim, RootMeanSquare(), Vmath::Vcopy(), and Vmath::Vvtvp().

Referenced by Nektar::MMFMaxwell::v_InitObject().

◆ ComputeDivCurlMF()

void Nektar::SolverUtils::MMFSystem::ComputeDivCurlMF ( )
protected

Definition at line 442 of file MMFSystem.cpp.

443{
444 int nq = m_fields[0]->GetNpoints();
445 int MMFdim = 3;
446
447 Array<OneD, NekDouble> tmp(nq);
448 Array<OneD, NekDouble> Dtmp(nq);
449
450 m_DivMF = Array<OneD, Array<OneD, NekDouble>>(MMFdim);
451 for (int j = 0; j < MMFdim; ++j)
452 {
453 m_DivMF[j] = Array<OneD, NekDouble>(nq, 0.0);
454 for (int k = 0; k < m_spacedim; ++k)
455 {
456 Vmath::Vcopy(nq, &m_movingframes[j][k * nq], 1, &tmp[0], 1);
457
458 m_fields[0]->PhysDeriv(k, tmp, Dtmp);
459 Vmath::Vadd(nq, &Dtmp[0], 1, &m_DivMF[j][0], 1, &m_DivMF[j][0], 1);
460 }
461 }
462
463 std::cout << "*** Divergence of MF1 = " << AvgInt(m_DivMF[0])
464 << ", MF2 = " << AvgInt(m_DivMF[1])
465 << ", MF3 = " << AvgInt(m_DivMF[2]) << std::endl;
466
467 // Compute Curl of MF: CurlMF[i][j] = (\nabla \times e^i) cdot e^j
468 m_CurlMF = Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MMFdim);
469 for (int i = 0; i < MMFdim; ++i)
470 {
471 m_CurlMF[i] = Array<OneD, Array<OneD, NekDouble>>(MMFdim);
472
473 for (int j = 0; j < MMFdim; ++j)
474 {
475 m_CurlMF[i][j] = Array<OneD, NekDouble>(nq, 0.0);
476 }
477 }
478
479 Array<OneD, Array<OneD, NekDouble>> MFtmp(m_spacedim);
480 Array<OneD, Array<OneD, NekDouble>> CurlMFtmp(m_spacedim);
481
482 for (int i = 0; i < m_spacedim; ++i)
483 {
484 MFtmp[i] = Array<OneD, NekDouble>(nq);
485 CurlMFtmp[i] = Array<OneD, NekDouble>(nq);
486 }
487
488 for (int dir = 0; dir < MMFdim; dir++)
489 {
490 for (int i = 0; i < m_spacedim; ++i)
491 {
492 Vmath::Vcopy(nq, &m_movingframes[dir][i * nq], 1, &MFtmp[i][0], 1);
493 }
494
495 ComputeCurl(MFtmp, CurlMFtmp);
496
497 for (int j = 0; j < MMFdim; ++j)
498 {
499 for (int i = 0; i < m_spacedim; ++i)
500 {
501 Vmath::Vvtvp(nq, &m_movingframes[j][i * nq], 1,
502 &CurlMFtmp[i][0], 1, &m_CurlMF[dir][j][0], 1,
503 &m_CurlMF[dir][j][0], 1);
504 }
505 }
506 }
507
508 std::cout << "*** Curl of MF1 = ( " << AvgInt(m_CurlMF[0][0]) << " , "
509 << AvgInt(m_CurlMF[0][1]) << " , " << AvgInt(m_CurlMF[0][2])
510 << " ) " << std::endl;
511 std::cout << "*** Curl of MF2 = ( " << AvgInt(m_CurlMF[1][0]) << " , "
512 << AvgInt(m_CurlMF[1][1]) << " , " << AvgInt(m_CurlMF[1][2])
513 << " ) " << std::endl;
514 std::cout << "*** Curl of MF3 = ( " << AvgInt(m_CurlMF[2][0]) << " , "
515 << AvgInt(m_CurlMF[2][1]) << " , " << AvgInt(m_CurlMF[2][2])
516 << " ) " << std::endl;
517}
SOLVER_UTILS_EXPORT void ComputeCurl(const Array< OneD, const Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray)
Definition: MMFSystem.cpp:731
Array< OneD, Array< OneD, NekDouble > > m_DivMF
Definition: MMFSystem.h:192
SOLVER_UTILS_EXPORT NekDouble AvgInt(const Array< OneD, const NekDouble > &inarray)
Definition: MMFSystem.cpp:2292
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_CurlMF
Definition: MMFSystem.h:193

References AvgInt(), ComputeCurl(), m_CurlMF, m_DivMF, Nektar::SolverUtils::EquationSystem::m_fields, m_movingframes, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vadd(), Vmath::Vcopy(), and Vmath::Vvtvp().

Referenced by MMFInitObject().

◆ ComputeMFtrace()

void Nektar::SolverUtils::MMFSystem::ComputeMFtrace ( )
protected

Definition at line 519 of file MMFSystem.cpp.

520{
521 int MFdim = 3;
522
523 int nq = m_fields[0]->GetNpoints();
524 int nTraceNumPoints = GetTraceTotPoints();
525
526 Array<OneD, NekDouble> tmp(nq);
527 Array<OneD, NekDouble> Fwdtmp(nq);
528 Array<OneD, NekDouble> Bwdtmp(nq);
529
530 m_MFtraceFwd = Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
531 m_MFtraceBwd = Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
532
533 for (int j = 0; j < MFdim; ++j)
534 {
535 m_MFtraceFwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
536 m_MFtraceBwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
537 }
538
539 // m_MFtraceFwd[0] = e^1_{Fwd}, m_MFtraceFwd[1] = e^2_{Fwd}
540 for (int j = 0; j < MFdim; ++j)
541 {
542 for (int i = 0; i < m_spacedim; ++i)
543 {
544 m_MFtraceFwd[j][i] = Array<OneD, NekDouble>(nTraceNumPoints);
545 m_MFtraceBwd[j][i] = Array<OneD, NekDouble>(nTraceNumPoints);
546
547 Vmath::Vcopy(nq, &m_movingframes[j][i * nq], 1, &tmp[0], 1);
548
549 m_fields[0]->GetFwdBwdTracePhys(tmp, Fwdtmp, Bwdtmp);
550
551 CopyBoundaryTrace(Fwdtmp, Bwdtmp, eFwdEQBwd);
552
553 Vmath::Vcopy(nTraceNumPoints, &Fwdtmp[0], 1, &m_MFtraceFwd[j][i][0],
554 1);
555 Vmath::Vcopy(nTraceNumPoints, &Bwdtmp[0], 1, &m_MFtraceBwd[j][i][0],
556 1);
557 }
558 }
559
560 std::cout << "*** MFtraceFwd = ( " << VectorAvgMagnitude(m_MFtraceFwd[0])
561 << " , " << VectorAvgMagnitude(m_MFtraceFwd[1]) << " , "
562 << VectorAvgMagnitude(m_MFtraceFwd[2]) << " ) " << std::endl;
563 std::cout << "*** MFtraceBwd = ( " << VectorAvgMagnitude(m_MFtraceBwd[0])
564 << " , " << VectorAvgMagnitude(m_MFtraceBwd[1]) << " , "
565 << VectorAvgMagnitude(m_MFtraceBwd[2]) << " ) " << std::endl;
566}
SOLVER_UTILS_EXPORT NekDouble VectorAvgMagnitude(const Array< OneD, const Array< OneD, NekDouble > > &inarray)
Definition: MMFSystem.cpp:2352
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_MFtraceFwd
Definition: MMFSystem.h:196
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_MFtraceBwd
Definition: MMFSystem.h:197

References CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_MFtraceBwd, m_MFtraceFwd, m_movingframes, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vcopy(), and VectorAvgMagnitude().

Referenced by MMFInitObject().

◆ ComputencdotMF()

void Nektar::SolverUtils::MMFSystem::ComputencdotMF ( )
protected

Definition at line 317 of file MMFSystem.cpp.

318{
319 int nq = m_fields[0]->GetNpoints();
320 int nTracePointsTot = GetTraceNpoints();
321
322 // Compute MFjFwd and MFjBwd
323 Array<OneD, NekDouble> tmp(nq);
324
325 Array<OneD, Array<OneD, Array<OneD, NekDouble>>> MFtraceFwd(m_shapedim);
326 Array<OneD, Array<OneD, Array<OneD, NekDouble>>> MFtraceBwd(m_shapedim);
327 Array<OneD, Array<OneD, NekDouble>> SurfaceNormalFwd;
328 Array<OneD, Array<OneD, NekDouble>> SurfaceNormalBwd;
329
330 for (int j = 0; j < m_shapedim; ++j)
331 {
332 MFtraceFwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
333 MFtraceBwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
334
335 SurfaceNormalFwd = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
336 SurfaceNormalBwd = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
337
338 for (int k = 0; k < m_spacedim; ++k)
339 {
340 MFtraceFwd[j][k] = Array<OneD, NekDouble>(nTracePointsTot);
341 MFtraceBwd[j][k] = Array<OneD, NekDouble>(nTracePointsTot);
342
343 SurfaceNormalFwd[k] = Array<OneD, NekDouble>(nTracePointsTot);
344 SurfaceNormalBwd[k] = Array<OneD, NekDouble>(nTracePointsTot);
345
346 Vmath::Vcopy(nq, &m_movingframes[j][k * nq], 1, &tmp[0], 1);
347
348 m_fields[0]->GetFwdBwdTracePhys(tmp, MFtraceFwd[j][k],
349 MFtraceBwd[j][k]);
350
351 CopyBoundaryTrace(MFtraceFwd[j][k], MFtraceBwd[j][k], eFwdEQBwd);
352 }
353 }
354
355 VectorCrossProd(MFtraceFwd[0], MFtraceFwd[1], SurfaceNormalFwd);
356 VectorCrossProd(MFtraceBwd[0], MFtraceBwd[1], SurfaceNormalBwd);
357
358 // Compute n \times e^i
359 m_ncdotMFFwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
360 m_ncdotMFBwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
361 for (int j = 0; j < m_shapedim; ++j)
362 {
363 m_ncdotMFFwd[j] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
364 m_ncdotMFBwd[j] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
365
366 VectorDotProd(m_traceNormals, MFtraceFwd[j], m_ncdotMFFwd[j]);
367 VectorDotProd(m_traceNormals, MFtraceBwd[j], m_ncdotMFBwd[j]);
368 }
369
370 // Compute n^{\perp} \times e^i
371 Array<OneD, Array<OneD, NekDouble>> SurfaceNormal(m_spacedim);
372 Array<OneD, Array<OneD, NekDouble>> Tracevector(m_spacedim);
373 for (int k = 0; k < m_spacedim; k++)
374 {
375 SurfaceNormal[k] = Array<OneD, NekDouble>(nTracePointsTot);
376 Tracevector[k] = Array<OneD, NekDouble>(nTracePointsTot);
377
378 Vmath::Vcopy(nTracePointsTot, &m_MFtraceFwd[2][k][0], 1,
379 &SurfaceNormal[k][0], 1);
380 }
381
382 VectorCrossProd(m_traceNormals, SurfaceNormal, Tracevector);
383
384 m_nperpcdotMFFwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
385 m_nperpcdotMFBwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
386 for (int j = 0; j < m_shapedim; ++j)
387 {
388 m_nperpcdotMFFwd[j] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
389 m_nperpcdotMFBwd[j] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
390
391 VectorDotProd(Tracevector, MFtraceFwd[j], m_nperpcdotMFFwd[j]);
392 VectorDotProd(Tracevector, MFtraceBwd[j], m_nperpcdotMFBwd[j]);
393 }
394
395 if (m_shapedim == 2)
396 {
397 std::cout << "*** m_ncdotMFFwd = ( " << RootMeanSquare(m_ncdotMFFwd[0])
398 << " , " << RootMeanSquare(m_ncdotMFFwd[1]) << " ) "
399 << std::endl;
400 std::cout << "*** m_ncdotMFBwd = ( " << RootMeanSquare(m_ncdotMFBwd[0])
401 << " , " << RootMeanSquare(m_ncdotMFBwd[1]) << " ) "
402 << std::endl;
403
404 std::cout << "*** m_nperpcdotMFFwd = ( "
405 << RootMeanSquare(m_nperpcdotMFFwd[0]) << " , "
406 << RootMeanSquare(m_nperpcdotMFFwd[1]) << " ) " << std::endl;
407 std::cout << "*** m_nperpcdotMFBwd = ( "
408 << RootMeanSquare(m_nperpcdotMFBwd[0]) << " , "
409 << RootMeanSquare(m_nperpcdotMFBwd[1]) << " ) " << std::endl;
410 }
411
412 else if (m_shapedim == 3)
413 {
414 std::cout << "*** m_ncdotMFFwd = ( "
415 << Vmath::Vsum(nTracePointsTot, m_ncdotMFFwd[0], 1) << " , "
416 << Vmath::Vsum(nTracePointsTot, m_ncdotMFFwd[1], 1) << " , "
417 << Vmath::Vsum(nTracePointsTot, m_ncdotMFFwd[2], 1) << " ) "
418 << std::endl;
419 std::cout << "*** m_ncdotMFBwd = ( "
420 << Vmath::Vsum(nTracePointsTot, m_ncdotMFBwd[0], 1) << " , "
421 << Vmath::Vsum(nTracePointsTot, m_ncdotMFBwd[1], 1) << " , "
422 << Vmath::Vsum(nTracePointsTot, m_ncdotMFBwd[2], 1) << " ) "
423 << std::endl;
424
425 std::cout << "*** m_nperpcdotMFFwd = ( "
426 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFFwd[0], 1)
427 << " , "
428 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFFwd[1], 1)
429 << " , "
430 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFFwd[2], 1)
431 << " ) " << std::endl;
432 std::cout << "*** m_nperpcdotMFBwd = ( "
433 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFBwd[0], 1)
434 << " , "
435 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFBwd[1], 1)
436 << " , "
437 << Vmath::Vsum(nTracePointsTot, m_nperpcdotMFBwd[2], 1)
438 << " ) " << std::endl;
439 }
440}
SOLVER_UTILS_EXPORT int GetTraceNpoints()
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)
Definition: MMFSystem.cpp:673
Array< OneD, Array< OneD, NekDouble > > m_nperpcdotMFFwd
Definition: MMFSystem.h:189
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
Array< OneD, Array< OneD, NekDouble > > m_nperpcdotMFBwd
Definition: MMFSystem.h:190
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFBwd
Definition: MMFSystem.h:187
Array< OneD, Array< OneD, NekDouble > > m_ncdotMFFwd
Definition: MMFSystem.h:186
T Vsum(int n, const T *x, const int incx)
Subtract return sum(x)
Definition: Vmath.hpp:608

References CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_MFtraceFwd, m_movingframes, m_ncdotMFBwd, m_ncdotMFFwd, m_nperpcdotMFBwd, m_nperpcdotMFFwd, m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, Nektar::SolverUtils::EquationSystem::m_traceNormals, RootMeanSquare(), Vmath::Vcopy(), VectorCrossProd(), VectorDotProd(), and Vmath::Vsum().

Referenced by MMFInitObject(), and Nektar::MMFAdvection::v_InitObject().

◆ ComputeNtimesF12()

void Nektar::SolverUtils::MMFSystem::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 
)
protected

Definition at line 945 of file MMFSystem.cpp.

953{
954 int nTraceNumPoints = GetTraceTotPoints();
955
956 NekDouble tmpFwd, tmpBwd, Aver1, Aver2, HFwdk, HBwdk;
957
958 Array<OneD, NekDouble> z1HAver(m_spacedim);
959 Array<OneD, NekDouble> z2HAver(m_spacedim);
960 Array<OneD, NekDouble> n1e1(m_spacedim);
961 Array<OneD, NekDouble> n2e2(m_spacedim);
962
963 Array<OneD, NekDouble> n1e1_times_z1HAver(m_spacedim);
964 Array<OneD, NekDouble> n2e2_times_z2HAver(m_spacedim);
965
966 for (int i = 0; i < nTraceNumPoints; ++i)
967 {
968 Aver1 = 0.5 * (im1Fwd[i] + im1Bwd[i]);
969 Aver2 = 0.5 * (im2Fwd[i] + im2Bwd[i]);
970
971 for (int k = 0; k < m_spacedim; k++)
972 {
973 // Compute \vec{HFwd} and \vec{HBwd}
974 HFwdk = Fwd[0][i] * m_MFtraceFwd[0][k][i] +
975 Fwd[1][i] * m_MFtraceFwd[1][k][i];
976 HBwdk = Bwd[0][i] * m_MFtraceBwd[0][k][i] +
977 Bwd[1][i] * m_MFtraceBwd[1][k][i];
978
979 // Compute z_i {{ \vec{H} }}
980 z1HAver[k] = 0.5 * (im1Fwd[i] * HFwdk + im1Bwd[i] * HBwdk);
981 z2HAver[k] = 0.5 * (im2Fwd[i] * HFwdk + im2Bwd[i] * HBwdk);
982
983 // Choose e^i for the one in anisotropy region
984 n1e1[k] = m_ncdotMFFwd[0][i] * m_MFtraceFwd[0][k][i];
985 n2e2[k] = m_ncdotMFFwd[1][i] * m_MFtraceFwd[1][k][i];
986 }
987
988 // Compute n1e1 \times z1HAver and n2e2 \times z2HAver
989 VectorCrossProd(n1e1, z1HAver, n1e1_times_z1HAver);
990 VectorCrossProd(n2e2, z2HAver, n2e2_times_z2HAver);
991
992 // e^3 \cdot ( n1e1 \times z1HAver + n2e2 \times z2HAver)
993 tmpFwd = 0.0;
994 tmpBwd = 0.0;
995 for (int k = 0; k < m_spacedim; k++)
996 {
997 tmpFwd += m_MFtraceFwd[2][k][i] * (n1e1_times_z1HAver[k] / Aver1 +
998 n2e2_times_z2HAver[k] / Aver2);
999 tmpBwd += m_MFtraceBwd[2][k][i] * (n1e1_times_z1HAver[k] / Aver1 +
1000 n2e2_times_z2HAver[k] / Aver2);
1001 }
1002
1003 outarrayFwd[i] = tmpFwd;
1004 outarrayBwd[i] = tmpBwd;
1005 }
1006}

References Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), m_MFtraceBwd, m_MFtraceFwd, m_ncdotMFFwd, Nektar::SolverUtils::EquationSystem::m_spacedim, and VectorCrossProd().

Referenced by NumericalMaxwellFluxTE(), and NumericalMaxwellFluxTM().

◆ ComputeNtimesFz()

void Nektar::SolverUtils::MMFSystem::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 
)
protected

Definition at line 909 of file MMFSystem.cpp.

917{
918 int nTraceNumPoints = GetTraceTotPoints();
919
920 NekDouble tmpFwd, tmpBwd;
921 NekDouble Aver, ntimesz;
922
923 for (int i = 0; i < nTraceNumPoints; ++i)
924 {
925 Aver = 0.5 * (imFwd[i] + imBwd[i]);
926
927 tmpFwd = 0.0;
928 tmpBwd = 0.0;
929 for (int k = 0; k < m_spacedim; ++k)
930 {
931 ntimesz = 0.5 * (imFwd[i] * Fwd[2][i] + imBwd[i] * Bwd[2][i]);
932
933 tmpFwd +=
934 m_MFtraceFwd[dir][k][i] * m_ntimesMFFwd[2][k][i] * ntimesz;
935 tmpBwd +=
936 m_MFtraceBwd[dir][k][i] * m_ntimesMFBwd[2][k][i] * ntimesz;
937 }
938
939 outarrayFwd[i] = tmpFwd / Aver;
940 outarrayBwd[i] = tmpBwd / Aver;
941 }
942}
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimesMFBwd
Definition: MMFSystem.h:200
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimesMFFwd
Definition: MMFSystem.h:199

References Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), m_MFtraceBwd, m_MFtraceFwd, m_ntimesMFBwd, m_ntimesMFFwd, and Nektar::SolverUtils::EquationSystem::m_spacedim.

Referenced by NumericalMaxwellFluxTE(), and NumericalMaxwellFluxTM().

◆ ComputeNtimesMF()

void Nektar::SolverUtils::MMFSystem::ComputeNtimesMF ( )
protected

Definition at line 615 of file MMFSystem.cpp.

616{
617 int MFdim = 3;
618 int nTracePointsTot = GetTraceNpoints();
619
620 m_ntimesMFFwd = Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
621 m_ntimesMFBwd = Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
623 Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
625 Array<OneD, Array<OneD, Array<OneD, NekDouble>>>(MFdim);
626 for (int j = 0; j < MFdim; ++j)
627 {
628 m_ntimesMFFwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
629 m_ntimesMFBwd[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
631 Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
633 Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
634
635 for (int k = 0; k < m_spacedim; ++k)
636 {
637 m_ntimesMFFwd[j][k] = Array<OneD, NekDouble>(nTracePointsTot);
638 m_ntimesMFBwd[j][k] = Array<OneD, NekDouble>(nTracePointsTot);
640 Array<OneD, NekDouble>(nTracePointsTot);
642 Array<OneD, NekDouble>(nTracePointsTot);
643 }
644
651 }
652
653 std::cout << "*** m_ntimesMFFwd = ( "
654 << VectorAvgMagnitude(m_ntimesMFFwd[0]) << " , "
655 << VectorAvgMagnitude(m_ntimesMFFwd[1]) << " , "
656 << VectorAvgMagnitude(m_ntimesMFFwd[2]) << " ) " << std::endl;
657 std::cout << "*** m_ntimesMFBwd = ( "
658 << VectorAvgMagnitude(m_ntimesMFBwd[0]) << " , "
659 << VectorAvgMagnitude(m_ntimesMFBwd[1]) << " , "
660 << VectorAvgMagnitude(m_ntimesMFBwd[2]) << " ) " << std::endl;
661 std::cout << "*** m_ntimes_ntimesMFFwd = ( "
665 << std::endl;
666 std::cout << "*** m_ntimes_ntimesMFBwd = ( "
670 << std::endl;
671}
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimes_ntimesMFBwd
Definition: MMFSystem.h:202
Array< OneD, Array< OneD, Array< OneD, NekDouble > > > m_ntimes_ntimesMFFwd
Definition: MMFSystem.h:201

References Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), m_MFtraceBwd, m_MFtraceFwd, m_ntimes_ntimesMFBwd, m_ntimes_ntimesMFFwd, m_ntimesMFBwd, m_ntimesMFFwd, Nektar::SolverUtils::EquationSystem::m_spacedim, Nektar::SolverUtils::EquationSystem::m_traceNormals, VectorAvgMagnitude(), and VectorCrossProd().

Referenced by Nektar::MMFMaxwell::v_InitObject().

◆ ComputeNtimestimesdF12()

void Nektar::SolverUtils::MMFSystem::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 
)
protected

Definition at line 1074 of file MMFSystem.cpp.

1082{
1083 int nTraceNumPoints = GetTraceTotPoints();
1084
1085 Array<OneD, NekDouble> directFwd(nTraceNumPoints);
1086 Array<OneD, NekDouble> directBwd(nTraceNumPoints);
1087
1088 NekDouble Aver1, Aver2;
1089 for (int i = 0; i < nTraceNumPoints; ++i)
1090 {
1091 Aver1 = im1Fwd[i] + im1Bwd[i];
1092 Aver2 = im2Fwd[i] + im2Bwd[i];
1093
1094 outarrayFwd[i] =
1095 -m_alpha * (Fwd[2][i] - Bwd[2][i]) * (1.0 / Aver1 + 1.0 / Aver2);
1096 outarrayBwd[i] =
1097 -m_alpha * (Fwd[2][i] - Bwd[2][i]) * (1.0 / Aver1 + 1.0 / Aver2);
1098 }
1099}

References Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), and m_alpha.

Referenced by NumericalMaxwellFluxTE(), and NumericalMaxwellFluxTM().

◆ ComputeNtimestimesdFz()

void Nektar::SolverUtils::MMFSystem::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 
)
protected

Definition at line 1008 of file MMFSystem.cpp.

1014{
1015 int nTraceNumPoints = GetTraceTotPoints();
1016
1017 Array<OneD, NekDouble> dH(m_spacedim);
1018 Array<OneD, NekDouble> nFwd(m_spacedim);
1019
1020 Array<OneD, NekDouble> eiFwd(m_spacedim);
1021 Array<OneD, NekDouble> eiBwd(m_spacedim);
1022
1023 Array<OneD, NekDouble> eitimesdHFwd(m_spacedim);
1024 Array<OneD, NekDouble> eitimesdHBwd(m_spacedim);
1025
1026 Array<OneD, NekDouble> ntimeseitimesdHFwd(m_spacedim);
1027 Array<OneD, NekDouble> ntimeseitimesdHBwd(m_spacedim);
1028
1029 NekDouble Aver, HFwdk, HBwdk, tmpFwd, tmpBwd;
1030 for (int i = 0; i < nTraceNumPoints; ++i)
1031 {
1032 Aver = 0.5 * (imFwd[i] + imBwd[i]);
1033
1034 // Get [H]
1035 for (int k = 0; k < m_spacedim; k++)
1036 {
1037 HFwdk = Fwd[0][i] * m_MFtraceFwd[0][k][i] +
1038 Fwd[1][i] * m_MFtraceFwd[1][k][i];
1039 HBwdk = Bwd[0][i] * m_MFtraceBwd[0][k][i] +
1040 Bwd[1][i] * m_MFtraceBwd[1][k][i];
1041 dH[k] = HFwdk - HBwdk;
1042
1043 eiFwd[k] = m_MFtraceFwd[dir][k][i];
1044 eiBwd[k] = m_MFtraceBwd[dir][k][i];
1045
1046 nFwd[k] = m_traceNormals[k][i];
1047 }
1048
1049 // MFtraceFwd (MFtraceBwd) \times [H]
1050 // VectorCrossProd(eiFwd, dH, eitimesdHFwd);
1051 // VectorCrossProd(eiBwd, dH, eitimesdHBwd);
1052 VectorCrossProd(nFwd, dH, eitimesdHFwd);
1053 VectorCrossProd(nFwd, dH, eitimesdHBwd);
1054
1055 // n times eitimesdH
1056 VectorCrossProd(nFwd, eitimesdHFwd, ntimeseitimesdHFwd);
1057 VectorCrossProd(nFwd, eitimesdHBwd, ntimeseitimesdHBwd);
1058
1059 // MFtraceFwd \cdot ntimeseitimesdH
1060 tmpFwd = 0.0;
1061 tmpBwd = 0.0;
1062 for (int k = 0; k < m_spacedim; k++)
1063 {
1064 tmpFwd += eiFwd[k] * ntimeseitimesdHFwd[k];
1065 tmpBwd += eiBwd[k] * ntimeseitimesdHBwd[k];
1066 }
1067
1068 outarrayFwd[i] = 0.5 * m_alpha * tmpFwd / Aver;
1069 outarrayBwd[i] = 0.5 * m_alpha * tmpBwd / Aver;
1070 }
1071}

References Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), m_alpha, m_MFtraceBwd, m_MFtraceFwd, Nektar::SolverUtils::EquationSystem::m_spacedim, Nektar::SolverUtils::EquationSystem::m_traceNormals, and VectorCrossProd().

Referenced by NumericalMaxwellFluxTE(), and NumericalMaxwellFluxTM().

◆ ComputeZimYim()

void Nektar::SolverUtils::MMFSystem::ComputeZimYim ( Array< OneD, Array< OneD, NekDouble > > &  epsvec,
Array< OneD, Array< OneD, NekDouble > > &  muvec 
)
protected

Definition at line 1102 of file MMFSystem.cpp.

1105{
1106 int nTraceNumPoints = GetTraceNpoints();
1107
1108 switch (m_TestMaxwellType)
1109 {
1110 case eMaxwell1D:
1111 case eScatField1D:
1112 {
1113 Array<OneD, NekDouble> Fwdeps(nTraceNumPoints, 1.0);
1114 Array<OneD, NekDouble> Bwdeps(nTraceNumPoints, 1.0);
1115 Array<OneD, NekDouble> Fwdmu(nTraceNumPoints, 1.0);
1116 Array<OneD, NekDouble> Bwdmu(nTraceNumPoints, 1.0);
1117 m_fields[0]->GetFwdBwdTracePhys(epsvec[0], Fwdeps, Bwdeps);
1118 m_fields[0]->GetFwdBwdTracePhys(muvec[0], Fwdeps, Bwdeps);
1119
1120 CopyBoundaryTrace(Fwdeps, Bwdeps, eFwdEQBwd, 0);
1121 CopyBoundaryTrace(Fwdmu, Bwdmu, eFwdEQBwd, 1);
1122
1123 m_ZimFwd = Array<OneD, Array<OneD, NekDouble>>(1);
1124 m_ZimBwd = Array<OneD, Array<OneD, NekDouble>>(1);
1125 m_YimFwd = Array<OneD, Array<OneD, NekDouble>>(1);
1126 m_YimBwd = Array<OneD, Array<OneD, NekDouble>>(1);
1127
1128 m_ZimFwd[0] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1129 m_ZimBwd[0] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1130 m_YimFwd[0] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1131 m_YimBwd[0] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1132
1133 // ZimFwd = sqrt( muFwd / epsFwd), ZimBwd = sqrt( muBwd / epsBwd)
1134 for (int i = 0; i < nTraceNumPoints; ++i)
1135 {
1136 m_ZimFwd[0][i] = sqrt(Fwdmu[i] / Fwdeps[i]);
1137 m_ZimBwd[0][i] = sqrt(Bwdmu[i] / Bwdeps[i]);
1138
1139 m_YimFwd[0][i] = 1.0 / m_ZimFwd[0][i];
1140 m_YimBwd[0][i] = 1.0 / m_ZimBwd[0][i];
1141 }
1142
1143 std::cout << "*** ZimFwd = " << RootMeanSquare(m_ZimFwd[0])
1144 << ", ZimBwd = " << RootMeanSquare(m_ZimBwd[0])
1145 << ", YimFwd = " << RootMeanSquare(m_YimFwd[0])
1146 << ", YimBwd = " << RootMeanSquare(m_YimBwd[0])
1147 << std::endl;
1148 }
1149 break;
1150
1151 case eTestMaxwell2DPEC:
1153 case eTestMaxwell2DPMC:
1154 case eScatField2D:
1155 case eTotField2D:
1156 case eMaxwellSphere:
1157 case eELF2DSurface:
1158 {
1159 m_ZimFwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1160 m_ZimBwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1161 m_YimFwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1162 m_YimBwd = Array<OneD, Array<OneD, NekDouble>>(m_shapedim);
1163
1164 for (int j = 0; j < m_shapedim; ++j)
1165 {
1166 m_ZimFwd[j] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1167 m_ZimBwd[j] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1168 m_YimFwd[j] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1169 m_YimBwd[j] = Array<OneD, NekDouble>(nTraceNumPoints, 1.0);
1170 }
1171
1172 switch (m_PolType)
1173 {
1174 case eTransMagnetic:
1175 {
1176 Array<OneD, NekDouble> Fwdmu1(nTraceNumPoints);
1177 Array<OneD, NekDouble> Bwdmu1(nTraceNumPoints);
1178 Array<OneD, NekDouble> Fwdmu2(nTraceNumPoints);
1179 Array<OneD, NekDouble> Bwdmu2(nTraceNumPoints);
1180 Array<OneD, NekDouble> Fwdeps3(nTraceNumPoints);
1181 Array<OneD, NekDouble> Bwdeps3(nTraceNumPoints);
1182
1183 m_fields[0]->GetFwdBwdTracePhys(muvec[0], Fwdmu1, Bwdmu1);
1184 m_fields[0]->GetFwdBwdTracePhys(muvec[1], Fwdmu2, Bwdmu2);
1185 m_fields[0]->GetFwdBwdTracePhys(epsvec[2], Fwdeps3,
1186 Bwdeps3);
1187
1188 CopyBoundaryTrace(Fwdmu1, Bwdmu1, eFwdEQBwd, 0);
1189 CopyBoundaryTrace(Fwdmu2, Bwdmu2, eFwdEQBwd, 1);
1190 CopyBoundaryTrace(Fwdeps3, Bwdeps3, eFwdEQBwd, 2);
1191
1192 // ZimFwd = sqrt( muFwd / epsFwd), ZimBwd = sqrt( muBwd /
1193 // epsBwd)
1194 for (int i = 0; i < nTraceNumPoints; ++i)
1195 {
1196 m_ZimFwd[0][i] = sqrt(Fwdmu2[i] / Fwdeps3[i]);
1197 m_ZimBwd[0][i] = sqrt(Bwdmu2[i] / Bwdeps3[i]);
1198
1199 m_YimFwd[0][i] = 1.0 / m_ZimFwd[0][i];
1200 m_YimBwd[0][i] = 1.0 / m_ZimBwd[0][i];
1201
1202 m_ZimFwd[1][i] = sqrt(Fwdmu1[i] / Fwdeps3[i]);
1203 m_ZimBwd[1][i] = sqrt(Bwdmu1[i] / Bwdeps3[i]);
1204
1205 m_YimFwd[1][i] = 1.0 / m_ZimFwd[1][i];
1206 m_YimBwd[1][i] = 1.0 / m_ZimBwd[1][i];
1207 }
1208 }
1209 break; // eTransMagnetic
1210
1211 case eTransElectric:
1212 {
1213 Array<OneD, NekDouble> Fwdeps1(nTraceNumPoints);
1214 Array<OneD, NekDouble> Bwdeps1(nTraceNumPoints);
1215 Array<OneD, NekDouble> Fwdeps2(nTraceNumPoints);
1216 Array<OneD, NekDouble> Bwdeps2(nTraceNumPoints);
1217 Array<OneD, NekDouble> Fwdmu3(nTraceNumPoints);
1218 Array<OneD, NekDouble> Bwdmu3(nTraceNumPoints);
1219
1220 m_fields[0]->GetFwdBwdTracePhys(epsvec[0], Fwdeps1,
1221 Bwdeps1);
1222 m_fields[0]->GetFwdBwdTracePhys(epsvec[1], Fwdeps2,
1223 Bwdeps2);
1224 m_fields[0]->GetFwdBwdTracePhys(muvec[2], Fwdmu3, Bwdmu3);
1225
1226 CopyBoundaryTrace(Fwdeps1, Bwdeps1, eFwdEQBwd, 0);
1227 CopyBoundaryTrace(Fwdeps2, Bwdeps2, eFwdEQBwd, 1);
1228 CopyBoundaryTrace(Fwdmu3, Bwdmu3, eFwdEQBwd, 2);
1229
1230 for (int i = 0; i < nTraceNumPoints; ++i)
1231 {
1232 m_ZimFwd[0][i] = sqrt(Fwdmu3[i] / Fwdeps2[i]);
1233 m_ZimBwd[0][i] = sqrt(Bwdmu3[i] / Bwdeps2[i]);
1234
1235 m_YimFwd[0][i] = 1.0 / m_ZimFwd[0][i];
1236 m_YimBwd[0][i] = 1.0 / m_ZimBwd[0][i];
1237
1238 m_ZimFwd[1][i] = sqrt(Fwdmu3[i] / Fwdeps1[i]);
1239 m_ZimBwd[1][i] = sqrt(Bwdmu3[i] / Bwdeps1[i]);
1240
1241 m_YimFwd[1][i] = 1.0 / m_ZimFwd[1][i];
1242 m_YimBwd[1][i] = 1.0 / m_ZimBwd[1][i];
1243 }
1244 }
1245 break; // eTransELectric
1246
1247 default:
1248 break;
1249 } // PolType
1250
1251 std::cout << "*** ZimFwd0 = [ "
1252 << Vmath::Vmin(nTraceNumPoints, m_ZimFwd[0], 1) << " , "
1253 << Vmath::Vmax(nTraceNumPoints, m_ZimFwd[0], 1)
1254 << " ], ZimBwd0 = [ "
1255 << Vmath::Vmin(nTraceNumPoints, m_ZimBwd[0], 1) << " , "
1256 << Vmath::Vmax(nTraceNumPoints, m_ZimBwd[0], 1) << " ] "
1257 << std::endl;
1258 std::cout << "*** ZimFwd1 = [ "
1259 << Vmath::Vmin(nTraceNumPoints, m_ZimFwd[1], 1) << " , "
1260 << Vmath::Vmax(nTraceNumPoints, m_ZimFwd[1], 1)
1261 << " ], ZimBwd1 = [ "
1262 << Vmath::Vmin(nTraceNumPoints, m_ZimBwd[1], 1) << " , "
1263 << Vmath::Vmax(nTraceNumPoints, m_ZimBwd[1], 1) << " ] "
1264 << std::endl;
1265 }
1266 break; // eMaxwell2D
1267
1268 default:
1269 break;
1270 } // TestMaxwellType
1271}
Array< OneD, Array< OneD, NekDouble > > m_ZimFwd
Definition: MMFSystem.h:204
Array< OneD, Array< OneD, NekDouble > > m_YimBwd
Definition: MMFSystem.h:207
Array< OneD, Array< OneD, NekDouble > > m_ZimBwd
Definition: MMFSystem.h:205
Array< OneD, Array< OneD, NekDouble > > m_YimFwd
Definition: MMFSystem.h:206
TestMaxwellType m_TestMaxwellType
Definition: MMFSystem.h:153
T Vmin(int n, const T *x, const int incx)
Return the minimum element in x - called vmin to avoid conflict with min.
Definition: Vmath.hpp:725
T Vmax(int n, const T *x, const int incx)
Return the maximum element in x – called vmax to avoid conflict with max.
Definition: Vmath.hpp:644

References CopyBoundaryTrace(), Nektar::SolverUtils::eELF2DSurface, Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eMaxwell1D, Nektar::SolverUtils::eMaxwellSphere, Nektar::SolverUtils::eScatField1D, Nektar::SolverUtils::eScatField2D, Nektar::SolverUtils::eTestMaxwell2DPEC, Nektar::SolverUtils::eTestMaxwell2DPECAVGFLUX, Nektar::SolverUtils::eTestMaxwell2DPMC, Nektar::SolverUtils::eTotField2D, Nektar::SolverUtils::eTransElectric, Nektar::SolverUtils::eTransMagnetic, Nektar::SolverUtils::EquationSystem::GetTraceNpoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_PolType, m_shapedim, m_TestMaxwellType, m_YimBwd, m_YimFwd, m_ZimBwd, m_ZimFwd, RootMeanSquare(), tinysimd::sqrt(), Vmath::Vmax(), and Vmath::Vmin().

Referenced by Nektar::MMFMaxwell::v_InitObject().

◆ CopyBoundaryTrace()

void Nektar::SolverUtils::MMFSystem::CopyBoundaryTrace ( const Array< OneD, const NekDouble > &  Fwd,
Array< OneD, NekDouble > &  Bwd,
const BoundaryCopyType  BDCopyType,
const int  var = 0,
const std::string  btype = "NoUserDefined" 
)

Definition at line 835 of file MMFSystem.cpp.

839{
840 int id1, id2, npts, nptselem, cnt = 0;
841 Array<OneD, NekDouble> Dirichlet, x0, x1, x2;
842
843 // loop over Boundary Regions
844 for (int n = 0; n < m_fields[var]->GetBndConditions().size(); ++n)
845 {
846 nptselem = m_fields[var]->GetBndCondExpansions()[n]->GetNpoints();
847
848 Dirichlet = Array<OneD, NekDouble>(nptselem);
849 x0 = Array<OneD, NekDouble>(nptselem);
850 x1 = Array<OneD, NekDouble>(nptselem);
851 x2 = Array<OneD, NekDouble>(nptselem);
852
853 if (BDCopyType == eDirichlet)
854 {
855 m_fields[var]->GetBndCondExpansions()[n]->GetCoords(x0, x1, x2);
857 m_session->GetFunction("BoundaryConditions", 0);
858 ifunc->Evaluate(x0, x1, x2, 0.0, Dirichlet);
859 }
860
861 for (int e = 0;
862 e < m_fields[var]->GetBndCondExpansions()[n]->GetExpSize(); ++e)
863 {
864 npts = m_fields[var]
865 ->GetBndCondExpansions()[n]
866 ->GetExp(e)
867 ->GetNumPoints(0);
868 id1 = m_fields[var]->GetBndCondExpansions()[n]->GetPhys_Offset(e);
869 id2 = m_fields[var]->GetTrace()->GetPhys_Offset(
870 m_fields[var]->GetTraceMap()->GetBndCondIDToGlobalTraceID(cnt +
871 e));
872
873 if (m_fields[var]->GetBndConditions()[n]->GetUserDefined() ==
874 BDtype ||
875 BDtype == "NoUserDefined")
876 {
877 switch (BDCopyType)
878 {
879 case eDirichlet:
880 {
881 Vmath::Vcopy(npts, &Dirichlet[id1], 1, &Bwd[id2], 1);
882 }
883 break;
884
885 case eFwdEQBwd:
886 {
887 Vmath::Vcopy(npts, &Fwd[id2], 1, &Bwd[id2], 1);
888 }
889 break;
890
891 case eFwdEQNegBwd:
892 {
893 Vmath::Vcopy(npts, &Fwd[id2], 1, &Bwd[id2], 1);
894 Vmath::Neg(npts, &Bwd[id2], 1);
895 }
896 break;
897
898 default:
899 break;
900 }
901 }
902 }
903
904 cnt += m_fields[var]->GetBndCondExpansions()[n]->GetExpSize();
905 }
906}
SOLVER_UTILS_EXPORT int GetExpSize()
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
std::shared_ptr< Equation > EquationSharedPtr
Definition: Equation.h:125

References Nektar::SolverUtils::eDirichlet, Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, Nektar::SolverUtils::EquationSystem::GetExpSize(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_session, Vmath::Neg(), and Vmath::Vcopy().

Referenced by AverageMaxwellFlux1D(), ComputeMFtrace(), ComputencdotMF(), ComputeZimYim(), LaxFriedrichMaxwellFlux1D(), NumericalMaxwellFluxTE(), NumericalMaxwellFluxTM(), Nektar::MMFSWE::NumericalSWEFlux(), and UpwindMaxwellFlux1D().

◆ DeriveCrossProductMF()

void Nektar::SolverUtils::MMFSystem::DeriveCrossProductMF ( Array< OneD, Array< OneD, NekDouble > > &  CrossProductMF)
protected

Definition at line 568 of file MMFSystem.cpp.

570{
571 int MFdim = 3;
572 int nq = GetTotPoints();
573
574 CrossProductMF = Array<OneD, Array<OneD, NekDouble>>(MFdim);
575
576 Array<OneD, Array<OneD, NekDouble>> MF1tmp(m_spacedim);
577 Array<OneD, Array<OneD, NekDouble>> MF2tmp(m_spacedim);
578 Array<OneD, Array<OneD, NekDouble>> MF3tmp(m_spacedim);
579 Array<OneD, Array<OneD, Array<OneD, NekDouble>>> MFtmpCurl(MFdim);
580 for (int j = 0; j < MFdim; ++j)
581 {
582 CrossProductMF[j] = Array<OneD, NekDouble>(nq * m_spacedim);
583 MFtmpCurl[j] = Array<OneD, Array<OneD, NekDouble>>(m_spacedim);
584 for (int k = 0; k < m_spacedim; ++k)
585 {
586 MFtmpCurl[j][k] = Array<OneD, NekDouble>(nq);
587 }
588 }
589
590 for (int k = 0; k < m_spacedim; ++k)
591 {
592 MF1tmp[k] = Array<OneD, NekDouble>(nq);
593 MF2tmp[k] = Array<OneD, NekDouble>(nq);
594 MF3tmp[k] = Array<OneD, NekDouble>(nq);
595
596 Vmath::Vcopy(nq, &m_movingframes[0][k * nq], 1, &MF1tmp[k][0], 1);
597 Vmath::Vcopy(nq, &m_movingframes[1][k * nq], 1, &MF2tmp[k][0], 1);
598 Vmath::Vcopy(nq, &m_movingframes[2][k * nq], 1, &MF3tmp[k][0], 1);
599 }
600
601 VectorCrossProd(MF3tmp, MF1tmp, MFtmpCurl[0]);
602 VectorCrossProd(MF2tmp, MF3tmp, MFtmpCurl[1]);
603 VectorCrossProd(MF1tmp, MF2tmp, MFtmpCurl[2]);
604
605 for (int j = 0; j < MFdim; ++j)
606 {
607 for (int k = 0; k < m_spacedim; ++k)
608 {
609 Vmath::Vcopy(nq, &MFtmpCurl[j][k][0], 1, &CrossProductMF[j][k * nq],
610 1);
611 }
612 }
613}

References Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_movingframes, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vcopy(), and VectorCrossProd().

Referenced by Nektar::MMFMaxwell::v_InitObject().

◆ GetIncidentField()

Array< OneD, NekDouble > Nektar::SolverUtils::MMFSystem::GetIncidentField ( const int  var,
const NekDouble  time 
)
protected

Definition at line 2005 of file MMFSystem.cpp.

2007{
2008 int nq = m_fields[0]->GetNpoints();
2009
2010 Array<OneD, NekDouble> x0(nq);
2011 Array<OneD, NekDouble> x1(nq);
2012 Array<OneD, NekDouble> x2(nq);
2013
2014 m_fields[0]->GetCoords(x0, x1, x2);
2015
2016 // GetSmoothFactor such that wave propages from the left to the object.
2017 // a = 0.1, ta = 1, f = 1.0./(1.0 + exp( -0.5.*(time-ta)/a ));
2018 Array<OneD, NekDouble> SmoothFactor(nq, 1.0);
2019 switch (m_SmoothFactor)
2020 {
2021 case 0:
2022 {
2023 for (int i = 0; i < nq; i++)
2024 {
2025 SmoothFactor[i] = 1.0 / (1.0 + exp(-1.0 * (time - 1.0) / 0.1));
2026 }
2027 }
2028 break;
2029
2030 case 1:
2031 {
2032 NekDouble xmin = Vmath::Vmin(nq, x0, 1);
2033 NekDouble xp;
2034 for (int i = 0; i < nq; i++)
2035 {
2036 xp = x0[i] - xmin - time - m_SFinit;
2037 if (xp > 0.0)
2038 {
2039 SmoothFactor[i] =
2040 2.0 / (1.0 + exp(0.5 * (sqrt(xp * xp) - 0.1)));
2041 }
2042
2043 else
2044 {
2045 SmoothFactor[i] = 1.0;
2046 }
2047 }
2048 }
2049 break;
2050
2051 default:
2052 break;
2053 }
2054
2055 // Generate a factor for smoothly increasing wave
2056 Array<OneD, NekDouble> F1(nq);
2057 Array<OneD, NekDouble> F2(nq);
2058 Array<OneD, NekDouble> F3(nq);
2059
2060 Array<OneD, NekDouble> dF1dt(nq);
2061 Array<OneD, NekDouble> dF2dt(nq);
2062 Array<OneD, NekDouble> dF3dt(nq);
2063
2064 Array<OneD, NekDouble> F1int(nq);
2065 Array<OneD, NekDouble> F2int(nq);
2066 Array<OneD, NekDouble> F3int(nq);
2067
2068 Array<OneD, NekDouble> outarray(nq);
2069
2070 switch (m_IncType)
2071 {
2072 case ePlaneWave:
2073 {
2074 NekDouble cs, sn;
2075 NekDouble e1y, e2y;
2076 switch (m_PolType)
2077 {
2078 case eTransMagnetic:
2079 {
2080 // H = 0 \hat{x} - e^{ikx} \hat{y}
2081 // E = e^{ikx} \hat{z}
2082 // outarray1 = Hr1inc
2083 // outarray2 = Hr2inc
2084 // outarray3 = Ezrinc
2085 for (int i = 0; i < nq; i++)
2086 {
2087 e1y = m_movingframes[0][nq + i];
2088 e2y = m_movingframes[1][nq + i];
2089
2090 cs = SmoothFactor[i] * cos(m_Incfreq * (x0[i] - time));
2091 sn = SmoothFactor[i] * sin(m_Incfreq * (x0[i] - time));
2092
2093 F1[i] = -1.0 * cs * e1y;
2094 F2[i] = -1.0 * cs * e2y;
2095 F3[i] = cs;
2096
2097 dF1dt[i] = -1.0 * m_Incfreq * sn * e1y;
2098 dF2dt[i] = -1.0 * m_Incfreq * sn * e2y;
2099 dF3dt[i] = 1.0 * m_Incfreq * sn;
2100
2101 F1int[i] = (1.0 / m_Incfreq) * sn * e1y;
2102 F2int[i] = (1.0 / m_Incfreq) * sn * e2y;
2103 F3int[i] = (-1.0 / m_Incfreq) * sn;
2104 }
2105 }
2106 break;
2107
2108 case eTransElectric:
2109 {
2110 // E = 0 \hat{x} + e^{ikx} \hat{y}
2111 // H = e^{ikx} \hat{z}
2112 // outarray1 = Er1inc
2113 // outarray2 = Er2inc
2114 // outarray3 = Hzrinc
2115 // outarray1 = Ei1inc
2116 // outarray2 = Ei2inc
2117 // outarray3 = Hziinc
2118 for (int i = 0; i < nq; i++)
2119 {
2120 e1y = m_movingframes[0][nq + i];
2121 e2y = m_movingframes[1][nq + i];
2122
2123 cs = SmoothFactor[i] * cos(m_Incfreq * (x0[i] - time));
2124 sn = SmoothFactor[i] * sin(m_Incfreq * (x0[i] - time));
2125
2126 F1[i] = cs * e1y;
2127 F2[i] = cs * e2y;
2128 F3[i] = cs;
2129
2130 dF1dt[i] = m_Incfreq * sn * e1y;
2131 dF2dt[i] = m_Incfreq * sn * e2y;
2132 dF3dt[i] = m_Incfreq * sn;
2133
2134 F1int[i] = (-1.0 / m_Incfreq) * sn * e1y;
2135 F2int[i] = (-1.0 / m_Incfreq) * sn * e2y;
2136 F3int[i] = (-1.0 / m_Incfreq) * sn;
2137 }
2138 }
2139 break;
2140
2141 default:
2142 break;
2143 }
2144 }
2145 break;
2146
2147 case ePlaneWaveImag:
2148 {
2149 NekDouble cs, sn;
2150 NekDouble e1y, e2y;
2151 switch (m_PolType)
2152 {
2153 case eTransMagnetic:
2154 {
2155 // H = 0 \hat{x} - e^{ikx} \hat{y}
2156 // E = e^{ikx} \hat{z}
2157 // outarray1 = Hr1inc
2158 // outarray2 = Hr2inc
2159 // outarray3 = Ezrinc
2160 for (int i = 0; i < nq; i++)
2161 {
2162 e1y = m_movingframes[0][nq + i];
2163 e2y = m_movingframes[1][nq + i];
2164
2165 cs = SmoothFactor[i] * cos(m_Incfreq * (x0[i] - time));
2166 sn = SmoothFactor[i] * sin(m_Incfreq * (x0[i] - time));
2167
2168 F1[i] = -1.0 * sn * e1y;
2169 F2[i] = -1.0 * sn * e2y;
2170 F3[i] = sn;
2171
2172 dF1dt[i] = m_Incfreq * cs * e1y;
2173 dF2dt[i] = m_Incfreq * cs * e2y;
2174 dF3dt[i] = -1.0 * m_Incfreq * cs;
2175
2176 F1int[i] = (-1.0 / m_Incfreq) * cs * e1y;
2177 F2int[i] = (-1.0 / m_Incfreq) * cs * e2y;
2178 F3int[i] = (1.0 / m_Incfreq) * cs;
2179 }
2180 }
2181 break;
2182
2183 case eTransElectric:
2184 {
2185 // E = 0 \hat{x} + e^{ikx} \hat{y}
2186 // H = e^{ikx} \hat{z}
2187 // outarray1 = Er1inc
2188 // outarray2 = Er2inc
2189 // outarray3 = Hzrinc
2190 // outarray1 = Ei1inc
2191 // outarray2 = Ei2inc
2192 // outarray3 = Hziinc
2193 for (int i = 0; i < nq; i++)
2194 {
2195 e1y = m_movingframes[0][nq + i];
2196 e2y = m_movingframes[1][nq + i];
2197
2198 cs = SmoothFactor[i] * cos(m_Incfreq * (x0[i] - time));
2199 sn = SmoothFactor[i] * sin(m_Incfreq * (x0[i] - time));
2200
2201 F1[i] = sn * e1y;
2202 F2[i] = sn * e2y;
2203 F3[i] = sn;
2204
2205 dF1dt[i] = -1.0 * m_Incfreq * cs * e1y;
2206 dF2dt[i] = -1.0 * m_Incfreq * cs * e2y;
2207 dF3dt[i] = -1.0 * m_Incfreq * cs;
2208
2209 F1int[i] = (1.0 / m_Incfreq) * cs * e1y;
2210 F2int[i] = (1.0 / m_Incfreq) * cs * e2y;
2211 F3int[i] = (1.0 / m_Incfreq) * cs;
2212 }
2213 }
2214 break;
2215
2216 default:
2217 break;
2218 }
2219 }
2220 break;
2221
2222 default:
2223 break;
2224 }
2225
2226 switch (var)
2227 {
2228 case 0:
2229 {
2230 outarray = F1;
2231 }
2232 break;
2233
2234 case 1:
2235 {
2236 outarray = F2;
2237 }
2238 break;
2239
2240 case 2:
2241 {
2242 outarray = F3;
2243 }
2244 break;
2245
2246 case 10:
2247 {
2248 outarray = dF1dt;
2249 }
2250 break;
2251
2252 case 11:
2253 {
2254 outarray = dF2dt;
2255 }
2256 break;
2257
2258 case 12:
2259 {
2260 outarray = dF3dt;
2261 }
2262 break;
2263
2264 case 20:
2265 {
2266 outarray = F1int;
2267 }
2268 break;
2269
2270 case 21:
2271 {
2272 outarray = F2int;
2273 }
2274 break;
2275
2276 case 22:
2277 {
2278 outarray = F3int;
2279 }
2280 break;
2281
2282 default:
2283 {
2284 Vmath::Zero(nq, outarray, 1);
2285 }
2286 break;
2287 }
2288
2289 return outarray;
2290}
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.hpp:273

References Nektar::SolverUtils::ePlaneWave, Nektar::SolverUtils::ePlaneWaveImag, Nektar::SolverUtils::eTransElectric, Nektar::SolverUtils::eTransMagnetic, Nektar::SolverUtils::EquationSystem::m_fields, m_Incfreq, m_IncType, m_movingframes, m_PolType, m_SFinit, m_SmoothFactor, tinysimd::sqrt(), Vmath::Vmin(), and Vmath::Zero().

Referenced by Nektar::MMFMaxwell::Checkpoint_TotalFieldOutput(), Nektar::MMFMaxwell::Checkpoint_TotPlotOutput(), Nektar::MMFMaxwell::ComputeSurfaceCurrent(), Nektar::MMFMaxwell::DoOdeRhs(), NumericalMaxwellFluxTE(), NumericalMaxwellFluxTM(), and Nektar::MMFMaxwell::v_DoSolve().

◆ GetMaxwellFlux1D()

void Nektar::SolverUtils::MMFSystem::GetMaxwellFlux1D ( const int  var,
const Array< OneD, const Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  flux 
)
protected

Definition at line 1633 of file MMFSystem.cpp.

1636{
1637 int nq = m_fields[0]->GetTotPoints();
1638
1639 switch (var)
1640 {
1641 case 0:
1642 {
1643 // H in flux 0
1644 Vmath::Vcopy(nq, physfield[1], 1, flux[0], 1);
1645
1646 // E in flux 1
1647 Vmath::Zero(nq, flux[1], 1);
1648 }
1649 break;
1650
1651 case 1:
1652 {
1653 // E in flux 0
1654 Vmath::Vcopy(nq, physfield[0], 1, flux[0], 1);
1655
1656 // H in flux 1
1657 Vmath::Zero(nq, flux[1], 1);
1658 }
1659 break;
1660 //----------------------------------------------------
1661
1662 default:
1663 break;
1664 }
1665}

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

Referenced by GetMaxwellFluxVector().

◆ GetMaxwellFlux2D()

void Nektar::SolverUtils::MMFSystem::GetMaxwellFlux2D ( const int  var,
const Array< OneD, const Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  flux 
)
protected

Definition at line 1667 of file MMFSystem.cpp.

1670{
1671 int nq = m_fields[0]->GetTotPoints();
1672
1673 NekDouble sign = 1.0;
1674 switch (m_PolType)
1675 {
1676 // TransMagnetic
1677 case 0:
1678 {
1679 sign = -1.0;
1680 }
1681 break;
1682
1683 // TransElectric
1684 case 1:
1685 {
1686 sign = 1.0;
1687 }
1688 break;
1689
1690 default:
1691 break;
1692 }
1693
1694 switch (var)
1695 {
1696 case 0:
1697 {
1698 // -Ez in flux 1
1699 Vmath::Smul(nq, sign, physfield[2], 1, flux[0], 1);
1700 Vmath::Zero(nq, flux[1], 1);
1701 }
1702 break;
1703
1704 case 1:
1705 {
1706 // Ez in flux 0
1707 Vmath::Zero(nq, flux[0], 1);
1708 Vmath::Smul(nq, -sign, physfield[2], 1, flux[1], 1);
1709 }
1710 break;
1711
1712 case 2:
1713 {
1714 Vmath::Smul(nq, sign, physfield[0], 1, flux[0], 1);
1715 Vmath::Smul(nq, -sign, physfield[1], 1, flux[1], 1);
1716 }
1717 break;
1718
1719 default:
1720 ASSERTL0(false, "GetFluxVector2D: illegal vector index");
1721 }
1722}
#define sign(a, b)
return the sign(b)*a
Definition: Polylib.cpp:47
void Smul(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Scalar multiply y = alpha*x.
Definition: Vmath.hpp:100

References ASSERTL0, Nektar::SolverUtils::EquationSystem::m_fields, m_PolType, sign, Vmath::Smul(), and Vmath::Zero().

Referenced by GetMaxwellFluxVector().

◆ GetMaxwellFluxVector()

void Nektar::SolverUtils::MMFSystem::GetMaxwellFluxVector ( const int  var,
const Array< OneD, const Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  flux 
)
protected

Definition at line 1603 of file MMFSystem.cpp.

1606{
1607 switch (m_TestMaxwellType)
1608 {
1609 case eMaxwell1D:
1610 case eScatField1D:
1611 {
1612 GetMaxwellFlux1D(var, physfield, flux);
1613 }
1614 break;
1615
1616 case eTestMaxwell2DPEC:
1618 case eTestMaxwell2DPMC:
1619 case eScatField2D:
1620 case eTotField2D:
1621 case eMaxwellSphere:
1622 case eELF2DSurface:
1623 {
1624 GetMaxwellFlux2D(var, physfield, flux);
1625 }
1626 break;
1627
1628 default:
1629 break;
1630 }
1631}
SOLVER_UTILS_EXPORT void GetMaxwellFlux1D(const int var, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
Definition: MMFSystem.cpp:1633
SOLVER_UTILS_EXPORT void GetMaxwellFlux2D(const int var, const Array< OneD, const Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &flux)
Definition: MMFSystem.cpp:1667

References Nektar::SolverUtils::eELF2DSurface, Nektar::SolverUtils::eMaxwell1D, Nektar::SolverUtils::eMaxwellSphere, Nektar::SolverUtils::eScatField1D, Nektar::SolverUtils::eScatField2D, Nektar::SolverUtils::eTestMaxwell2DPEC, Nektar::SolverUtils::eTestMaxwell2DPECAVGFLUX, Nektar::SolverUtils::eTestMaxwell2DPMC, Nektar::SolverUtils::eTotField2D, GetMaxwellFlux1D(), GetMaxwellFlux2D(), and m_TestMaxwellType.

Referenced by Nektar::MMFMaxwell::AddGreenDerivCompensate(), and Nektar::MMFMaxwell::WeakDGMaxwellDirDeriv().

◆ GramSchumitz()

void Nektar::SolverUtils::MMFSystem::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 
)
protected

Definition at line 2399 of file MMFSystem.cpp.

2403{
2404
2405 int nq = v1[0].size();
2406 Array<OneD, NekDouble> tmp(nq, 0.0);
2407 Array<OneD, NekDouble> mag(nq, 0.0);
2408
2409 for (int i = 0; i < m_spacedim; ++i)
2410 {
2411 // u2 = v2 - < u1 , v2 > ( u1 / < u1, u1 > )
2412 Vmath::Vvtvp(nq, &v1[i][0], 1, &v2[i][0], 1, &tmp[0], 1, &tmp[0], 1);
2413 Vmath::Vvtvp(nq, &v1[i][0], 1, &v1[i][0], 1, &mag[0], 1, &mag[0], 1);
2414 }
2415 Vmath::Vdiv(nq, &tmp[0], 1, &mag[0], 1, &tmp[0], 1);
2416 Vmath::Neg(nq, &tmp[0], 1);
2417
2418 // outarray = Array<OneD, Array<OneD, NekDouble> >(m_spacedim);
2419
2420 // u2 = v2 - < u1 , v2 > ( u1 / < u1, u1 > )
2421 for (int i = 0; i < m_spacedim; ++i)
2422 {
2423 outarray[i] = Array<OneD, NekDouble>(nq, 0.0);
2424 Vmath::Vvtvp(nq, &tmp[0], 1, &v1[i][0], 1, &v2[i][0], 1,
2425 &outarray[i][0], 1);
2426 }
2427
2428 if (KeepTheMagnitude)
2429 {
2430 Array<OneD, NekDouble> magorig(nq, 0.0);
2431 Array<OneD, NekDouble> magnew(nq, 0.0);
2432
2433 for (int i = 0; i < m_spacedim; ++i)
2434 {
2435 Vmath::Vmul(nq, &v2[0][0], 1, &v2[0][0], 1, &magorig[0], 1);
2436 Vmath::Vvtvp(nq, &v2[1][0], 1, &v2[1][0], 1, &magorig[0], 1,
2437 &magorig[0], 1);
2438 Vmath::Vvtvp(nq, &v2[2][0], 1, &v2[2][0], 1, &magorig[0], 1,
2439 &magorig[0], 1);
2440
2441 Vmath::Vmul(nq, &outarray[0][0], 1, &outarray[0][0], 1, &magnew[0],
2442 1);
2443 Vmath::Vvtvp(nq, &outarray[1][0], 1, &outarray[1][0], 1, &magnew[0],
2444 1, &magnew[0], 1);
2445 Vmath::Vvtvp(nq, &outarray[2][0], 1, &outarray[2][0], 1, &magnew[0],
2446 1, &magnew[0], 1);
2447 }
2448
2449 for (int i = 0; i < m_spacedim; ++i)
2450 {
2451 for (int j = 0; j < nq; ++j)
2452 {
2453 if (fabs(magnew[j]) > 0.000000001)
2454 {
2455 outarray[i][j] =
2456 outarray[i][j] * sqrt(magorig[j] / magnew[j]);
2457 }
2458 }
2459 }
2460 }
2461}
void Vdiv(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:126

References Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Neg(), tinysimd::sqrt(), Vmath::Vdiv(), Vmath::Vmul(), and Vmath::Vvtvp().

Referenced by Nektar::MMFAdvection::EvaluateAdvectionVelocity().

◆ LaxFriedrichMaxwellFlux1D()

void Nektar::SolverUtils::MMFSystem::LaxFriedrichMaxwellFlux1D ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd 
)
protected

Definition at line 1509 of file MMFSystem.cpp.

1513{
1514 int i;
1515 int nTraceNumPoints = GetTraceTotPoints();
1516 int nvar = 2;
1517
1518 // get temporary arrays
1519 Array<OneD, Array<OneD, NekDouble>> Fwd(nvar);
1520 Array<OneD, Array<OneD, NekDouble>> Bwd(nvar);
1521
1522 for (i = 0; i < nvar; ++i)
1523 {
1524 Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1525 Bwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1526 }
1527
1528 // get the physical values at the trace from the dependent variables
1529 for (i = 0; i < nvar; ++i)
1530 {
1531 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd[i], Bwd[i]);
1532 }
1533
1534 // E = 0 at the boundaries
1535 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0, "PEC");
1536
1537 // d H / d n = 0 at the boundaries
1538 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PEC");
1539
1540 Array<OneD, NekDouble> dE(nTraceNumPoints);
1541 Array<OneD, NekDouble> dH(nTraceNumPoints);
1542
1543 Vmath::Vsub(nTraceNumPoints, &Fwd[0][0], 1, &Bwd[0][0], 1, &dE[0], 1);
1544 Vmath::Vsub(nTraceNumPoints, &Fwd[1][0], 1, &Bwd[1][0], 1, &dH[0], 1);
1545
1546 NekDouble nx;
1547 for (i = 0; i < nTraceNumPoints; ++i)
1548 {
1549 nx = m_traceNormals[0][i];
1550 numfluxFwd[0][i] =
1551 0.5 * nx * ((Fwd[1][i] + Bwd[1][i]) - nx * m_YimFwd[0][i] * dE[i]);
1552 numfluxFwd[1][i] =
1553 0.5 * nx * ((Fwd[0][i] + Bwd[0][i]) - nx * m_ZimFwd[0][i] * dH[i]);
1554
1555 numfluxBwd[0][i] =
1556 0.5 * nx * ((Fwd[1][i] + Bwd[1][i]) - nx * m_YimFwd[0][i] * dE[i]);
1557 numfluxBwd[1][i] =
1558 0.5 * nx * ((Fwd[0][i] + Bwd[0][i]) - nx * m_ZimFwd[0][i] * dH[i]);
1559 }
1560}

References CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_traceNormals, m_YimFwd, m_ZimFwd, and Vmath::Vsub().

Referenced by NumericalMaxwellFlux().

◆ MMFInitObject()

void Nektar::SolverUtils::MMFSystem::MMFInitObject ( const Array< OneD, const Array< OneD, NekDouble > > &  Anisotropy,
const int  TangentXelem = -1 
)

Definition at line 51 of file MMFSystem.cpp.

54{
55 m_pi = 3.14159265358979323846;
56 m_shapedim = m_fields[0]->GetShapeDimension();
57
58 ASSERTL0(m_spacedim == 3, "Space Dimension should be 3.");
59
60 // Factor for Numerical Flux
61 m_session->LoadParameter("alpha", m_alpha, 1.0);
62
63 // Factor for Numerical Flux
64 m_session->LoadParameter("Incfreq", m_Incfreq, 1.0);
65
66 // SmoothFactor
67 m_session->LoadParameter("SmoothFactor", m_SmoothFactor, 1);
68 m_session->LoadParameter("SFinit", m_SFinit, 0.0);
69
70 // Define SurfaceType
71 if (m_session->DefinesSolverInfo("SURFACETYPE"))
72 {
73 std::string SurfaceTypeStr = m_session->GetSolverInfo("SURFACETYPE");
74 int i;
75 for (i = 0; i < (int)SIZE_SurfaceType; ++i)
76 {
77 if (boost::iequals(SurfaceTypeMap[i], SurfaceTypeStr))
78 {
80 break;
81 }
82 }
83 }
84 else
85 {
87 }
88
89 // if discontinuous Galerkin determine numerical flux to use
90 for (int i = 0; i < (int)SIZE_UpwindType; ++i)
91 {
92 bool match;
93 m_session->MatchSolverInfo("UPWINDTYPE", UpwindTypeMap[i], match,
94 false);
95 if (match)
96 {
98 break;
99 }
100 }
101
102 // SetUpMovingFrames: To generate m_movingframes
103 SetUpMovingFrames(Anisotropy, TangentXelem);
104
105 // Derive movingfrmaes at interfaces
106 // Get: m_MFtraceFwd and m_MFtraceBwd
108
109 // Check Movingframes and surfraceNormal
110 // Get: m_ncdotMFFwd,m_ncdotMFBwd,m_nperpcdotMFFwd,m_nperpcdotMFBwd
112
113 // Compute nabla cdot m_movingframes
114 // Get: m_DivMF, m_CurlMF
116}
void SetUpMovingFrames(const Array< OneD, const Array< OneD, NekDouble > > &Anisotropy, const int TangentXelem)
Definition: MMFSystem.cpp:118
SOLVER_UTILS_EXPORT void ComputeDivCurlMF()
Definition: MMFSystem.cpp:442
SOLVER_UTILS_EXPORT void ComputeMFtrace()
Definition: MMFSystem.cpp:519
SOLVER_UTILS_EXPORT void ComputencdotMF()
Definition: MMFSystem.cpp:317
const char *const UpwindTypeMap[]
Definition: MMFSystem.h:86
@ SIZE_UpwindType
Length of enum list.
Definition: MMFSystem.h:83
const char *const SurfaceTypeMap[]
Definition: MMFSystem.h:54

References ASSERTL0, ComputeDivCurlMF(), ComputeMFtrace(), ComputencdotMF(), m_alpha, Nektar::SolverUtils::EquationSystem::m_fields, m_Incfreq, m_pi, Nektar::SolverUtils::EquationSystem::m_session, m_SFinit, m_shapedim, m_SmoothFactor, Nektar::SolverUtils::EquationSystem::m_spacedim, m_surfaceType, m_upwindType, SetUpMovingFrames(), Nektar::SolverUtils::SIZE_SurfaceType, Nektar::SolverUtils::SIZE_UpwindType, Nektar::SolverUtils::SurfaceTypeMap, and Nektar::SolverUtils::UpwindTypeMap.

Referenced by Nektar::MMFDiffusion::v_InitObject(), Nektar::MMFAdvection::v_InitObject(), Nektar::MMFMaxwell::v_InitObject(), and Nektar::MMFSWE::v_InitObject().

◆ NumericalMaxwellFlux()

void Nektar::SolverUtils::MMFSystem::NumericalMaxwellFlux ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd,
const NekDouble  time = 0.0 
)
protected

Definition at line 1724 of file MMFSystem.cpp.

1728{
1729
1730 switch (m_TestMaxwellType)
1731 {
1732 case eMaxwell1D:
1733 case eScatField1D:
1734 {
1735 switch (m_upwindType)
1736 {
1737 case eAverage:
1738 {
1739 AverageMaxwellFlux1D(physfield, numfluxFwd, numfluxBwd);
1740 }
1741 break;
1742
1743 case eLaxFriedrich:
1744 {
1745 LaxFriedrichMaxwellFlux1D(physfield, numfluxFwd,
1746 numfluxBwd);
1747 }
1748 break;
1749
1750 case eUpwind:
1751 {
1752 UpwindMaxwellFlux1D(physfield, numfluxFwd, numfluxBwd);
1753 }
1754 break;
1755
1756 default:
1757 {
1758 ASSERTL0(false,
1759 "populate switch statement for upwind flux");
1760 }
1761 break; // upwindType
1762 }
1763 }
1764 break; // eMaxwell1D
1765
1766 case eTestMaxwell2DPEC:
1768 case eTestMaxwell2DPMC:
1769 case eScatField2D:
1770 case eTotField2D:
1771 case eMaxwellSphere:
1772 case eELF2DSurface:
1773 {
1774 switch (m_PolType)
1775 {
1776 case eTransMagnetic:
1777 {
1778 NumericalMaxwellFluxTM(physfield, numfluxFwd, numfluxBwd,
1779 time);
1780 }
1781 break;
1782
1783 case eTransElectric:
1784 {
1785 NumericalMaxwellFluxTE(physfield, numfluxFwd, numfluxBwd,
1786 time);
1787 }
1788 break;
1789
1790 default:
1791 break;
1792 }
1793 }
1794 break; // eMaxwell2D
1795
1796 default:
1797 break;
1798 } // m_TestMaxwellType
1799}
SOLVER_UTILS_EXPORT void LaxFriedrichMaxwellFlux1D(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numfluxFwd, Array< OneD, Array< OneD, NekDouble > > &numfluxBwd)
Definition: MMFSystem.cpp:1509
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)
Definition: MMFSystem.cpp:1801
SOLVER_UTILS_EXPORT void AverageMaxwellFlux1D(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numfluxFwd, Array< OneD, Array< OneD, NekDouble > > &numfluxBwd)
Definition: MMFSystem.cpp:1562
SOLVER_UTILS_EXPORT void UpwindMaxwellFlux1D(Array< OneD, Array< OneD, NekDouble > > &physfield, Array< OneD, Array< OneD, NekDouble > > &numfluxFwd, Array< OneD, Array< OneD, NekDouble > > &numfluxBwd)
Definition: MMFSystem.cpp:1455
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)
Definition: MMFSystem.cpp:1902
@ eAverage
averaged (or centred) flux
Definition: MMFSystem.h:77
@ eLaxFriedrich
Lax-Friedrich flux.
Definition: MMFSystem.h:78

References ASSERTL0, AverageMaxwellFlux1D(), Nektar::SolverUtils::eAverage, Nektar::SolverUtils::eELF2DSurface, Nektar::SolverUtils::eLaxFriedrich, Nektar::SolverUtils::eMaxwell1D, Nektar::SolverUtils::eMaxwellSphere, Nektar::SolverUtils::eScatField1D, Nektar::SolverUtils::eScatField2D, Nektar::SolverUtils::eTestMaxwell2DPEC, Nektar::SolverUtils::eTestMaxwell2DPECAVGFLUX, Nektar::SolverUtils::eTestMaxwell2DPMC, Nektar::SolverUtils::eTotField2D, Nektar::SolverUtils::eTransElectric, Nektar::SolverUtils::eTransMagnetic, Nektar::SolverUtils::eUpwind, LaxFriedrichMaxwellFlux1D(), m_PolType, m_TestMaxwellType, m_upwindType, NumericalMaxwellFluxTE(), NumericalMaxwellFluxTM(), and UpwindMaxwellFlux1D().

Referenced by Nektar::MMFMaxwell::WeakDGMaxwellDirDeriv().

◆ NumericalMaxwellFluxTE()

void Nektar::SolverUtils::MMFSystem::NumericalMaxwellFluxTE ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd,
const NekDouble  time 
)
protected

Definition at line 1902 of file MMFSystem.cpp.

1907{
1908 int nq = m_fields[0]->GetNpoints();
1909 int nTraceNumPoints = GetTraceTotPoints();
1910 int nvar = physfield.size();
1911
1912 // Get temporary arrays
1913 Array<OneD, Array<OneD, NekDouble>> Fwd(nvar);
1914 Array<OneD, Array<OneD, NekDouble>> Bwd(nvar);
1915 for (int i = 0; i < nvar; ++i)
1916 {
1917 Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1918 Bwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1919
1920 // get the physical values at the trace
1921 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd[i], Bwd[i]);
1922 }
1923
1924 // E = 0 at the PEC boundaries:
1925 Array<OneD, NekDouble> IncField(nq, 0.0);
1926 Array<OneD, NekDouble> IncFieldBwd(nTraceNumPoints, 0.0);
1927 Array<OneD, Array<OneD, NekDouble>> IncFieldFwd(m_spacedim);
1928 for (int i = 0; i < m_spacedim; ++i)
1929 {
1930 IncFieldFwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1931
1932 IncField = GetIncidentField(i, time);
1933 m_fields[0]->GetFwdBwdTracePhys(IncField, IncFieldFwd[i], IncFieldBwd);
1934
1935 Vmath::Svtvp(nTraceNumPoints, 2.0, &IncFieldFwd[i][0], 1, &Fwd[i][0], 1,
1936 &IncFieldFwd[i][0], 1);
1937 }
1938
1939 CopyBoundaryTrace(IncFieldFwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0,
1940 "PEC_Forces");
1941 CopyBoundaryTrace(IncFieldFwd[1], Bwd[1], SolverUtils::eFwdEQNegBwd, 1,
1942 "PEC_Forces");
1943 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQBwd, 2, "PEC_Forces");
1944
1945 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQBwd, 0, "PMC_Forces");
1946 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PMC_Forces");
1947 CopyBoundaryTrace(IncFieldFwd[2], Bwd[2], SolverUtils::eFwdEQNegBwd, 2,
1948 "PMC_Forces");
1949
1950 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0, "PEC");
1951 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQNegBwd, 1, "PEC");
1952 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQBwd, 2, "PEC");
1953
1954 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQBwd, 0, "PMC");
1955 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PMC");
1956 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQNegBwd, 2, "PMC");
1957
1958 Array<OneD, NekDouble> e1Fwd_cdot_ncrossdE(nTraceNumPoints);
1959 Array<OneD, NekDouble> e1Bwd_cdot_ncrossdE(nTraceNumPoints);
1960 Array<OneD, NekDouble> e2Fwd_cdot_ncrossdE(nTraceNumPoints);
1961 Array<OneD, NekDouble> e2Bwd_cdot_ncrossdE(nTraceNumPoints);
1962 Array<OneD, NekDouble> e3Fwd_cdot_dHe3(nTraceNumPoints);
1963 Array<OneD, NekDouble> e3Bwd_cdot_dHe3(nTraceNumPoints);
1964
1965 // Compute numfluxFwd[dir] = (eFwd^[dir] \cdot n \times e^3) * (ZimFwd *
1966 // HFwd^3 + ZimBwd * HBwd^3 )
1967 // Compute numfluxBwd[dir] = (eBwd^[dir] \cdot n \times e^3) * (ZimFwd *
1968 // HFwd^3 + ZimBwd * HBwd^3 )
1969 ComputeNtimesFz(0, Fwd, Bwd, m_ZimFwd[0], m_ZimBwd[0], numfluxFwd[0],
1970 numfluxBwd[0]);
1971 ComputeNtimesFz(1, Fwd, Bwd, m_ZimFwd[1], m_ZimBwd[1], numfluxFwd[1],
1972 numfluxBwd[1]);
1973
1974 // Compute numfluxFwd[2] = eFwd^3 \cdot ( n1e1 \times ( imFwd EFwd + imBwd
1975 // EBwd ) ) / 2 {{Y_i}}
1976 ComputeNtimesF12(Fwd, Bwd, m_YimFwd[0], m_YimBwd[0], m_YimFwd[1],
1977 m_YimBwd[1], numfluxFwd[2], numfluxBwd[2]);
1978
1979 // Compute e1Fwd_cdot_ncrossdE = eFwd[dir] \cdot \alpha n \times n \times
1980 // [E] / 2 {{ZimFwd}}
1981 ComputeNtimestimesdFz(0, Fwd, Bwd, m_ZimFwd[0], m_ZimBwd[0],
1982 e1Fwd_cdot_ncrossdE, e1Bwd_cdot_ncrossdE);
1983 ComputeNtimestimesdFz(1, Fwd, Bwd, m_ZimFwd[1], m_ZimBwd[1],
1984 e2Fwd_cdot_ncrossdE, e2Bwd_cdot_ncrossdE);
1985
1986 // Compute - \alpha [H3] * ( 1/2{{Yim1}} + 1/2{{Yim2}} )
1987 ComputeNtimestimesdF12(Fwd, Bwd, m_YimFwd[0], m_YimBwd[0], m_YimFwd[1],
1988 m_YimBwd[1], e3Fwd_cdot_dHe3, e3Bwd_cdot_dHe3);
1989
1990 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxFwd[0], 1, e1Fwd_cdot_ncrossdE, 1,
1991 numfluxFwd[0], 1);
1992 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxFwd[1], 1, e2Fwd_cdot_ncrossdE, 1,
1993 numfluxFwd[1], 1);
1994 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxFwd[2], 1, e3Fwd_cdot_dHe3, 1,
1995 numfluxFwd[2], 1);
1996
1997 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxBwd[0], 1, e1Bwd_cdot_ncrossdE, 1,
1998 numfluxBwd[0], 1);
1999 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxBwd[1], 1, e2Bwd_cdot_ncrossdE, 1,
2000 numfluxBwd[1], 1);
2001 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxBwd[2], 1, e3Bwd_cdot_dHe3, 1,
2002 numfluxBwd[2], 1);
2003}
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)
Definition: MMFSystem.cpp:945
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)
Definition: MMFSystem.cpp:909
SOLVER_UTILS_EXPORT Array< OneD, NekDouble > GetIncidentField(const int var, const NekDouble time)
Definition: MMFSystem.cpp:2005
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)
Definition: MMFSystem.cpp:1074
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)
Definition: MMFSystem.cpp:1008
void Svtvp(int n, const T alpha, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Svtvp (scalar times vector plus vector): z = alpha*x + y.
Definition: Vmath.hpp:396

References ComputeNtimesF12(), ComputeNtimesFz(), ComputeNtimestimesdF12(), ComputeNtimestimesdFz(), CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, GetIncidentField(), Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_spacedim, m_YimBwd, m_YimFwd, m_ZimBwd, m_ZimFwd, and Vmath::Svtvp().

Referenced by NumericalMaxwellFlux().

◆ NumericalMaxwellFluxTM()

void Nektar::SolverUtils::MMFSystem::NumericalMaxwellFluxTM ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd,
const NekDouble  time 
)
protected

Definition at line 1801 of file MMFSystem.cpp.

1805{
1806 int nq = m_fields[0]->GetNpoints();
1807 int nTraceNumPoints = GetTraceTotPoints();
1808 int nvar = physfield.size();
1809
1810 // get temporary arrays
1811 Array<OneD, Array<OneD, NekDouble>> Fwd(nvar);
1812 Array<OneD, Array<OneD, NekDouble>> Bwd(nvar);
1813 for (int i = 0; i < nvar; ++i)
1814 {
1815 Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1816 Bwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1817
1818 // get the physical values at the trace
1819 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd[i], Bwd[i]);
1820 }
1821
1822 // E^|| = 0 at the PEC boundaries vs. H^|| = 0 at PMC boundaries
1823 Array<OneD, NekDouble> IncField(nq, 0.0);
1824 Array<OneD, NekDouble> IncFieldBwd(nTraceNumPoints, 0.0);
1825 Array<OneD, Array<OneD, NekDouble>> IncFieldFwd(m_spacedim);
1826 for (int i = 0; i < m_spacedim; ++i)
1827 {
1828 IncFieldFwd[i] = Array<OneD, NekDouble>(nTraceNumPoints, 0.0);
1829
1830 IncField = GetIncidentField(i, time);
1831 m_fields[0]->GetFwdBwdTracePhys(IncField, IncFieldFwd[i], IncFieldBwd);
1832
1833 Vmath::Svtvp(nTraceNumPoints, 2.0, &IncFieldFwd[i][0], 1, &Fwd[i][0], 1,
1834 &IncFieldFwd[i][0], 1);
1835 }
1836
1837 // Total Field Formulation
1838 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQBwd, 0, "PEC_Forces");
1839 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PEC_Forces");
1840 CopyBoundaryTrace(IncFieldFwd[2], Bwd[2], SolverUtils::eFwdEQNegBwd, 2,
1841 "PEC_Forces");
1842
1843 CopyBoundaryTrace(IncFieldFwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0,
1844 "PMC_Forces");
1845 CopyBoundaryTrace(IncFieldFwd[1], Bwd[1], SolverUtils::eFwdEQNegBwd, 1,
1846 "PMC_Forces");
1847 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQBwd, 2, "PMC_Forces");
1848
1849 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQBwd, 0, "PEC");
1850 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PEC");
1851 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQNegBwd, 2, "PEC");
1852
1853 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0, "PMC");
1854 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQNegBwd, 1, "PMC");
1855 CopyBoundaryTrace(Fwd[2], Bwd[2], SolverUtils::eFwdEQBwd, 2, "PMC");
1856
1857 Array<OneD, NekDouble> e1Fwd_cdot_ncrossdH(nTraceNumPoints, 0.0);
1858 Array<OneD, NekDouble> e1Bwd_cdot_ncrossdH(nTraceNumPoints, 0.0);
1859 Array<OneD, NekDouble> e2Fwd_cdot_ncrossdH(nTraceNumPoints, 0.0);
1860 Array<OneD, NekDouble> e2Bwd_cdot_ncrossdH(nTraceNumPoints, 0.0);
1861 Array<OneD, NekDouble> e3Fwd_cdot_dEe3(nTraceNumPoints, 0.0);
1862 Array<OneD, NekDouble> e3Bwd_cdot_dEe3(nTraceNumPoints, 0.0);
1863
1864 // Compute numfluxFwd[dir] = (eFwd^[dir] \cdot n \times e^3) * (YimFwd *
1865 // EFwd^3 + YimBwd * EBwd^3 )
1866 ComputeNtimesFz(0, Fwd, Bwd, m_YimFwd[0], m_YimBwd[0], numfluxFwd[0],
1867 numfluxBwd[0]);
1868 ComputeNtimesFz(1, Fwd, Bwd, m_YimFwd[1], m_YimBwd[1], numfluxFwd[1],
1869 numfluxBwd[1]);
1870
1871 // Compute numfluxFwd[2] = eFwd^3 \cdot ( n1e1 \times ( ZimFwd HFwd + ZimBwd
1872 // HBwd ) ) / 2 {{Z_i}}
1873 ComputeNtimesF12(Fwd, Bwd, m_ZimFwd[0], m_ZimBwd[0], m_ZimFwd[1],
1874 m_ZimBwd[1], numfluxFwd[2], numfluxBwd[2]);
1875
1876 // Compute e1Fwd_cdot_ncrossdE = eFwd[dir] \cdot \alpha n \times n \times
1877 // [H] / 2 {{YimFwd}}
1878 ComputeNtimestimesdFz(0, Fwd, Bwd, m_YimFwd[0], m_YimBwd[0],
1879 e1Fwd_cdot_ncrossdH, e1Bwd_cdot_ncrossdH);
1880 ComputeNtimestimesdFz(1, Fwd, Bwd, m_YimFwd[1], m_YimBwd[1],
1881 e2Fwd_cdot_ncrossdH, e2Bwd_cdot_ncrossdH);
1882
1883 // Compute \alpha [E3] * ( 1/2{{Zim1}} + 1/2{{Zim2}} )
1884 ComputeNtimestimesdF12(Fwd, Bwd, m_ZimFwd[0], m_ZimBwd[0], m_ZimFwd[1],
1885 m_ZimBwd[1], e3Fwd_cdot_dEe3, e3Bwd_cdot_dEe3);
1886
1887 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxFwd[0], 1, e1Fwd_cdot_ncrossdH,
1888 1, numfluxFwd[0], 1);
1889 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxFwd[1], 1, e2Fwd_cdot_ncrossdH,
1890 1, numfluxFwd[1], 1);
1891 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxFwd[2], 1, e3Fwd_cdot_dEe3, 1,
1892 numfluxFwd[2], 1);
1893
1894 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxBwd[0], 1, e1Bwd_cdot_ncrossdH,
1895 1, numfluxBwd[0], 1);
1896 Vmath::Svtvp(nTraceNumPoints, -1.0, numfluxBwd[1], 1, e2Bwd_cdot_ncrossdH,
1897 1, numfluxBwd[1], 1);
1898 Vmath::Svtvp(nTraceNumPoints, 1.0, numfluxBwd[2], 1, e3Bwd_cdot_dEe3, 1,
1899 numfluxBwd[2], 1);
1900}

References ComputeNtimesF12(), ComputeNtimesFz(), ComputeNtimestimesdF12(), ComputeNtimestimesdFz(), CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, GetIncidentField(), Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_spacedim, m_YimBwd, m_YimFwd, m_ZimBwd, m_ZimFwd, and Vmath::Svtvp().

Referenced by NumericalMaxwellFlux().

◆ RootMeanSquare()

NekDouble Nektar::SolverUtils::MMFSystem::RootMeanSquare ( const Array< OneD, const NekDouble > &  inarray)
protected

Definition at line 2338 of file MMFSystem.cpp.

2339{
2340 int Ntot = inarray.size();
2341
2342 NekDouble reval = 0.0;
2343 for (int i = 0; i < Ntot; ++i)
2344 {
2345 reval = reval + inarray[i] * inarray[i];
2346 }
2347 reval = sqrt(reval / Ntot);
2348
2349 return reval;
2350}

References tinysimd::sqrt().

Referenced by Computedemdxicdote(), ComputencdotMF(), ComputeZimYim(), Nektar::MMFAdvection::EvaluateAdvectionVelocity(), Nektar::MMFMaxwell::v_DoSolve(), Nektar::MMFAdvection::v_InitObject(), Nektar::MMFMaxwell::v_InitObject(), and VectorAvgMagnitude().

◆ SetUpMovingFrames()

void Nektar::SolverUtils::MMFSystem::SetUpMovingFrames ( const Array< OneD, const Array< OneD, NekDouble > > &  Anisotropy,
const int  TangentXelem 
)
protected

Definition at line 118 of file MMFSystem.cpp.

121{
122
123 int nq = m_fields[0]->GetNpoints();
124 int MFdim = 3;
125
126 // Construct The Moving Frames
127 m_movingframes = Array<OneD, Array<OneD, NekDouble>>(MFdim);
128
129 for (int j = 0; j < MFdim; ++j)
130 {
131 m_movingframes[j] = Array<OneD, NekDouble>(m_spacedim * nq, 0.0);
132 }
133
134 // Read MMF Geom Info
135 std::string conn = "TangentX";
136 m_MMFfactors = Array<OneD, NekDouble>(4);
137
138 m_session->LoadSolverInfo("MMFDir", conn, "LOCAL");
139
140 // (x-x_0)^2/a^2 + (y-y_0)^2/b^2 = 1
141 // factors[0] = a
142 // factors[1] = b
143 // factors[2] = x_0
144 // factors[3] = y_0
145 m_session->LoadParameter("MMFCircAxisX", m_MMFfactors[0], 1.0);
146 m_session->LoadParameter("MMFCircAxisY", m_MMFfactors[1], 1.0);
147 m_session->LoadParameter("MMFCircCentreX", m_MMFfactors[2], 0.0);
148 m_session->LoadParameter("MMFCircCentreY", m_MMFfactors[3], 0.0);
149
150 if (conn == "TangentX")
151 {
153 }
154 if (conn == "TangentY")
155 {
157 }
158 if (conn == "TangentXY")
159 {
161 }
162 if (conn == "TangentZ")
163 {
165 }
166 if (conn == "TangentCircular")
167 {
169 }
170 if (conn == "TangentIrregular")
171 {
173 }
174 if (conn == "TangentNonconvex")
175 {
177 }
178 if (conn == "LOCAL")
179 {
181 }
182
183 // Get Tangetn vectors from GeomFactors2D, Orthonormalized = true
184 m_fields[0]->GetMovingFrames(m_MMFdir, m_MMFfactors, m_movingframes);
185
186 // Align the tangentX direction after TangentXelem
187 if (TangentXelem > 0)
188 {
189 Array<OneD, NekDouble> tmp(nq);
190 m_fields[0]->GenerateElementVector(TangentXelem, 1.0, 0.0, tmp);
191 for (int j = 0; j < m_shapedim; j++)
192 {
193 for (int k = 0; k < m_spacedim; k++)
194 {
195 Vmath::Vmul(nq, &tmp[0], 1, &m_movingframes[j][k * nq], 1,
196 &m_movingframes[j][k * nq], 1);
197 }
198 }
199
200 m_fields[0]->GenerateElementVector(TangentXelem, 0.0, 1.0, tmp);
201 Vmath::Vadd(nq, &tmp[0], 1, &m_movingframes[0][0 * nq], 1,
202 &m_movingframes[0][0 * nq], 1);
203 Vmath::Vadd(nq, &tmp[0], 1, &m_movingframes[1][1 * nq], 1,
204 &m_movingframes[1][1 * nq], 1);
205
206 int indxtmp = Vmath::Imax(nq, tmp, 1);
207 std::cout << "*** MF in PML Region is aligned as MF1 = ( "
208 << m_movingframes[0][indxtmp] << " , "
209 << m_movingframes[0][nq + indxtmp] << " , "
210 << m_movingframes[0][2 * nq + indxtmp] << " ) "
211 << ", MF2 = ( " << m_movingframes[1][indxtmp] << " , "
212 << m_movingframes[1][nq + indxtmp] << " , "
213 << m_movingframes[1][2 * nq + indxtmp] << " ) " << std::endl;
214 }
215
216 // Multiply Anisotropy to movingframes
217 for (int j = 0; j < m_shapedim; ++j)
218 {
219 for (int k = 0; k < m_spacedim; ++k)
220 {
221 Vmath::Vmul(nq, &Anisotropy[j][0], 1, &m_movingframes[j][k * nq], 1,
222 &m_movingframes[j][k * nq], 1);
223 }
224 }
225
226 // Test the moving frames
228}
void CheckMovingFrames(const Array< OneD, const Array< OneD, NekDouble > > &movingframes)
Definition: MMFSystem.cpp:230
Array< OneD, NekDouble > m_MMFfactors
Definition: MMFSystem.h:157
@ eLOCAL
No Principal direction.
@ eTangentIrregular
Circular around the centre of domain.
@ eTangentX
X coordinate direction.
@ eTangentCircular
Circular around the centre of domain.
@ eTangentNonconvex
Circular around the centre of domain.
@ eTangentXY
XY direction.
@ eTangentZ
Z coordinate direction.
@ eTangentY
Y coordinate direction.
int Imax(int n, const T *x, const int incx)
Return the index of the maximum element in x.
Definition: Vmath.hpp:623

References CheckMovingFrames(), Nektar::SpatialDomains::eLOCAL, Nektar::SpatialDomains::eTangentCircular, Nektar::SpatialDomains::eTangentIrregular, Nektar::SpatialDomains::eTangentNonconvex, Nektar::SpatialDomains::eTangentX, Nektar::SpatialDomains::eTangentXY, Nektar::SpatialDomains::eTangentY, Nektar::SpatialDomains::eTangentZ, Vmath::Imax(), Nektar::SolverUtils::EquationSystem::m_fields, m_MMFdir, m_MMFfactors, m_movingframes, Nektar::SolverUtils::EquationSystem::m_session, m_shapedim, Nektar::SolverUtils::EquationSystem::m_spacedim, Vmath::Vadd(), and Vmath::Vmul().

Referenced by MMFInitObject().

◆ UpwindMaxwellFlux1D()

void Nektar::SolverUtils::MMFSystem::UpwindMaxwellFlux1D ( Array< OneD, Array< OneD, NekDouble > > &  physfield,
Array< OneD, Array< OneD, NekDouble > > &  numfluxFwd,
Array< OneD, Array< OneD, NekDouble > > &  numfluxBwd 
)
protected

Definition at line 1455 of file MMFSystem.cpp.

1459{
1460 int i;
1461 int nTraceNumPoints = GetTraceTotPoints();
1462 int nvar = 2;
1463
1464 // get temporary arrays
1465 Array<OneD, Array<OneD, NekDouble>> Fwd(nvar);
1466 Array<OneD, Array<OneD, NekDouble>> Bwd(nvar);
1467
1468 for (i = 0; i < nvar; ++i)
1469 {
1470 Fwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1471 Bwd[i] = Array<OneD, NekDouble>(nTraceNumPoints);
1472 }
1473
1474 // get the physical values at the trace from the dependent variables
1475 for (i = 0; i < nvar; ++i)
1476 {
1477 m_fields[i]->GetFwdBwdTracePhys(physfield[i], Fwd[i], Bwd[i]);
1478 }
1479
1480 // E = 0 at the boundaries
1481 CopyBoundaryTrace(Fwd[0], Bwd[0], SolverUtils::eFwdEQNegBwd, 0, "PEC");
1482
1483 // d H / d n = 0 at the boundaries
1484 CopyBoundaryTrace(Fwd[1], Bwd[1], SolverUtils::eFwdEQBwd, 1, "PEC");
1485
1486 Array<OneD, NekDouble> dE(nTraceNumPoints);
1487 Array<OneD, NekDouble> dH(nTraceNumPoints);
1488
1489 Vmath::Vsub(nTraceNumPoints, &Fwd[0][0], 1, &Bwd[0][0], 1, &dE[0], 1);
1490 Vmath::Vsub(nTraceNumPoints, &Fwd[1][0], 1, &Bwd[1][0], 1, &dH[0], 1);
1491
1492 NekDouble nx, AverZ, AverY, AverZH, AverYE;
1493 for (i = 0; i < nTraceNumPoints; ++i)
1494 {
1495 nx = m_traceNormals[0][i];
1496 AverZ = m_ZimFwd[0][i] + m_ZimBwd[0][i];
1497 AverY = m_YimFwd[0][i] + m_YimBwd[0][i];
1498 AverZH = m_ZimFwd[0][i] * Fwd[1][i] + m_ZimBwd[0][i] * Bwd[1][i];
1499 AverYE = m_YimFwd[0][i] * Fwd[0][i] + m_YimBwd[0][i] * Bwd[0][i];
1500
1501 numfluxFwd[0][i] = nx / AverZ * (AverZH - nx * dE[i]);
1502 numfluxFwd[1][i] = nx / AverY * (AverYE - nx * dH[i]);
1503
1504 numfluxBwd[0][i] = nx / AverZ * (AverZH - nx * dE[i]);
1505 numfluxBwd[1][i] = nx / AverY * (AverYE - nx * dH[i]);
1506 }
1507}

References CopyBoundaryTrace(), Nektar::SolverUtils::eFwdEQBwd, Nektar::SolverUtils::eFwdEQNegBwd, Nektar::SolverUtils::EquationSystem::GetTraceTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_traceNormals, m_YimBwd, m_YimFwd, m_ZimBwd, m_ZimFwd, and Vmath::Vsub().

Referenced by NumericalMaxwellFlux().

◆ v_GenerateSummary()

void Nektar::SolverUtils::MMFSystem::v_GenerateSummary ( SummaryList l)
overridevirtual

Virtual function for generating summary information.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Reimplemented in Nektar::MMFAdvection, Nektar::MMFDiffusion, Nektar::MMFMaxwell, and Nektar::MMFSWE.

Definition at line 2463 of file MMFSystem.cpp.

2464{
2465 int nq = m_fields[0]->GetNpoints();
2467
2468 AddSummaryItem(s, "Total grids", nq);
2469 AddSummaryItem(s, "Shape Dimension", m_shapedim);
2471 if (m_surfaceType == eSphere)
2472 {
2473 NekDouble MeshError;
2474
2475 Array<OneD, NekDouble> x(nq);
2476 Array<OneD, NekDouble> y(nq);
2477 Array<OneD, NekDouble> z(nq);
2478 Array<OneD, NekDouble> rad(nq, 0.0);
2479
2480 m_fields[0]->GetCoords(x, y, z);
2481
2482 Vmath::Vvtvp(nq, x, 1, x, 1, rad, 1, rad, 1);
2483 Vmath::Vvtvp(nq, y, 1, y, 1, rad, 1, rad, 1);
2484 Vmath::Vvtvp(nq, z, 1, z, 1, rad, 1, rad, 1);
2485 Vmath::Vsqrt(nq, rad, 1, rad, 1);
2486
2487 Vmath::Sadd(nq, -1.0, rad, 1, rad, 1);
2488 Vmath::Vabs(nq, rad, 1, rad, 1);
2489
2490 MeshError = m_fields[0]->Integral(rad);
2491 SolverUtils::AddSummaryItem(s, "Mesh Error", MeshError);
2492 }
2493
2495 AddSummaryItem(s, "SmoothFactor", m_SmoothFactor);
2496 AddSummaryItem(s, "SFinit", m_SFinit);
2497
2498 if (fabs(m_alpha - 1.0) > 0.001)
2499 {
2500 AddSummaryItem(s, "Alpha", m_alpha);
2501 }
2502
2503 if (m_Incfreq > 0.0)
2504 {
2505 AddSummaryItem(s, "Incfreq", m_Incfreq);
2506 }
2507
2508 if (m_MMFdir == 4)
2509 {
2510 AddSummaryItem(s, "MMFCircAxisX", m_MMFfactors[0]);
2511 AddSummaryItem(s, "MMFCircAxisY", m_MMFfactors[1]);
2512 AddSummaryItem(s, "MMFCircCentreX", m_MMFfactors[2]);
2513 AddSummaryItem(s, "MMFCircCentreY", m_MMFfactors[3]);
2514 }
2515}
SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Print a summary of time stepping parameters.
static NekDouble rad(NekDouble x, NekDouble y)
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
std::vector< double > z(NPUPPER)
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add vector y = alpha + x.
Definition: Vmath.hpp:194

References Nektar::SolverUtils::AddSummaryItem(), Nektar::SolverUtils::eSphere, Nektar::SpatialDomains::GeomMMFMap, m_alpha, Nektar::SolverUtils::EquationSystem::m_fields, m_Incfreq, m_MMFdir, m_MMFfactors, m_SFinit, m_shapedim, m_SmoothFactor, m_surfaceType, Nektar::LibUtilities::rad(), Vmath::Sadd(), Nektar::SolverUtils::SurfaceTypeMap, Nektar::SolverUtils::UnsteadySystem::v_GenerateSummary(), Vmath::Vabs(), Vmath::Vsqrt(), Vmath::Vvtvp(), and Nektar::UnitTests::z().

Referenced by Nektar::MMFAdvection::v_GenerateSummary(), Nektar::MMFDiffusion::v_GenerateSummary(), Nektar::MMFMaxwell::v_GenerateSummary(), and Nektar::MMFSWE::v_GenerateSummary().

◆ VectorAvgMagnitude()

NekDouble Nektar::SolverUtils::MMFSystem::VectorAvgMagnitude ( const Array< OneD, const Array< OneD, NekDouble > > &  inarray)
protected

Definition at line 2352 of file MMFSystem.cpp.

2354{
2355 int nq = inarray[0].size();
2356
2357 Array<OneD, NekDouble> tmp(nq, 0.0);
2358 for (int k = 0; k < m_spacedim; k++)
2359 {
2360 Vmath::Vvtvp(nq, &inarray[k][0], 1, &inarray[k][0], 1, &tmp[0], 1,
2361 &tmp[0], 1);
2362 }
2363 Vmath::Vsqrt(nq, tmp, 1, tmp, 1);
2364
2365 return RootMeanSquare(tmp);
2366}

References Nektar::SolverUtils::EquationSystem::m_spacedim, RootMeanSquare(), Vmath::Vsqrt(), and Vmath::Vvtvp().

Referenced by ComputeMFtrace(), and ComputeNtimesMF().

◆ VectorCrossProd() [1/2]

void Nektar::SolverUtils::MMFSystem::VectorCrossProd ( const Array< OneD, const Array< OneD, NekDouble > > &  v1,
const Array< OneD, const Array< OneD, NekDouble > > &  v2,
Array< OneD, Array< OneD, NekDouble > > &  v3 
)
protected

Computes the vector cross-product in 3D of v1 and v2, storing the result in v3.

Parameters
v1First input vector.
v2Second input vector.
v3Output vector computed to be orthogonal to both v1 and v2.

Definition at line 696 of file MMFSystem.cpp.

700{
701 ASSERTL0(v1.size() == 3, "Input 1 has dimension not equal to 3.");
702 ASSERTL0(v2.size() == 3, "Input 2 has dimension not equal to 3.");
703 ASSERTL0(v3.size() == 3, "Output vector has dimension not equal to 3.");
704
705 int nq = v1[0].size();
706 Array<OneD, NekDouble> temp(nq);
707
708 Vmath::Vmul(nq, v1[2], 1, v2[1], 1, temp, 1);
709 Vmath::Vvtvm(nq, v1[1], 1, v2[2], 1, temp, 1, v3[0], 1);
710
711 Vmath::Vmul(nq, v1[0], 1, v2[2], 1, temp, 1);
712 Vmath::Vvtvm(nq, v1[2], 1, v2[0], 1, temp, 1, v3[1], 1);
713
714 Vmath::Vmul(nq, v1[1], 1, v2[0], 1, temp, 1);
715 Vmath::Vvtvm(nq, v1[0], 1, v2[1], 1, temp, 1, v3[2], 1);
716}
void Vvtvm(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)
vvtvm (vector times vector minus vector): z = w*x - y
Definition: Vmath.hpp:381

References ASSERTL0, Vmath::Vmul(), and Vmath::Vvtvm().

Referenced by ComputencdotMF(), ComputeNtimesF12(), ComputeNtimesMF(), ComputeNtimestimesdFz(), DeriveCrossProductMF(), and Nektar::MMFAdvection::EvaluateAdvectionVelocity().

◆ VectorCrossProd() [2/2]

void Nektar::SolverUtils::MMFSystem::VectorCrossProd ( const Array< OneD, NekDouble > &  v1,
const Array< OneD, NekDouble > &  v2,
Array< OneD, NekDouble > &  v3 
)
protected

Definition at line 718 of file MMFSystem.cpp.

721{
722 ASSERTL0(v1.size() == 3, "Input 1 has dimension not equal to 3.");
723 ASSERTL0(v2.size() == 3, "Input 2 has dimension not equal to 3.");
724 ASSERTL0(v3.size() == 3, "Output vector has dimension not equal to 3.");
725
726 v3[0] = v1[1] * v2[2] - v1[2] * v2[1];
727 v3[1] = v1[2] * v2[0] - v1[0] * v2[2];
728 v3[2] = v1[0] * v2[1] - v1[1] * v2[0];
729}

References ASSERTL0.

◆ VectorDotProd()

void Nektar::SolverUtils::MMFSystem::VectorDotProd ( const Array< OneD, const Array< OneD, NekDouble > > &  v1,
const Array< OneD, const Array< OneD, NekDouble > > &  v2,
Array< OneD, NekDouble > &  v3 
)
protected

Definition at line 673 of file MMFSystem.cpp.

677{
678 int coordim = v1.size();
679 int nq = v1[0].size();
680
681 v3 = Array<OneD, NekDouble>(nq, 0.0);
682 for (int i = 0; i < coordim; ++i)
683 {
684 Vmath::Vvtvp(nq, &v1[i][0], 1, &v2[i][0], 1, &v3[0], 1, &v3[0], 1);
685 }
686}

References Vmath::Vvtvp().

Referenced by ComputencdotMF().

Member Data Documentation

◆ m_alpha

NekDouble Nektar::SolverUtils::MMFSystem::m_alpha
protected

◆ m_CurlMF

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_CurlMF
protected

◆ m_dedxi_cdot_e

Array<OneD, Array<OneD, Array<OneD, Array<OneD, NekDouble> > > > Nektar::SolverUtils::MMFSystem::m_dedxi_cdot_e
protected

Definition at line 217 of file MMFSystem.h.

Referenced by AdddedtMaxwell(), and Computedemdxicdote().

◆ m_DivMF

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_DivMF
protected

◆ m_epsvec

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_epsvec
protected

◆ m_Incfreq

NekDouble Nektar::SolverUtils::MMFSystem::m_Incfreq
protected

◆ m_IncType

IncType Nektar::SolverUtils::MMFSystem::m_IncType

◆ m_MFtraceBwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_MFtraceBwd
protected

◆ m_MFtraceFwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_MFtraceFwd
protected

◆ m_MMFdir

SpatialDomains::GeomMMF Nektar::SolverUtils::MMFSystem::m_MMFdir
protected

◆ m_MMFfactors

Array<OneD, NekDouble> Nektar::SolverUtils::MMFSystem::m_MMFfactors

◆ m_movingframes

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_movingframes
protected

◆ m_muvec

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_muvec
protected

◆ m_ncdotMFBwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_ncdotMFBwd
protected

◆ m_ncdotMFFwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_ncdotMFFwd
protected

◆ m_negepsvecminus1

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_negepsvecminus1
protected

Definition at line 212 of file MMFSystem.h.

Referenced by Nektar::MMFMaxwell::DoOdeRhs(), and Nektar::MMFMaxwell::v_InitObject().

◆ m_negmuvecminus1

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_negmuvecminus1
protected

Definition at line 213 of file MMFSystem.h.

Referenced by Nektar::MMFMaxwell::DoOdeRhs(), and Nektar::MMFMaxwell::v_InitObject().

◆ m_nperpcdotMFBwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_nperpcdotMFBwd
protected

Definition at line 190 of file MMFSystem.h.

Referenced by ComputencdotMF(), and Nektar::MMFSWE::NumericalSWEFlux().

◆ m_nperpcdotMFFwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_nperpcdotMFFwd
protected

◆ m_ntimes_ntimesMFBwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_ntimes_ntimesMFBwd
protected

Definition at line 202 of file MMFSystem.h.

Referenced by ComputeNtimesMF().

◆ m_ntimes_ntimesMFFwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_ntimes_ntimesMFFwd
protected

Definition at line 201 of file MMFSystem.h.

Referenced by ComputeNtimesMF().

◆ m_ntimesMFBwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_ntimesMFBwd
protected

Definition at line 200 of file MMFSystem.h.

Referenced by ComputeNtimesFz(), and ComputeNtimesMF().

◆ m_ntimesMFFwd

Array<OneD, Array<OneD, Array<OneD, NekDouble> > > Nektar::SolverUtils::MMFSystem::m_ntimesMFFwd
protected

◆ m_pi

NekDouble Nektar::SolverUtils::MMFSystem::m_pi

◆ m_PolType

PolType Nektar::SolverUtils::MMFSystem::m_PolType

◆ m_SFinit

NekDouble Nektar::SolverUtils::MMFSystem::m_SFinit
protected

Definition at line 181 of file MMFSystem.h.

Referenced by GetIncidentField(), MMFInitObject(), and v_GenerateSummary().

◆ m_shapedim

int Nektar::SolverUtils::MMFSystem::m_shapedim

◆ m_SmoothFactor

int Nektar::SolverUtils::MMFSystem::m_SmoothFactor
protected

Definition at line 180 of file MMFSystem.h.

Referenced by GetIncidentField(), MMFInitObject(), and v_GenerateSummary().

◆ m_surfaceNormal

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_surfaceNormal
protected

Definition at line 184 of file MMFSystem.h.

◆ m_surfaceType

SurfaceType Nektar::SolverUtils::MMFSystem::m_surfaceType

◆ m_TestMaxwellType

TestMaxwellType Nektar::SolverUtils::MMFSystem::m_TestMaxwellType

◆ m_upwindType

UpwindType Nektar::SolverUtils::MMFSystem::m_upwindType

◆ m_YimBwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_YimBwd
protected

◆ m_YimFwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_YimFwd
protected

◆ m_ZimBwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_ZimBwd
protected

◆ m_ZimFwd

Array<OneD, Array<OneD, NekDouble> > Nektar::SolverUtils::MMFSystem::m_ZimFwd
protected