Nektar++
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Public Member Functions | Static Public Member Functions | Public Attributes | Static Public Attributes | Protected Member Functions | Protected Attributes | Private Attributes | Friends | List of all members
Nektar::MMFDiffusion Class Reference

A model for cardiac conduction. More...

#include <MMFDiffusion.h>

Inheritance diagram for Nektar::MMFDiffusion:
[legend]

Public Member Functions

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

Static Public Member Functions

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

Public Attributes

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

Static Public Attributes

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

Protected Member Functions

 MMFDiffusion (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Constructor. More...
 
virtual void v_InitObject (bool DeclareField=true) override
 Init object for UnsteadySystem class. More...
 
void DoImplicitSolve (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, NekDouble time, NekDouble lambda)
 Solve for the diffusion term. More...
 
void DoOdeRhs (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
 Computes the reaction terms \(f(u,v)\) and \(g(u,v)\). More...
 
void TestPlaneProblem (const NekDouble time, Array< OneD, NekDouble > &outfield)
 
void TestCubeProblem (const NekDouble time, Array< OneD, NekDouble > &outfield)
 
void Morphogenesis (const NekDouble time, unsigned int field, Array< OneD, NekDouble > &outfield)
 
Array< OneD, NekDoublePlanePhiWave ()
 
virtual void v_SetInitialConditions (NekDouble initialtime, bool dumpInitialConditions, const int domain) override
 Sets a custom initial condition. More...
 
virtual void v_GenerateSummary (SolverUtils::SummaryList &s) override
 Prints a summary of the model parameters. More...
 
virtual void v_EvaluateExactSolution (unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time) override
 
- Protected Member Functions inherited from Nektar::SolverUtils::MMFSystem
void SetUpMovingFrames (const Array< OneD, const Array< OneD, NekDouble >> &Anisotropy, const int TangentXelem)
 
void CheckMovingFrames (const Array< OneD, const Array< OneD, NekDouble >> &movingframes)
 
SOLVER_UTILS_EXPORT void ComputencdotMF ()
 
SOLVER_UTILS_EXPORT void ComputeDivCurlMF ()
 
SOLVER_UTILS_EXPORT void ComputeMFtrace ()
 
SOLVER_UTILS_EXPORT void VectorDotProd (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, NekDouble > &v3)
 
SOLVER_UTILS_EXPORT void VectorCrossProd (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &v3)
 
SOLVER_UTILS_EXPORT void VectorCrossProd (const Array< OneD, NekDouble > &v1, const Array< OneD, NekDouble > &v2, Array< OneD, NekDouble > &v3)
 
SOLVER_UTILS_EXPORT void ComputeCurl (const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleCartesianToMovingframes (const Array< OneD, const Array< OneD, NekDouble >> &uvec, unsigned int field)
 
SOLVER_UTILS_EXPORT void DeriveCrossProductMF (Array< OneD, Array< OneD, NekDouble >> &CrossProductMF)
 
SOLVER_UTILS_EXPORT void ComputeNtimesMF ()
 
SOLVER_UTILS_EXPORT void ComputeNtimesFz (const int dir, const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &imFwd, const Array< OneD, const NekDouble > &imBwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimesF12 (const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &im1Fwd, const Array< OneD, const NekDouble > &im1Bwd, const Array< OneD, const NekDouble > &im2Fwd, const Array< OneD, const NekDouble > &im2Bwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimestimesdFz (const int dir, const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &imFwd, const Array< OneD, const NekDouble > &imBwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void ComputeNtimestimesdF12 (const Array< OneD, Array< OneD, NekDouble >> &Fwd, const Array< OneD, Array< OneD, NekDouble >> &Bwd, const Array< OneD, const NekDouble > &im1Fwd, const Array< OneD, const NekDouble > &im1Bwd, const Array< OneD, const NekDouble > &im2Fwd, const Array< OneD, const NekDouble > &im2Bwd, Array< OneD, NekDouble > &outarrayFwd, Array< OneD, NekDouble > &outarrayBwd)
 
SOLVER_UTILS_EXPORT void CartesianToSpherical (const NekDouble x0j, const NekDouble x1j, const NekDouble x2j, NekDouble &sin_varphi, NekDouble &cos_varphi, NekDouble &sin_theta, NekDouble &cos_theta)
 
SOLVER_UTILS_EXPORT void ComputeZimYim (Array< OneD, Array< OneD, NekDouble >> &epsvec, Array< OneD, Array< OneD, NekDouble >> &muvec)
 
SOLVER_UTILS_EXPORT void AdddedtMaxwell (const Array< OneD, const Array< OneD, NekDouble >> &physarray, Array< OneD, Array< OneD, NekDouble >> &outarray)
 
SOLVER_UTILS_EXPORT void GetMaxwellFluxVector (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void GetMaxwellFlux1D (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void GetMaxwellFlux2D (const int var, const Array< OneD, const Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &flux)
 
SOLVER_UTILS_EXPORT void LaxFriedrichMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void UpwindMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void AverageMaxwellFlux1D (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFlux (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time=0.0)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFluxTM (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time)
 
SOLVER_UTILS_EXPORT void NumericalMaxwellFluxTE (Array< OneD, Array< OneD, NekDouble >> &physfield, Array< OneD, Array< OneD, NekDouble >> &numfluxFwd, Array< OneD, Array< OneD, NekDouble >> &numfluxBwd, const NekDouble time)
 
SOLVER_UTILS_EXPORT Array< OneD, NekDoubleGetIncidentField (const int var, const NekDouble time)
 
SOLVER_UTILS_EXPORT void Computedemdxicdote ()
 
SOLVER_UTILS_EXPORT NekDouble AvgInt (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble AvgAbsInt (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble AbsIntegral (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble RootMeanSquare (const Array< OneD, const NekDouble > &inarray)
 
SOLVER_UTILS_EXPORT NekDouble VectorAvgMagnitude (const Array< OneD, const Array< OneD, NekDouble >> &inarray)
 
SOLVER_UTILS_EXPORT void GramSchumitz (const Array< OneD, const Array< OneD, NekDouble >> &v1, const Array< OneD, const Array< OneD, NekDouble >> &v2, Array< OneD, Array< OneD, NekDouble >> &outarray, bool KeepTheMagnitude=true)
 
SOLVER_UTILS_EXPORT void BubbleSort (Array< OneD, NekDouble > &refarray, Array< OneD, NekDouble > &sortarray)
 
- Protected Member Functions inherited from Nektar::SolverUtils::UnsteadySystem
SOLVER_UTILS_EXPORT UnsteadySystem (const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
 Initialises UnsteadySystem class members. More...
 
SOLVER_UTILS_EXPORT NekDouble MaxTimeStepEstimator ()
 Get the maximum timestep estimator for cfl control. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoSolve () override
 Solves an unsteady problem. More...
 
virtual SOLVER_UTILS_EXPORT void v_DoInitialise () override
 Sets up initial conditions. 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)
 
SOLVER_UTILS_EXPORT void CheckForRestartTime (NekDouble &time, int &nchk)
 
SOLVER_UTILS_EXPORT void SVVVarDiffCoeff (const Array< OneD, Array< OneD, NekDouble >> vel, StdRegions::VarCoeffMap &varCoeffMap)
 Evaluate the SVV diffusion coefficient according to Moura's paper where it should proportional to h time velocity. More...
 
virtual SOLVER_UTILS_EXPORT bool v_UpdateTimeStepCheck ()
 
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 NekDouble v_LinfError (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray)
 Virtual function for the L_inf error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT NekDouble v_L2Error (unsigned int field, const Array< OneD, NekDouble > &exactsoln=NullNekDouble1DArray, bool Normalised=false)
 Virtual function for the L_2 error computation between fields and a given exact solution. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransCoeffToPhys ()
 Virtual function for transformation to physical space. More...
 
virtual SOLVER_UTILS_EXPORT void v_TransPhysToCoeff ()
 Virtual function for transformation to coefficient space. More...
 
virtual SOLVER_UTILS_EXPORT void v_Output (void)
 
virtual SOLVER_UTILS_EXPORT MultiRegions::ExpListSharedPtr v_GetPressure (void)
 
virtual SOLVER_UTILS_EXPORT void v_ExtraFldOutput (std::vector< Array< OneD, NekDouble >> &fieldcoeffs, std::vector< std::string > &variables)
 

Protected Attributes

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

Private Attributes

StdRegions::VarCoeffMap m_varcoeff
 Variable diffusivity. More...
 
Array< OneD, NekDoublem_epsilon
 
Array< OneD, NekDoublem_epsu
 

Friends

class MemoryManager< MMFDiffusion >
 

Additional Inherited Members

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

Detailed Description

A model for cardiac conduction.

Definition at line 79 of file MMFDiffusion.h.

Constructor & Destructor Documentation

◆ ~MMFDiffusion()

Nektar::MMFDiffusion::~MMFDiffusion ( )
virtual

Desctructor.

Definition at line 181 of file MMFDiffusion.cpp.

182 {
183 }

◆ MMFDiffusion()

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

Constructor.

Definition at line 57 of file MMFDiffusion.cpp.

59  : UnsteadySystem(pSession, pGraph), MMFSystem(pSession, pGraph)
60 {
61 }
SOLVER_UTILS_EXPORT MMFSystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Definition: MMFSystem.cpp:43
SOLVER_UTILS_EXPORT UnsteadySystem(const LibUtilities::SessionReaderSharedPtr &pSession, const SpatialDomains::MeshGraphSharedPtr &pGraph)
Initialises UnsteadySystem class members.

Member Function Documentation

◆ create()

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

Creates an instance of this class.

Definition at line 85 of file MMFDiffusion.h.

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

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

◆ DoImplicitSolve()

void Nektar::MMFDiffusion::DoImplicitSolve ( const Array< OneD, const Array< OneD, NekDouble >> &  inarray,
Array< OneD, Array< OneD, NekDouble >> &  outarray,
NekDouble  time,
NekDouble  lambda 
)
protected

Solve for the diffusion term.

OdeRhs

Parameters
inarrayInput array.
outarrayOutput array.
timeCurrent simulation time.
lambdaTimestep.

Definition at line 191 of file MMFDiffusion.cpp.

195 {
196  int nvariables = inarray.size();
197  int nq = m_fields[0]->GetNpoints();
198 
200  factors[StdRegions::eFactorTau] = 1.0;
201 
202  Array<OneD, Array<OneD, NekDouble>> F(nvariables);
203  factors[StdRegions::eFactorLambda] = 1.0 / lambda;
204  F[0] = Array<OneD, NekDouble>(nq * nvariables);
205 
206  for (int n = 1; n < nvariables; ++n)
207  {
208  F[n] = F[n - 1] + nq;
209  // cout << "F["<< n<<"=" << F[n][1] <<endl;
210  }
211 
212  // We solve ( \nabla^2 - HHlambda ) Y[i] = rhs [i]
213  // inarray = input: \hat{rhs} -> output: \hat{Y}
214  // outarray = output: nabla^2 \hat{Y}
215  // where \hat = modal coeffs
216  SetBoundaryConditions(time);
217 
218  for (int i = 0; i < nvariables; ++i)
219  {
220  factors[StdRegions::eFactorLambda] = 1.0 / lambda / m_epsu[i];
221 
222  // Multiply 1.0/timestep
223  Vmath::Smul(nq, -factors[StdRegions::eFactorLambda], inarray[i], 1,
224  F[i], 1);
225 
226  /* for (int k = 0; k < 15; ++k)
227  cout << "inarray["<<i << "]"<< k<<"=" << inarray[i][k]<<endl;*/
228  // Solve a system of equations with Helmholtz solver and transform
229  // back into physical space.
230  m_fields[i]->HelmSolve(F[i], m_fields[i]->UpdateCoeffs(), factors,
231  m_varcoeff);
232 
233  m_fields[i]->BwdTrans(m_fields[i]->GetCoeffs(), outarray[i]);
234  /* Array<OneD, NekDouble> coefarray = m_fields[i]->GetCoeffs();
235  for (int k = 0; k < 15; ++k)
236  cout << "inarray["<< k<<"=" << coefarray[k]<<endl;*/
237  }
238  /* for (int kk = 0; kk < 15; ++kk)
239  cout << "inarray["<< kk<<"=" <<
240  m_varcoeff[StdRegions::eVarCoeffMF3Mag][kk]<<endl;*/
241 }
StdRegions::VarCoeffMap m_varcoeff
Variable diffusivity.
Definition: MMFDiffusion.h:150
Array< OneD, NekDouble > m_epsu
Definition: MMFDiffusion.h:153
Array< OneD, MultiRegions::ExpListSharedPtr > m_fields
Array holding all dependent variables.
SOLVER_UTILS_EXPORT void SetBoundaryConditions(NekDouble time)
Evaluates the boundary conditions at the given time.
std::map< ConstFactorType, NekDouble > ConstFactorMap
Definition: StdRegions.hpp:399
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.cpp:248

References Nektar::StdRegions::eFactorLambda, Nektar::StdRegions::eFactorTau, m_epsu, Nektar::SolverUtils::EquationSystem::m_fields, m_varcoeff, Nektar::SolverUtils::EquationSystem::SetBoundaryConditions(), and Vmath::Smul().

Referenced by v_InitObject().

◆ DoOdeRhs()

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

Computes the reaction terms \(f(u,v)\) and \(g(u,v)\).

Definition at line 246 of file MMFDiffusion.cpp.

249 {
250  int nq = GetTotPoints();
251 
252  switch (m_TestType)
253  {
254  case eTestPlane:
255  {
256 
260 
261  m_fields[0]->GetCoords(x, y, z);
262 
263  for (int k = 0; k < nq; k++)
264  {
265  outarray[0][k] = (m_epsilon[0] + m_epsilon[1] - 1.0) * m_pi *
266  m_pi * exp(-1.0 * m_pi * m_pi * time) *
267  sin(m_pi * x[k]) * cos(m_pi * y[k]);
268  }
269  }
270  break;
271 
272  case eTestCube:
273  {
274 
278 
279  m_fields[0]->GetCoords(x, y, z);
280 
281  for (int k = 0; k < nq; k++)
282  {
283  outarray[0][k] =
284  (m_epsilon[0] + m_epsilon[1] + m_epsilon[2] - 1.0) * m_pi *
285  m_pi * exp(-1.0 * m_pi * m_pi * time) * sin(m_pi * x[k]) *
286  sin(m_pi * y[k]) * sin(m_pi * z[k]);
287  }
288  }
289  break;
290 
291  case eTestLinearSphere:
292  {
293  Array<OneD, NekDouble> temp(nq);
294 
295  NekDouble A = 2.0;
296  NekDouble B = 5.0;
297 
298  NekDouble m_a, m_b, m_c, m_d;
299  m_a = B - 1.0;
300  m_b = A * A;
301  m_c = -1.0 * B;
302  m_d = -1.0 * A * A;
303 
304  temp = Array<OneD, NekDouble>(nq, 0.0);
305  Vmath::Svtvp(nq, m_a, &inarray[0][0], 1, &temp[0], 1, &temp[0], 1);
306  Vmath::Svtvp(nq, m_b, &inarray[1][0], 1, &temp[0], 1,
307  &outarray[0][0], 1);
308 
309  temp = Array<OneD, NekDouble>(nq, 0.0);
310  Vmath::Svtvp(nq, m_c, &inarray[0][0], 1, &temp[0], 1, &temp[0], 1);
311  Vmath::Svtvp(nq, m_d, &inarray[1][0], 1, &temp[0], 1,
312  &outarray[1][0], 1);
313  }
314  break;
315 
317  {
318  NekDouble A = 2.0;
319  NekDouble B = 5.0;
320 
321  Array<OneD, NekDouble> Aonevec(nq, A);
322 
323  // cube = phys0*phys0*phy1
324  Array<OneD, NekDouble> cube(nq);
325  Vmath::Vmul(nq, &inarray[0][0], 1, &inarray[0][0], 1, &cube[0], 1);
326  Vmath::Vmul(nq, &inarray[1][0], 1, &cube[0], 1, &cube[0], 1);
327 
328  // outarray[0] = A - B*phy0 + phy0*phy0*phy1 - phy0
329  NekDouble coeff = -1.0 * B - 1.0;
330  Array<OneD, NekDouble> tmp(nq);
331  Vmath::Svtvp(nq, coeff, &inarray[0][0], 1, &cube[0], 1, &tmp[0], 1);
332  Vmath::Vadd(nq, &Aonevec[0], 1, &tmp[0], 1, &outarray[0][0], 1);
333 
334  // outarray[1] = B*phys0 - phy0*phy0*phy1
335  Vmath::Svtvm(nq, B, &inarray[0][0], 1, &cube[0], 1, &outarray[1][0],
336  1);
337  }
338  break;
339 
340  case eFHNStandard:
341  {
342  // \phi - \phi^3/3 - \psi
343  NekDouble a = 0.12;
344  NekDouble b = 0.011;
345  NekDouble c1 = 0.175;
346  NekDouble c2 = 0.03;
347  NekDouble d = 0.55;
348 
349  Array<OneD, NekDouble> tmp(nq);
350 
351  // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 v
352  Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
353  Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
354  Vmath::Vmul(nq, tmp, 1, inarray[0], 1, outarray[0], 1);
355  Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
356  Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
357 
358  Vmath::Smul(nq, -1.0 * c2, inarray[1], 1, tmp, 1);
359  Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
360 
361  // Reaction for \psi = b (\phi - d \psi )
362  Vmath::Svtvp(nq, -1.0 * d, inarray[1], 1, inarray[0], 1,
363  outarray[1], 1);
364  Vmath::Smul(nq, b, outarray[1], 1, outarray[1], 1);
365  }
366  break;
367 
368  case eFHNRogers:
369  {
370  NekDouble a = 0.13;
371  NekDouble b = 0.013;
372  NekDouble c1 = 0.26;
373  NekDouble c2 = 0.1;
374  NekDouble d = 1.0;
375 
376  Array<OneD, NekDouble> tmp(nq);
377 
378  // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 u v
379  Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
380  Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
381  Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
382  Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
383  Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
384 
385  Vmath::Vmul(nq, inarray[0], 1, inarray[1], 1, tmp, 1);
386  Vmath::Smul(nq, -1.0 * c2, tmp, 1, tmp, 1);
387  Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
388 
389  // Reaction for \psi = b (\phi - d \psi )
390  Vmath::Svtvp(nq, -1.0 * d, inarray[1], 1, inarray[0], 1,
391  outarray[1], 1);
392  Vmath::Smul(nq, b, outarray[1], 1, outarray[1], 1);
393  }
394  break;
395 
396  case eFHNAlievPanf:
397  {
398 
399  NekDouble a = 0.15;
400  NekDouble c1 = 8.0;
401  NekDouble c2 = 1.0;
402  NekDouble c0 = 0.002;
403  NekDouble mu1 = 0.2;
404  NekDouble mu2 = 0.3;
405 
406  Array<OneD, NekDouble> tmp(nq);
407 
408  // Reaction for \phi = c1 \phi ( \phi - a)*(1 - \phi) - c2 u v
409  Vmath::Smul(nq, -1.0 * c1, inarray[0], 1, outarray[0], 1);
410  Vmath::Sadd(nq, -1.0 * a, inarray[0], 1, tmp, 1);
411  Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
412  Vmath::Sadd(nq, -1.0, inarray[0], 1, tmp, 1);
413  Vmath::Vmul(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
414 
415  Vmath::Vmul(nq, inarray[0], 1, inarray[1], 1, tmp, 1);
416  Vmath::Smul(nq, -1.0 * c2, tmp, 1, tmp, 1);
417  Vmath::Vadd(nq, tmp, 1, outarray[0], 1, outarray[0], 1);
418 
419  // Reaction for \psi = (c0 + (\mu1 \psi/(\mu2+\phi) ) )*(-\psi - c1
420  // * \phi*(\phi - a - 1) )
421 
422  Vmath::Smul(nq, mu1, inarray[1], 1, outarray[1], 1);
423  Vmath::Sadd(nq, mu2, inarray[0], 1, tmp, 1);
424  Vmath::Vdiv(nq, outarray[1], 1, tmp, 1, outarray[1], 1);
425  Vmath::Sadd(nq, c0, outarray[1], 1, outarray[1], 1);
426 
427  Vmath::Sadd(nq, (-a - 1.0), inarray[0], 1, tmp, 1);
428  Vmath::Vmul(nq, inarray[0], 1, tmp, 1, tmp, 1);
429  Vmath::Smul(nq, c1, tmp, 1, tmp, 1);
430  Vmath::Vadd(nq, inarray[1], 1, tmp, 1, tmp, 1);
431  Vmath::Neg(nq, tmp, 1);
432 
433  Vmath::Vmul(nq, tmp, 1, outarray[1], 1, outarray[1], 1);
434  }
435  break;
436 
437  default:
438  break;
439  }
440 }
Array< OneD, NekDouble > m_epsilon
Definition: MMFDiffusion.h:152
SOLVER_UTILS_EXPORT int GetTotPoints()
@ eFHNStandard
Definition: MMFDiffusion.h:52
@ eTestLinearSphere
Definition: MMFDiffusion.h:50
@ eTestPlane
Definition: MMFDiffusion.h:48
@ eTestCube
Definition: MMFDiffusion.h:49
@ eTestNonlinearSphere
Definition: MMFDiffusion.h:51
@ eFHNRogers
Definition: MMFDiffusion.h:53
@ eFHNAlievPanf
Definition: MMFDiffusion.h:54
double NekDouble
void Vmul(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Multiply vector z = x*y.
Definition: Vmath.cpp:209
void 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.cpp:622
void Neg(int n, T *x, const int incx)
Negate x = -x.
Definition: Vmath.cpp:518
void Svtvm(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 minus vector): z = alpha*x - y
Definition: Vmath.cpp:664
void Vadd(int n, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
Add vector z = x+y.
Definition: Vmath.cpp:359
void 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.cpp:284
void Sadd(int n, const T alpha, const T *x, const int incx, T *y, const int incy)
Add scalar y = alpha + x.
Definition: Vmath.cpp:384

References Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_epsilon, Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::MMFSystem::m_pi, m_TestType, Vmath::Neg(), Vmath::Sadd(), Vmath::Smul(), Vmath::Svtvm(), Vmath::Svtvp(), Vmath::Vadd(), Vmath::Vdiv(), and Vmath::Vmul().

Referenced by v_InitObject().

◆ Morphogenesis()

void Nektar::MMFDiffusion::Morphogenesis ( const NekDouble  time,
unsigned int  field,
Array< OneD, NekDouble > &  outfield 
)
protected

Definition at line 566 of file MMFDiffusion.cpp.

568 {
569  int nq = GetTotPoints();
570 
571  int i, m, n, ind;
572  NekDouble a_n, d_n, gamma_n;
573  NekDouble A_mn, C_mn, theta, phi, radius;
574 
575  std::complex<double> Spericharmonic, delta_n, temp;
576  std::complex<double> varphi0, varphi1;
577  std::complex<double> B_mn, D_mn;
578 
579  // Set some parameter values
580  int Maxn = 6;
581  int Maxm = 2 * Maxn - 1;
582 
583  NekDouble A = 2.0;
584  NekDouble B = 5.0;
585 
586  NekDouble m_mu = 0.001;
587  NekDouble m_nu = 0.002;
588 
589  NekDouble m_a, m_b, m_c, m_d;
590 
591  m_a = B - 1.0;
592  m_b = A * A;
593  m_c = -1.0 * B;
594  m_d = -1.0 * A * A;
595 
598 
599  for (i = 0; i < Maxn; ++i)
600  {
601  Ainit[i] = Array<OneD, NekDouble>(Maxm, 0.0);
602  Binit[i] = Array<OneD, NekDouble>(Maxm, 0.0);
603  }
604 
605  Ainit[5][0] = -0.5839;
606  Ainit[5][1] = -0.8436;
607  Ainit[5][2] = -0.4764;
608  Ainit[5][3] = 0.6475;
609  Ainit[5][4] = 0.1886;
610  Ainit[5][5] = 0.8709;
611  Ainit[5][6] = -0.8338;
612  Ainit[5][7] = 0.1795;
613  Ainit[5][8] = -0.7873;
614  Ainit[5][9] = 0.8842;
615  Ainit[5][10] = 0.2943;
616 
617  Binit[5][0] = -0.6263;
618  Binit[5][1] = 0.9803;
619  Binit[5][2] = 0.7222;
620  Binit[5][3] = 0.5945;
621  Binit[5][4] = 0.6026;
622  Binit[5][5] = -0.2076;
623  Binit[5][6] = 0.4556;
624  Binit[5][7] = 0.6024;
625  Binit[5][8] = 0.9695;
626  Binit[5][9] = -0.4936;
627  Binit[5][10] = 0.1098;
628 
634 
635  m_fields[0]->GetCoords(x, y, z);
636  for (int i = 0; i < nq; ++i)
637  {
638  radius = sqrt(x[i] * x[i] + y[i] * y[i] + z[i] * z[i]);
639 
640  // theta is in [0, pi]
641  theta = asin(z[i] / radius) + 0.5 * m_pi;
642 
643  // phi is in [0, 2*pi]
644  phi = atan2(y[i], x[i]) + m_pi;
645 
646  varphi0 = 0.0 * varphi0;
647  varphi1 = 0.0 * varphi1;
648  for (n = 0; n < Maxn; ++n)
649  {
650  // Set up parameters
651  a_n = m_a - m_mu * (n * (n + 1) / radius / radius);
652  d_n = m_d - m_nu * (n * (n + 1) / radius / radius);
653 
654  gamma_n = 0.5 * (a_n + d_n);
655 
656  temp = (a_n + d_n) * (a_n + d_n) - 4.0 * (a_n * d_n - m_b * m_c);
657  delta_n = 0.5 * sqrt(temp);
658 
659  for (m = -n; m <= n; ++m)
660  {
661  ind = m + n;
662  A_mn = Ainit[n][ind];
663  C_mn = Binit[n][ind];
664 
665  B_mn = ((a_n - gamma_n) * Ainit[n][ind] + m_b * Binit[n][ind]) /
666  delta_n;
667  D_mn = (m_c * Ainit[n][ind] + (d_n - gamma_n) * Binit[n][ind]) /
668  delta_n;
669 
670  Spericharmonic =
671  boost::math::spherical_harmonic(n, m, theta, phi);
672  varphi0 += exp(gamma_n * time) *
673  (A_mn * cosh(delta_n * time) +
674  B_mn * sinh(delta_n * time)) *
675  Spericharmonic;
676  varphi1 += exp(gamma_n * time) *
677  (C_mn * cosh(delta_n * time) +
678  D_mn * sinh(delta_n * time)) *
679  Spericharmonic;
680  }
681  }
682 
683  u[i] = varphi0.real();
684  v[i] = varphi1.real();
685  }
686 
687  switch (field)
688  {
689  case 0:
690  {
691  outfield = u;
692  }
693  break;
694 
695  case 1:
696  {
697  outfield = v;
698  }
699  break;
700  }
701 }
NekDouble m_mu
scalarT< T > sqrt(scalarT< T > in)
Definition: scalar.hpp:294

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

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ PlanePhiWave()

Array< OneD, NekDouble > Nektar::MMFDiffusion::PlanePhiWave ( )
protected

Definition at line 703 of file MMFDiffusion.cpp.

704 {
705  int nq = GetTotPoints();
706  Array<OneD, NekDouble> outarray(nq, 0.0);
707 
711 
712  m_fields[0]->GetCoords(x, y, z);
713 
714  NekDouble xmin, ymin, xmax;
715 
716  xmin = Vmath::Vmin(nq, x, 1);
717  xmax = Vmath::Vmax(nq, x, 1);
718  ymin = Vmath::Vmin(nq, y, 1);
719 
720  NekDouble xp, yp, xp2;
721  for (int i = 0; i < nq; i++)
722  {
723  switch (m_InitWaveType)
724  {
725  case eLeft:
726  {
727  NekDouble radiusofinit = 4.0;
728  NekDouble frontstiff = 0.1;
729 
730  xp = x[i] - xmin;
731  outarray[i] =
732  1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff));
733  }
734  break;
735 
736  case eBothEnds:
737  {
738  NekDouble radiusofinit = 3.0;
739  NekDouble frontstiff = 0.1;
740 
741  xp = x[i] - xmin;
742  xp2 = x[i] - xmax;
743 
744  outarray[i] =
745  1.0 / (1.0 +
746  exp((sqrt(xp * xp) - radiusofinit) / frontstiff)) +
747  1.0 / (1.0 +
748  exp((sqrt(xp2 * xp2) - radiusofinit) / frontstiff));
749  }
750  break;
751 
752  case eCenter:
753  {
754  NekDouble radiusofinit = 6.0;
755  NekDouble frontstiff = 0.1;
756 
757  // NekDouble xc = 0.5*(Vmath::Vmax(nq, x, 1) + Vmath::Vmin(nq,
758  // x, 1));
759 
760  xp = x[i] - xmin;
761  outarray[i] =
762  1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff));
763  }
764  break;
765 
766  case eLeftBottomCorner:
767  {
768  NekDouble radiusofinit = 6.0;
769  NekDouble frontstiff = 0.1;
770  NekDouble bs = 2.0;
771 
772  xp = x[i] - xmin;
773  yp = y[i] - ymin;
774  outarray[i] =
775  1.0 /
776  (1.0 + exp((sqrt(xp * xp + yp * yp) / bs - radiusofinit) /
777  frontstiff));
778  }
779  break;
780 
781  case ePoint:
782  {
783  NekDouble xloc, yloc, zloc, rad;
784  NekDouble radiusofinit = 10.0;
785 
786  xloc = x[i] - m_InitPtx;
787  yloc = y[i] - m_InitPty;
788  zloc = z[i] - m_InitPtz;
789 
790  rad = sqrt(xloc * xloc + yloc * yloc + zloc * zloc);
791 
792  xloc = xloc / radiusofinit;
793  yloc = yloc / radiusofinit;
794  zloc = zloc / radiusofinit;
795 
796  if (rad < radiusofinit)
797  {
798  outarray[i] =
799  exp(-(1.0 / 2.0) *
800  (xloc * xloc + yloc * yloc + zloc * zloc));
801  }
802 
803  else
804  {
805  outarray[i] = 0.0;
806  }
807  }
808  break;
809 
810  case eSpiralDock:
811  {
812  NekDouble radiusofinit = 3.0;
813  NekDouble frontstiff = 0.1;
814  xp = x[i] - 4.0;
815  yp = y[i];
816  outarray[i] =
817  (1.0 / (1.0 + exp(2.0 * yp))) *
818  (1.0 / (1.0 + exp(-2.0 * xp))) *
819  (1.0 / (1.0 + exp((xp - radiusofinit) / frontstiff)));
820  }
821  break;
822 
823  default:
824  break;
825  }
826  }
827 
828  return outarray;
829 }
InitWaveType m_InitWaveType
Definition: MMFDiffusion.h:108
static NekDouble rad(NekDouble x, NekDouble y)
@ eLeftBottomCorner
Definition: MMFDiffusion.h:68
@ eBothEnds
Definition: MMFDiffusion.h:66
@ eSpiralDock
Definition: MMFDiffusion.h:70
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.cpp:1050
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.cpp:945

References Nektar::eBothEnds, Nektar::eCenter, Nektar::eLeft, Nektar::eLeftBottomCorner, Nektar::ePoint, Nektar::eSpiralDock, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, m_InitPtx, m_InitPty, m_InitPtz, m_InitWaveType, Nektar::LibUtilities::rad(), tinysimd::sqrt(), Vmath::Vmax(), and Vmath::Vmin().

Referenced by v_SetInitialConditions().

◆ TestCubeProblem()

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

Definition at line 546 of file MMFDiffusion.cpp.

549 {
550  int nq = GetTotPoints();
551 
555 
556  m_fields[0]->GetCoords(x, y, z);
557 
558  outfield = Array<OneD, NekDouble>(nq);
559  for (int k = 0; k < nq; k++)
560  {
561  outfield[k] = exp(-1.0 * m_pi * m_pi * time) * sin(m_pi * x[k]) *
562  sin(m_pi * y[k]) * sin(m_pi * z[k]);
563  }
564 }

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

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ TestPlaneProblem()

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

Definition at line 526 of file MMFDiffusion.cpp.

529 {
530  int nq = GetTotPoints();
531 
535 
536  m_fields[0]->GetCoords(x, y, z);
537 
538  outfield = Array<OneD, NekDouble>(nq);
539  for (int k = 0; k < nq; k++)
540  {
541  outfield[k] = exp(-1.0 * m_pi * m_pi * time) * sin(m_pi * x[k]) *
542  cos(m_pi * y[k]);
543  }
544 }

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

Referenced by v_EvaluateExactSolution(), and v_SetInitialConditions().

◆ v_EvaluateExactSolution()

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

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 831 of file MMFDiffusion.cpp.

834 {
835  switch (m_TestType)
836  {
837  case eTestPlane:
838  {
839  TestPlaneProblem(time, outfield);
840  }
841  break;
842 
843  case eTestCube:
844  {
845  TestCubeProblem(time, outfield);
846  }
847  break;
848 
849  case eTestLinearSphere:
851  {
852  Morphogenesis(time, field, outfield);
853  }
854  break;
855 
856  case eFHNStandard:
857  case eFHNRogers:
858  case eFHNAlievPanf:
859  {
860  int nq = GetTotPoints();
861  outfield = Array<OneD, NekDouble>(nq, 0.0);
862  }
863  /* Falls through. */
864  default:
865  {
866  EquationSystem::v_EvaluateExactSolution(field, outfield, time);
867  }
868  break;
869  }
870 }
void Morphogenesis(const NekDouble time, unsigned int field, Array< OneD, NekDouble > &outfield)
void TestCubeProblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
void TestPlaneProblem(const NekDouble time, Array< OneD, NekDouble > &outfield)
virtual SOLVER_UTILS_EXPORT void v_EvaluateExactSolution(unsigned int field, Array< OneD, NekDouble > &outfield, const NekDouble time)

References Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, Nektar::SolverUtils::EquationSystem::GetTotPoints(), m_TestType, Morphogenesis(), TestCubeProblem(), TestPlaneProblem(), and Nektar::SolverUtils::EquationSystem::v_EvaluateExactSolution().

◆ v_GenerateSummary()

void Nektar::MMFDiffusion::v_GenerateSummary ( SolverUtils::SummaryList s)
overrideprotectedvirtual

Prints a summary of the model parameters.

Reimplemented from Nektar::SolverUtils::MMFSystem.

Definition at line 872 of file MMFDiffusion.cpp.

873 {
876  SolverUtils::AddSummaryItem(s, "epsilon0", m_epsilon[0]);
877  SolverUtils::AddSummaryItem(s, "epsilon1", m_epsilon[1]);
878  SolverUtils::AddSummaryItem(s, "epsilon2", m_epsilon[2]);
880  {
881  SolverUtils::AddSummaryItem(s, "epsilon for u", m_epsu[0]);
882  SolverUtils::AddSummaryItem(s, "epsilon for v", m_epsu[1]);
883  }
884 }
virtual SOLVER_UTILS_EXPORT void v_GenerateSummary(SummaryList &s) override
Print a summary of time stepping parameters.
Definition: MMFSystem.cpp:2469
void AddSummaryItem(SummaryList &l, const std::string &name, const std::string &value)
Adds a summary item to the summary info list.
Definition: Misc.cpp:49
const char *const TestTypeMap[]
Definition: MMFDiffusion.h:58

References Nektar::SolverUtils::AddSummaryItem(), Nektar::eTestLinearSphere, m_epsilon, m_epsu, m_TestType, Nektar::TestTypeMap, and Nektar::SolverUtils::MMFSystem::v_GenerateSummary().

◆ v_InitObject()

void Nektar::MMFDiffusion::v_InitObject ( bool  DeclareField = true)
overrideprotectedvirtual

Init object for UnsteadySystem class.

Initialization object for UnsteadySystem class.

Reimplemented from Nektar::SolverUtils::UnsteadySystem.

Definition at line 63 of file MMFDiffusion.cpp.

64 {
65  UnsteadySystem::v_InitObject(DeclareFields);
66 
67  int nq = m_fields[0]->GetNpoints();
68  int nvar = m_fields.size();
69  int MFdim = 3;
70 
71  // Diffusivity coefficient for e^j
73  m_session->LoadParameter("epsilon0", m_epsilon[0], 1.0);
74  m_session->LoadParameter("epsilon1", m_epsilon[1], 1.0);
75  m_session->LoadParameter("epsilon2", m_epsilon[2], 1.0);
76 
77  // Diffusivity coefficient for u^j
78  m_epsu = Array<OneD, NekDouble>(nvar + 1);
79  m_session->LoadParameter("epsu0", m_epsu[0], 1.0);
80  m_session->LoadParameter("epsu1", m_epsu[1], 1.0);
81 
82  m_session->LoadParameter("InitPtx", m_InitPtx, 0.0);
83  m_session->LoadParameter("InitPty", m_InitPty, 0.0);
84  m_session->LoadParameter("InitPtz", m_InitPtz, 0.0);
85 
86  int shapedim = m_fields[0]->GetShapeDimension();
87  Array<OneD, Array<OneD, NekDouble>> Anisotropy(shapedim);
88  for (int j = 0; j < shapedim; ++j)
89  {
90  Anisotropy[j] = Array<OneD, NekDouble>(nq, 1.0);
91  Vmath::Fill(nq, sqrt(m_epsilon[j]), &Anisotropy[j][0], 1);
92  }
93 
94  MMFSystem::MMFInitObject(Anisotropy);
95 
96  // Define ProblemType
97  if (m_session->DefinesSolverInfo("TESTTYPE"))
98  {
99  std::string TestTypeStr = m_session->GetSolverInfo("TESTTYPE");
100  int i;
101  for (i = 0; i < (int)SIZE_TestType; ++i)
102  {
103  if (boost::iequals(TestTypeMap[i], TestTypeStr))
104  {
105  m_TestType = (TestType)i;
106  break;
107  }
108  }
109  }
110  else
111  {
112  m_TestType = (TestType)0;
113  }
114 
115  if (m_session->DefinesSolverInfo("INITWAVETYPE"))
116  {
117  std::string InitWaveTypeStr = m_session->GetSolverInfo("INITWAVETYPE");
118  for (int i = 0; i < (int)SIZE_TestType; ++i)
119  {
120  if (boost::iequals(InitWaveTypeMap[i], InitWaveTypeStr))
121  {
123  break;
124  }
125  }
126  }
127  else
128  {
130  }
131 
132  StdRegions::VarCoeffType MMFCoeffs[15] = {
141 
142  int indx;
143  Array<OneD, NekDouble> tmp(nq);
144  for (int k = 0; k < MFdim; ++k)
145  {
146  // For Moving Frames
147  indx = 5 * k;
148 
149  for (int j = 0; j < m_spacedim; ++j)
150  {
151  Vmath::Vcopy(nq, &m_movingframes[k][j * nq], 1, &tmp[0], 1);
152  m_varcoeff[MMFCoeffs[indx + j]] = tmp;
153  }
154 
155  // m_DivMF
156  Vmath::Vcopy(nq, &m_DivMF[k][0], 1, &tmp[0], 1);
157  m_varcoeff[MMFCoeffs[indx + 3]] = tmp;
158 
159  // \| e^k \|
160  tmp = Array<OneD, NekDouble>(nq, 0.0);
161  for (int i = 0; i < m_spacedim; ++i)
162  {
163  Vmath::Vvtvp(nq, &m_movingframes[k][i * nq], 1,
164  &m_movingframes[k][i * nq], 1, &tmp[0], 1, &tmp[0], 1);
165  }
166 
167  m_varcoeff[MMFCoeffs[indx + 4]] = tmp;
168  }
169 
170  if (!m_explicitDiffusion)
171  {
173  }
174 
176 }
void DefineOdeRhs(FuncPointerT func, ObjectPointerT obj)
void DefineImplicitSolve(FuncPointerT func, ObjectPointerT obj)
void DoOdeRhs(const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, const NekDouble time)
Computes the reaction terms and .
void DoImplicitSolve(const Array< OneD, const Array< OneD, NekDouble >> &inarray, Array< OneD, Array< OneD, NekDouble >> &outarray, NekDouble time, NekDouble lambda)
Solve for the diffusion term.
int m_spacedim
Spatial dimension (>= expansion dim).
LibUtilities::SessionReaderSharedPtr m_session
The session reader.
Array< OneD, Array< OneD, NekDouble > > m_DivMF
Definition: MMFSystem.h:195
Array< OneD, Array< OneD, NekDouble > > m_movingframes
Definition: MMFSystem.h:186
SOLVER_UTILS_EXPORT void MMFInitObject(const Array< OneD, const Array< OneD, NekDouble >> &Anisotropy, const int TangentXelem=-1)
Definition: MMFSystem.cpp:53
LibUtilities::TimeIntegrationSchemeOperators m_ode
The time integration scheme operators to use.
bool m_explicitDiffusion
Indicates if explicit or implicit treatment of diffusion is used.
virtual SOLVER_UTILS_EXPORT void v_InitObject(bool DeclareField=true) override
Init object for UnsteadySystem class.
const char *const InitWaveTypeMap[]
Definition: MMFDiffusion.h:74
@ SIZE_TestType
Length of enum list.
Definition: MMFDiffusion.h:55
void Vvtvp(int n, const T *w, const int incw, const T *x, const int incx, const T *y, const int incy, T *z, const int incz)
vvtvp (vector times vector plus vector): z = w*x + y
Definition: Vmath.cpp:574
void Fill(int n, const T alpha, T *x, const int incx)
Fill a vector with a constant value.
Definition: Vmath.cpp:45
void Vcopy(int n, const T *x, const int incx, T *y, const int incy)
Definition: Vmath.cpp:1255

References Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineImplicitSolve(), Nektar::LibUtilities::TimeIntegrationSchemeOperators::DefineOdeRhs(), DoImplicitSolve(), DoOdeRhs(), Nektar::StdRegions::eVarCoeffMF1Div, Nektar::StdRegions::eVarCoeffMF1Mag, Nektar::StdRegions::eVarCoeffMF1x, Nektar::StdRegions::eVarCoeffMF1y, Nektar::StdRegions::eVarCoeffMF1z, Nektar::StdRegions::eVarCoeffMF2Div, Nektar::StdRegions::eVarCoeffMF2Mag, Nektar::StdRegions::eVarCoeffMF2x, Nektar::StdRegions::eVarCoeffMF2y, Nektar::StdRegions::eVarCoeffMF2z, Nektar::StdRegions::eVarCoeffMF3Div, Nektar::StdRegions::eVarCoeffMF3Mag, Nektar::StdRegions::eVarCoeffMF3x, Nektar::StdRegions::eVarCoeffMF3y, Nektar::StdRegions::eVarCoeffMF3z, Vmath::Fill(), Nektar::InitWaveTypeMap, Nektar::SolverUtils::MMFSystem::m_DivMF, m_epsilon, m_epsu, Nektar::SolverUtils::UnsteadySystem::m_explicitDiffusion, Nektar::SolverUtils::EquationSystem::m_fields, m_InitPtx, m_InitPty, m_InitPtz, m_InitWaveType, Nektar::SolverUtils::MMFSystem::m_movingframes, Nektar::SolverUtils::UnsteadySystem::m_ode, Nektar::SolverUtils::EquationSystem::m_session, Nektar::SolverUtils::EquationSystem::m_spacedim, m_TestType, m_varcoeff, Nektar::SolverUtils::MMFSystem::MMFInitObject(), Nektar::SIZE_TestType, tinysimd::sqrt(), Nektar::TestTypeMap, Nektar::SolverUtils::UnsteadySystem::v_InitObject(), Vmath::Vcopy(), and Vmath::Vvtvp().

◆ v_SetInitialConditions()

void Nektar::MMFDiffusion::v_SetInitialConditions ( NekDouble  initialtime,
bool  dumpInitialConditions,
const int  domain 
)
overrideprotectedvirtual

Sets a custom initial condition.

Reimplemented from Nektar::SolverUtils::EquationSystem.

Definition at line 445 of file MMFDiffusion.cpp.

448 {
449  boost::ignore_unused(domain);
450 
451  int nq = GetTotPoints();
452 
453  switch (m_TestType)
454  {
455  case eTestPlane:
456  {
458 
459  TestPlaneProblem(initialtime, u);
460  m_fields[0]->SetPhys(u);
461  }
462  break;
463 
464  case eTestCube:
465  {
467 
468  TestCubeProblem(initialtime, u);
469  m_fields[0]->SetPhys(u);
470  /*for (int k=0; k<nq; ++k)
471  {
472  //for (int j=0; j<m_spacedim; ++j)
473  //{
474  cout << "_varcoeff" << u[k] <<endl;
475  // }
476  }*/
477  }
478  break;
479 
480  case eTestLinearSphere:
482  {
485 
486  Morphogenesis(initialtime, 0, u);
487  Morphogenesis(initialtime, 1, v);
488 
489  m_fields[0]->SetPhys(u);
490  m_fields[1]->SetPhys(v);
491  }
492  break;
493 
494  case eFHNStandard:
495  case eFHNRogers:
496  case eFHNAlievPanf:
497  {
498  Array<OneD, NekDouble> Zero(nq, 0.0);
499  m_fields[0]->SetPhys(PlanePhiWave());
500  m_fields[1]->SetPhys(Zero);
501  }
502  break;
503 
504  default:
505  {
506  EquationSystem::v_SetInitialConditions(initialtime, false);
507  }
508  break;
509  }
510 
511  // forward transform to fill the modal coeffs
512  for (int i = 0; i < m_fields.size(); ++i)
513  {
514  m_fields[i]->SetPhysState(true);
515  m_fields[i]->FwdTrans(m_fields[i]->GetPhys(),
516  m_fields[i]->UpdateCoeffs());
517  }
518 
519  if (dumpInitialConditions)
520  {
521  std::string outname = m_sessionName + "_initial.chk";
522  WriteFld(outname);
523  }
524 }
Array< OneD, NekDouble > PlanePhiWave()
virtual SOLVER_UTILS_EXPORT void v_SetInitialConditions(NekDouble initialtime=0.0, bool dumpInitialConditions=true, const int domain=0)
SOLVER_UTILS_EXPORT void WriteFld(const std::string &outname)
Write field data to the given filename.
std::string m_sessionName
Name of the session.
void Zero(int n, T *x, const int incx)
Zero vector.
Definition: Vmath.cpp:492

References Nektar::eFHNAlievPanf, Nektar::eFHNRogers, Nektar::eFHNStandard, Nektar::eTestCube, Nektar::eTestLinearSphere, Nektar::eTestNonlinearSphere, Nektar::eTestPlane, Nektar::SolverUtils::EquationSystem::GetTotPoints(), Nektar::SolverUtils::EquationSystem::m_fields, Nektar::SolverUtils::EquationSystem::m_sessionName, m_TestType, Morphogenesis(), PlanePhiWave(), TestCubeProblem(), TestPlaneProblem(), Nektar::SolverUtils::EquationSystem::v_SetInitialConditions(), Nektar::SolverUtils::EquationSystem::WriteFld(), and Vmath::Zero().

Friends And Related Function Documentation

◆ MemoryManager< MMFDiffusion >

friend class MemoryManager< MMFDiffusion >
friend

Definition at line 74 of file MMFDiffusion.h.

Member Data Documentation

◆ className

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

Name of class.

Definition at line 96 of file MMFDiffusion.h.

◆ m_epsilon

Array<OneD, NekDouble> Nektar::MMFDiffusion::m_epsilon
private

Definition at line 152 of file MMFDiffusion.h.

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

◆ m_epsu

Array<OneD, NekDouble> Nektar::MMFDiffusion::m_epsu
private

Definition at line 153 of file MMFDiffusion.h.

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

◆ m_InitPtx

NekDouble Nektar::MMFDiffusion::m_InitPtx
protected

Definition at line 146 of file MMFDiffusion.h.

Referenced by PlanePhiWave(), and v_InitObject().

◆ m_InitPty

NekDouble Nektar::MMFDiffusion::m_InitPty
protected

Definition at line 146 of file MMFDiffusion.h.

Referenced by PlanePhiWave(), and v_InitObject().

◆ m_InitPtz

NekDouble Nektar::MMFDiffusion::m_InitPtz
protected

Definition at line 146 of file MMFDiffusion.h.

Referenced by PlanePhiWave(), and v_InitObject().

◆ m_InitWaveType

InitWaveType Nektar::MMFDiffusion::m_InitWaveType
protected

Definition at line 108 of file MMFDiffusion.h.

Referenced by PlanePhiWave(), and v_InitObject().

◆ m_TestType

TestType Nektar::MMFDiffusion::m_TestType

◆ m_varcoeff

StdRegions::VarCoeffMap Nektar::MMFDiffusion::m_varcoeff
private

Variable diffusivity.

Definition at line 150 of file MMFDiffusion.h.

Referenced by DoImplicitSolve(), and v_InitObject().