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::AdvectionNonConservative Class Reference

#include <AdvectionNonConservative.h>

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

Static Public Member Functions

static AdvectionSharedPtr create (std::string advType)
 

Static Public Attributes

static std::string type
 

Protected Member Functions

 AdvectionNonConservative ()
 
virtual void v_InitObject (LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
 Initialise AdvectionNonConservative objects and store them before starting the time-stepping. More...
 
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)
 Advects a vector field. More...
 
- 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. More...
 

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. More...
 
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. More...
 
template<typename FuncPointerT , typename ObjectPointerT >
void SetFluxVector (FuncPointerT func, ObjectPointerT obj)
 Set the flux vector callback function. More...
 
void SetRiemannSolver (RiemannSolverSharedPtr riemann)
 Set a Riemann solver object for this advection object. More...
 
void SetFluxVector (AdvectionFluxVecCB fluxVector)
 Set the flux vector callback function. More...
 
void SetBaseFlow (const Array< OneD, Array< OneD, NekDouble > > &inarray)
 Set the base flow used for linearised advection objects. More...
 
- Protected Attributes inherited from Nektar::SolverUtils::Advection
AdvectionFluxVecCB m_fluxVector
 Callback function to the flux vector (set when advection is in conservative form). More...
 
RiemannSolverSharedPtr m_riemann
 Riemann solver for DG-type schemes. More...
 
int m_spaceDim
 Storage for space dimension. Used for homogeneous extension. More...
 

Detailed Description

Definition at line 45 of file AdvectionNonConservative.h.

Constructor & Destructor Documentation

Nektar::SolverUtils::AdvectionNonConservative::AdvectionNonConservative ( )
protected

Definition at line 46 of file AdvectionNonConservative.cpp.

Referenced by create().

47  {
48 
49  }

Member Function Documentation

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

Definition at line 48 of file AdvectionNonConservative.h.

References AdvectionNonConservative().

49  {
51  }
boost::shared_ptr< Advection > AdvectionSharedPtr
A shared pointer to an Advection object.
Definition: Advection.h:158
void Nektar::SolverUtils::AdvectionNonConservative::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

Advects a vector field.

Implements Nektar::SolverUtils::Advection.

Definition at line 65 of file AdvectionNonConservative.cpp.

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

72  {
73  int nDim = advVel.num_elements();
74  int nPointsTot = fields[0]->GetNpoints();
75  Array<OneD, NekDouble> grad0,grad1,grad2;
76 
77  grad0 = Array<OneD, NekDouble> (nPointsTot);
78 
79  if (nDim > 1)
80  {
81  grad1 = Array<OneD,NekDouble>(nPointsTot);
82  }
83 
84  if (nDim > 2)
85  {
86  grad2 = Array<OneD,NekDouble>(nPointsTot);
87  }
88 
89 
90  for (int i = 0; i < nConvectiveFields; ++i)
91  {
92  // Evaluate V \cdot Grad(u)
93  switch(nDim)
94  {
95  case 1:
96  fields[0]->PhysDeriv(inarray[i], grad0);
97 
98  Vmath::Vmul(nPointsTot,
99  grad0, 1,
100  advVel[0], 1,
101  outarray[i], 1);
102  break;
103  case 2:
104  fields[0]->PhysDeriv(inarray[i], grad0, grad1);
105 
106 
107  // Calculate advection terms
108  Vmath::Vmul (nPointsTot,
109  grad0, 1,
110  advVel[0], 1,
111  outarray[i], 1);
112 
113  Vmath::Vvtvp(nPointsTot,
114  grad1, 1,
115  advVel[1], 1,
116  outarray[i], 1,
117  outarray[i], 1);
118 
119  break;
120  case 3:
121  fields[0]->PhysDeriv(inarray[i], grad0, grad1, grad2);
122 
123  // Calculate advection terms
124  Vmath::Vmul (nPointsTot,
125  grad0, 1,
126  advVel[0], 1,
127  outarray[i], 1);
128 
129  Vmath::Vvtvp(nPointsTot,
130  grad1, 1,
131  advVel[1], 1,
132  outarray[i], 1,
133  outarray[i], 1);
134 
135  Vmath::Vvtvp(nPointsTot,
136  grad2, 1,
137  advVel[2], 1,
138  outarray[i], 1,
139  outarray[i], 1);
140  break;
141  default:
142  ASSERTL0(false,"dimension unknown");
143  }
144  }
145  }
#define ASSERTL0(condition, msg)
Definition: ErrorUtil.hpp:188
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:428
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:169
void Nektar::SolverUtils::AdvectionNonConservative::v_InitObject ( LibUtilities::SessionReaderSharedPtr  pSession,
Array< OneD, MultiRegions::ExpListSharedPtr pFields 
)
protectedvirtual

Initialise AdvectionNonConservative 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 AdvectionNonConservative.cpp.

References Nektar::SolverUtils::Advection::v_InitObject().

61  {
62  Advection::v_InitObject(pSession, pFields);
63  }
virtual SOLVER_UTILS_EXPORT void v_InitObject(LibUtilities::SessionReaderSharedPtr pSession, Array< OneD, MultiRegions::ExpListSharedPtr > pFields)
Initialises the advection object.
Definition: Advection.cpp:97

Member Data Documentation

std::string Nektar::SolverUtils::AdvectionNonConservative::type
static
Initial value:
RegisterCreatorFunction("NonConservative",

Definition at line 53 of file AdvectionNonConservative.h.