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

#include <AdvectionWeakDG.h>

Inheritance diagram for Nektar::SolverUtils::AdvectionWeakDG:
Inheritance graph
[legend]
Collaboration diagram for Nektar::SolverUtils::AdvectionWeakDG:
Collaboration graph
[legend]

Static Public Member Functions

static AdvectionSharedPtr create (std::string advType)

Static Public Attributes

static std::string type

Protected Member Functions

 AdvectionWeakDG ()
virtual void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 Initialise AdvectionWeakDG objects and store them before starting the time-stepping.
virtual void v_Advect (const int nConvective, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &advVel, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 Compute the advection term at each time-step using the Discontinuous Glaerkin approach (DG).
- Protected Member Functions inherited from Nektar::SolverUtils::Advection
virtual SOLVER_UTILS_EXPORT void v_SetBaseFlow (const Array< OneD, Array< OneD, NekDouble > > &inarray)
 Overrides the base flow used during linearised advection.

Additional Inherited Members

- Public Member Functions inherited from Nektar::SolverUtils::Advection
SOLVER_UTILS_EXPORT void InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 Interface function to initialise the advection object.
SOLVER_UTILS_EXPORT void Advect (const int nConvectiveFields, const Array< OneD, MultiRegions::ExpListSharedPtr > &fields, const Array< OneD, Array< OneD, NekDouble > > &advVel, const Array< OneD, Array< OneD, NekDouble > > &inarray, Array< OneD, Array< OneD, NekDouble > > &outarray, const NekDouble &time)
 Interface function to advect the vector field.
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxVector (FuncPointerT func, ObjectPointerT obj)
 Set the flux vector callback function.
void SetRiemannSolver (RiemannSolverSharedPtr riemann)
 Set a Riemann solver object for this advection object.
void SetFluxVector (AdvectionFluxVecCB fluxVector)
 Set the flux vector callback function.
void SetBaseFlow (const Array< OneD, Array< OneD, NekDouble > > &inarray)
 Set the base flow used for linearised advection objects.
- Protected Attributes inherited from Nektar::SolverUtils::Advection
AdvectionFluxVecCB m_fluxVector
 Callback function to the flux vector (set when advection is in conservative form).
RiemannSolverSharedPtr m_riemann
 Riemann solver for DG-type schemes.
int m_spaceDim
 Storage for space dimension. Used for homogeneous extension.

Detailed Description

Definition at line 45 of file AdvectionWeakDG.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::AdvectionWeakDG::AdvectionWeakDG ( )
protected

Definition at line 47 of file AdvectionWeakDG.cpp.

Referenced by create().

{
}

Member Function Documentation

static AdvectionSharedPtr Nektar::SolverUtils::AdvectionWeakDG::create ( std::string  advType)
inlinestatic

Definition at line 48 of file AdvectionWeakDG.h.

References AdvectionWeakDG().

{
}
void Nektar::SolverUtils::AdvectionWeakDG::v_Advect ( const int  nConvectiveFields,
const Array< OneD, MultiRegions::ExpListSharedPtr > &  fields,
const Array< OneD, Array< OneD, NekDouble > > &  advVel,
const Array< OneD, Array< OneD, NekDouble > > &  inarray,
Array< OneD, Array< OneD, NekDouble > > &  outarray,
const NekDouble time 
)
protectedvirtual

Compute the advection term at each time-step using the Discontinuous Glaerkin approach (DG).

Parameters
nConvectiveFieldsNumber of fields.
fieldsPointer to fields.
advVelAdvection velocities.
inarraySolution at the previous time-step.
outarrayAdvection term to be passed at the time integration class.

Implements Nektar::SolverUtils::Advection.

Definition at line 76 of file AdvectionWeakDG.cpp.

References ASSERTL1, Nektar::SolverUtils::Advection::m_fluxVector, Nektar::SolverUtils::Advection::m_riemann, Nektar::SolverUtils::Advection::m_spaceDim, Vmath::Neg(), and Vmath::Vadd().

{
int nDim = fields[0]->GetCoordim(0);
int nPointsTot = fields[0]->GetTotPoints();
int nCoeffs = fields[0]->GetNcoeffs();
int nTracePointsTot = fields[0]->GetTrace()->GetTotPoints();
int i, j;
Array<OneD, Array<OneD, NekDouble> > tmp(nConvectiveFields);
Array<OneD, Array<OneD, Array<OneD, NekDouble> > > fluxvector(
nConvectiveFields);
// Allocate storage for flux vector F(u).
for (i = 0; i < nConvectiveFields; ++i)
{
fluxvector[i] =
Array<OneD, Array<OneD, NekDouble> >(m_spaceDim);
for (j = 0; j < m_spaceDim; ++j)
{
fluxvector[i][j] = Array<OneD, NekDouble>(nPointsTot);
}
}
"Riemann solver must be provided for AdvectionWeakDG.");
m_fluxVector(inarray, fluxvector);
// Get the advection part (without numerical flux)
for(i = 0; i < nConvectiveFields; ++i)
{
tmp[i] = Array<OneD, NekDouble>(nCoeffs, 0.0);
for (j = 0; j < nDim; ++j)
{
fields[i]->IProductWRTDerivBase(j, fluxvector[i][j],
outarray[i]);
Vmath::Vadd(nCoeffs, outarray[i], 1, tmp[i], 1, tmp[i], 1);
}
}
// Store forwards/backwards space along trace space
Array<OneD, Array<OneD, NekDouble> > Fwd (nConvectiveFields);
Array<OneD, Array<OneD, NekDouble> > Bwd (nConvectiveFields);
Array<OneD, Array<OneD, NekDouble> > numflux(nConvectiveFields);
for(i = 0; i < nConvectiveFields; ++i)
{
Fwd[i] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
Bwd[i] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
numflux[i] = Array<OneD, NekDouble>(nTracePointsTot, 0.0);
fields[i]->GetFwdBwdTracePhys(inarray[i], Fwd[i], Bwd[i]);
}
m_riemann->Solve(m_spaceDim, Fwd, Bwd, numflux);
// Evaulate <\phi, \hat{F}\cdot n> - OutField[i]
for(i = 0; i < nConvectiveFields; ++i)
{
Vmath::Neg (nCoeffs, tmp[i], 1);
fields[i]->AddTraceIntegral (numflux[i], tmp[i]);
fields[i]->MultiplyByElmtInvMass(tmp[i], tmp[i]);
fields[i]->BwdTrans (tmp[i], outarray[i]);
}
}
void Nektar::SolverUtils::AdvectionWeakDG::v_InitObject ( LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields 
)
protectedvirtual

Initialise AdvectionWeakDG objects and store them before starting the time-stepping.

Parameters
pSessionPointer to session reader.
pFieldsPointer to fields.

Reimplemented from Nektar::SolverUtils::Advection.

Definition at line 58 of file AdvectionWeakDG.cpp.

{
Advection::v_InitObject(pSession, pFields);
}

Member Data Documentation

std::string Nektar::SolverUtils::AdvectionWeakDG::type
static
Initial value:
RegisterCreatorFunction("WeakDG", AdvectionWeakDG::create)

Definition at line 53 of file AdvectionWeakDG.h.