#include <dgcartesiansyseqnacousticperturb.h>

Inheritance diagram for DgSysEqnAcousticPerturb< nDim >:

Collaboration diagram for DgSysEqnAcousticPerturb< nDim >:

Classes
struct	CV

Public Member Functions
	DgSysEqnAcousticPerturb (MInt solverId)
	Constructor calls parent constructor & loads all necessary properties for this equation. More...

void	calcInitialCondition (const MFloat t, const MFloat x, MFloat const nodeVars, MFloat *u) const
	Calculate the initial condition for a certain point in space. More...

void	calcFlux (const MFloat const nodeVars, const MFloat const q, const MInt noNodes1D, MFloat *const flux) const
	Calculates the physical fluxes in all dimensions for all integrations points within a cell. More...

void	calcSource (const MFloat const nodeVars, const MFloat const u, const MInt noNodes1D, const MFloat t, const MFloat const x, MFloat const src) const
	Calculates the source terms for all integration points within a cell. More...

void	calcSpongeSource (const MFloat nodeVars, const MFloat u, const MInt noNodes1D, const MFloat eta, MFloat src) const
	Calculates the sponge source terms for all integration points within a cell. More...

MFloat	getTimeStep (const MFloat nodeVars, const MFloat const u, const MInt noNodes1D, const MFloat invJacobian, const MInt sbpMode) const
	Calculate the time step for an explicit time stepping scheme for a given element. More...

void	calcRiemann (const MFloat nodeVarsL, const MFloat nodeVarsR, const MFloat stateL, const MFloat stateR, const MInt noNodes1D, const MInt dirId, MFloat *flux) const
	Calculates the numerical flux at a surface given two states (left and right). More...

void	calcRiemannLaxFriedich (const MFloat nodeVarsL, const MFloat nodeVarsR, const MFloat stateL, const MFloat stateR, const MInt noNodes1D, const MInt dirId, MFloat *flux) const
	Calculate Riemann flux using the local Lax-Friedrichs flux scheme. More...

void	calcRiemannRoe (const MFloat nodeVarsL, const MFloat nodeVarsR, const MFloat stateL, const MFloat stateR, const MInt noNodes1D, const MInt dirId, MFloat *flux) const

void	primToCons (const MFloat prim, MFloat cons) const
	Calculates a set of conservative variables from a set of primitive variables. More...

void	consToPrim (const MFloat cons, MFloat prim) const
	Calculates a set of primitive variables from a set of conservative variables. More...

void	getDefaultNodeVars (MFloat *const nodeVars) const
	Return the default node variables. More...

void	getDefaultNodeVarsBody (MFloat *const nodeVars) const
	Return the default node variables inside a body. More...

MBool	extendNodeVar (const MInt varId) const
	Return if the given node variable should be extended. More...

void	calcRiemannRoe (const MFloat nodeVarsL, const MFloat nodeVarsR, const MFloat stateL, const MFloat stateR, const MInt noNodes1D, const MInt dirId, MFloat *flux) const

void	calcRiemannRoe (const MFloat nodeVarsL, const MFloat nodeVarsR, const MFloat stateL, const MFloat stateR, const MInt noNodes1D, const MInt dirId, MFloat *flux) const

Public Member Functions inherited from DgSysEqn< nDim, DgSysEqnAcousticPerturb< nDim > >
MFloat	cfl () const

MFloat	cflScaled (const MInt polyDeg) const

Private Member Functions
void	calcFlux1D (const MFloat const nodeVars, const MFloat const q, const MInt noNodes1D, const MInt dirId, MFloat *const flux) const
	Calculates the physical fluxes in dimension `dirId` for all integrations points on an element face. More...

MFloat	cflFactor (const MInt polyDeg) const

const MString	s_consVarNames []

const MString	s_primVarNames []

const MString	s_nodeVarNames []

const MString	s_consVarNames []

const MString	s_primVarNames []

const MString	s_nodeVarNames []

Private Attributes
std::array< MFloat, nDim >	m_meanVelocity

MFloat	m_meanDensity

MFloat	m_meanSpeedOfSound

MBool	m_constantSpeedOfSound = false

MFloat	m_spongePressureInfy

MFloat	m_spongeSigma

MBool	m_compressibleSourceTerm = true

MString	m_dgIntegrationMethod

MInt	m_dgTimeIntegrationScheme

std::array< std::array< MFloat, s_maxPolyDeg+1 >, 2 >	m_cflFactor

Static Private Attributes
static const MString	s_sysEqnName = "DG_SYSEQN_ACOUSTICPERTURB"

static const MInt	s_noVariables = nDim + 1

static const MString	s_consVarNames [s_noVariables]

static const MString	s_primVarNames [s_noVariables]

static const MInt	s_noNodeVars = 2 * nDim + 2

static const MBool	s_hasTimeDependentNodeVars = false

static const MString	s_nodeVarNames [s_noNodeVars]

static const MInt	s_maxPolyDeg = 31

Friends
class	DgSysEqn< nDim, DgSysEqnAcousticPerturb >

Additional Inherited Members
Static Public Member Functions inherited from DgSysEqn< nDim, DgSysEqnAcousticPerturb< nDim > >
static constexpr MInt	noVars ()

static constexpr MInt	noNodeVars ()

static constexpr MBool	hasTimeDependentNodeVars ()

static const MString &	sysEqnName ()

static const MString &	consVarNames (MInt i)

static const MString &	primVarNames (MInt i)

static const MString &	nodeVarNames (MInt i)

Public Attributes inherited from DgSysEqn< nDim, DgSysEqnAcousticPerturb< nDim > >
MInt	m_initialCondition

MFloat	m_initialNumberWaves

MInt	m_sourceTerm

MInt	m_riemannSolver

Protected Member Functions inherited from DgSysEqn< nDim, DgSysEqnAcousticPerturb< nDim > >
	DgSysEqn (MInt solverId)

Protected Attributes inherited from DgSysEqn< nDim, DgSysEqnAcousticPerturb< nDim > >
const MInt	m_solverId

Detailed Description

template<MInt nDim>
class DgSysEqnAcousticPerturb< nDim >

Definition at line 27 of file dgcartesiansyseqnacousticperturb.h.

Constructor & Destructor Documentation

◆ DgSysEqnAcousticPerturb()

template<MInt nDim>

DgSysEqnAcousticPerturb< nDim >::DgSysEqnAcousticPerturb ( MInt solverId )

explicit

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in] solverId Current solver id.

Definition at line 164 of file dgcartesiansyseqnacousticperturb.h.

Member Function Documentation

◆ calcFlux()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcFlux	(	const MFloat *const	nodeVars,
		const MFloat *const	q,
		const MInt	noNodes1D,
		MFloat *const	flux
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	q	Pointer to solution variables (size: noNodes1D^nDim*noVars)).
[in]	noNodes1D	Number of nodes 1D in the cell.
[out]	flux	Calculated flux (size: noNodes1D^nDimnoVarsnDim)).

Definition at line 827 of file dgcartesiansyseqnacousticperturb.h.

                                                                       {
  // TRACE();
 
  const MInt noNodes1D3 = (nDim == 3) ? noNodes1D : 1;
  const MFloatTensor U(const_cast<MFloat*>(q), noNodes1D, noNodes1D, noNodes1D3, s_noVariables);
  const MFloatTensor U0(const_cast<MFloat*>(nodeVars), noNodes1D, noNodes1D, noNodes1D3, s_noNodeVars);
  MFloatTensor f(flux, noNodes1D, noNodes1D, noNodes1D3, nDim, s_noVariables);
 
  // The following flux equations are calculated here (if 2D, consider all
  // z-components to be zero; for an incompressible case: cBar=1, rhoBar=1):
 
  // Flux in x-direction:
  // f_x(CV[u]) = u*uBar + v*vBar + w*wBar + p/rhoBar
  // f_x(CV[v]) = 0
  // f_x(CV[w]) = 0
  // f_x(CV[p]) = rhoBar*cBar^2*u + uBar*p
 
  // Flux in y-direction:
  // f_y(CV[u]) = 0
  // f_y(CV[v]) = u*uBar + v*vBar + w*wBar + p/rhoBar
  // f_y(CV[w]) = 0
  // f_y(CV[p]) = rhoBar*cBar^2*v + vBar*p
 
  // Flux in z-direction:
  // f_z(CV[u]) = 0
  // f_z(CV[v]) = 0
  // f_z(CV[w]) = u*uBar + v*vBar + w*wBar + p/rhoBar
  // f_z(CV[p]) = rhoBar*cBar^2*w + wBar*p
 
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D; j++) {
      for(MInt k = 0; k < noNodes1D3; k++) {
        for(MInt d = 0; d < nDim; d++) {
          // Momentum fluxes
          for(MInt fd = 0; fd < nDim; fd++) {
            f(i, j, k, d, CV::UU[fd]) = 0.0;
          }
          f(i, j, k, d, CV::UU[d]) = std::inner_product(&U0(i, j, k, CV::UU0[0]),
                                                        &U0(i, j, k, CV::UU0[0]) + nDim,
                                                        &U(i, j, k, CV::UU[0]),
                                                        U(i, j, k, CV::P) / U0(i, j, k, CV::RHO0));
 
          // Energy flux
          f(i, j, k, d, CV::P) = U0(i, j, k, CV::RHO0) * POW2(U0(i, j, k, CV::C0)) * U(i, j, k, CV::UU[d])
                                 + U0(i, j, k, CV::UU0[d]) * U(i, j, k, CV::P);
        }
      }
    }
  }
}

◆ calcFlux1D()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcFlux1D	(	const MFloat *const	nodeVars,
		const MFloat *const	q,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *const	flux
	)		const

private

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	q	Pointer to solution variables (size: noNodes1D^(nDim-1)*noVars)).
[in]	noNodes1D	Number of nodes 1D in the cell.
[in]	dirId	Orientation of the face (0 - x, 1 - y, 2 - z).
[out]	flux	Calculated flux (size: noNodes1D^(nDim-1)*noVars)).

Definition at line 896 of file dgcartesiansyseqnacousticperturb.h.

                                                                         {
  // TRACE();
 
  const MInt noNodes1D3 = (nDim == 3) ? noNodes1D : 1;
  const MFloatTensor U(const_cast<MFloat*>(q), noNodes1D, noNodes1D3, s_noVariables);
  const MFloatTensor U0(const_cast<MFloat*>(nodeVars), noNodes1D, noNodes1D3, s_noNodeVars);
  MFloatTensor f(flux, noNodes1D, noNodes1D3, s_noVariables);
 
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D3; j++) {
      // Momentum fluxes
      for(MInt d = 0; d < nDim; d++) {
        if(d == dirId) {
          f(i, j, CV::UU[d]) = std::inner_product(&U0(i, j, CV::UU0[0]),
                                                  &U0(i, j, CV::UU0[0]) + nDim,
                                                  &U(i, j, CV::UU[0]),
                                                  U(i, j, CV::P) / U0(i, j, CV::RHO0));
        } else {
          f(i, j, CV::UU[d]) = 0.0;
        }
      }
 
      // Energy flux
      f(i, j, CV::P) = U0(i, j, CV::RHO0) * POW2(U0(i, j, CV::C0)) * U(i, j, CV::UU[dirId])
                       + U0(i, j, CV::UU0[dirId]) * U(i, j, CV::P);
    }
  }
}

◆ calcInitialCondition()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcInitialCondition	(	const MFloat	t,
		const MFloat *	x,
		MFloat *const	nodeVars,
		MFloat *	u
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2012-11-28

Parameters

[in]	t	Time at which the initial condition should be applied.
[in]	x	A pointer to the node coordinates.
[out]	u	A pointer to the solution storage.

Note: If you (for some reason) use nodeVariables information to calculate the initial condition for the conservative variables, you also must set the nodeVariables here, otherwise the calcErrorNorms() method of the DG solver procudes faulty results.

Definition at line 397 of file dgcartesiansyseqnacousticperturb.h.

◆ calcRiemann()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcRiemann	(	const MFloat *	nodeVarsL,
		const MFloat *	nodeVarsR,
		const MFloat *	stateL,
		const MFloat *	stateR,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *	flux
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	nodeVarsL	Surface node variables on the left side of the surface (i.e. in -x/-y/-z-direction).
[in]	nodeVarsR	Surface node variables on the right side of the surface (i.e. in +x/+y/+z-direction).
[in]	stateL	Surface variables on the left side of the surface (i.e. in -x/-y/-z-direction).
[in]	stateR	Surface variables on the right side of the surface (i.e. in +x/+y/+z-direction).
[in]	noNodes1D	Number of nodes 1D.
[in]	dirId	Direction (0 = x-direction, 1 = y-direction, 2 = z-direction).
[out]	flux	Calculated Riemann flux (size: noNodes1D^(nDim-1)*noVars).

Definition at line 1423 of file dgcartesiansyseqnacousticperturb.h.

                                                                    {
  // TRACE();
 
  // Solve Riemann problem
  switch(this->m_riemannSolver) {
    case 0: // local Lax-Friedrichs flux
    {
      calcRiemannLaxFriedich(nodeVarsL, nodeVarsR, stateL, stateR, noNodes1D, dirId, flux);
      break;
    }
 
    case 1: // Roe's flux
    {
      TERMM(1, "Implementation does not work for internal fluxes, see comments");
      calcRiemannRoe(nodeVarsL, nodeVarsR, stateL, stateR, noNodes1D, dirId, flux);
      break;
    }
 
    default:
      mTerm(1, AT_, "The specified Riemann solver is not valid!");
      break;
  }
}

◆ calcRiemannLaxFriedich()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcRiemannLaxFriedich	(	const MFloat *	nodeVarsL,
		const MFloat *	nodeVarsR,
		const MFloat *	stateL,
		const MFloat *	stateR,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *	flux
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2015-03-13

Note: for argument description, see calcRieman(...).

Definition at line 1461 of file dgcartesiansyseqnacousticperturb.h.

                                                                               {
  // TRACE();
 
  const MInt noNodes1D3 = (nDim == 3) ? noNodes1D : 1;
  const MFloatTensor uL(const_cast<MFloat*>(stateL), noNodes1D, noNodes1D3, s_noVariables);
  const MFloatTensor uR(const_cast<MFloat*>(stateR), noNodes1D, noNodes1D3, s_noVariables);
 
  const MFloatTensor nL(const_cast<MFloat*>(nodeVarsL), noNodes1D, noNodes1D3, s_noNodeVars);
  const MFloatTensor nR(const_cast<MFloat*>(nodeVarsR), noNodes1D, noNodes1D3, s_noNodeVars);
 
#ifndef _OPENMP
  MFloatScratchSpace maxLambda(noNodes1D, noNodes1D3, AT_, "maxLambda");
#else
  MFloatTensor maxLambda(noNodes1D, noNodes1D3);
#endif
 
  // Calculate maximum eigenvalue
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D3; j++) {
      maxLambda(i, j) = std::max(fabs(nL(i, j, CV::UU0[dirId])) + nL(i, j, CV::C0),
                                 fabs(nR(i, j, CV::UU0[dirId])) + nR(i, j, CV::C0));
    }
  }
 
#ifndef _OPENMP
  MFloatScratchSpace fluxL(noNodes1D, noNodes1D3, s_noVariables, AT_, "fluxL");
  MFloatScratchSpace fluxR(noNodes1D, noNodes1D3, s_noVariables, AT_, "fluxR");
#else
  MFloatTensor fluxL(noNodes1D, noNodes1D3, s_noVariables);
  MFloatTensor fluxR(noNodes1D, noNodes1D3, s_noVariables);
#endif
 
  // Calculate flux from left and right state
  calcFlux1D(nodeVarsL, stateL, noNodes1D, dirId, &fluxL[0]);
  calcFlux1D(nodeVarsR, stateR, noNodes1D, dirId, &fluxR[0]);
 
  // Solve Riemann problem
  MFloatTensor riemann(flux, noNodes1D, noNodes1D3, s_noVariables);
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D3; j++) {
      for(MInt n = 0; n < s_noVariables; n++) {
        riemann(i, j, n) = 0.5 * ((fluxL(i, j, n) + fluxR(i, j, n)) - maxLambda(i, j) * (uR(i, j, n) - uL(i, j, n)));
      }
    }
  }
}

◆ calcRiemannRoe() [1/3]

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcRiemannRoe	(	const MFloat *	nodeVarsL,
		const MFloat *	nodeVarsR,
		const MFloat *	stateL,
		const MFloat *	stateR,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *	flux
	)		const

◆ calcRiemannRoe() [2/3]

void DgSysEqnAcousticPerturb< 2 >::calcRiemannRoe	(	const MFloat *	nodeVarsL,
		const MFloat *	nodeVarsR,
		const MFloat *	stateL,
		const MFloat *	stateR,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *	flux
	)		const

inline

Definition at line 1532 of file dgcartesiansyseqnacousticperturb.h.

                                                                           {
  // TRACE();
 
  const MInt nDim = 2;
  const MFloatTensor uL(const_cast<MFloat*>(stateL), noNodes1D, noVars());
  const MFloatTensor uR(const_cast<MFloat*>(stateR), noNodes1D, noVars());
  const MFloatTensor nL(const_cast<MFloat*>(nodeVarsL), noNodes1D, noNodeVars());
  const MFloatTensor nR(const_cast<MFloat*>(nodeVarsR), noNodes1D, noNodeVars());
  MFloatTensor f(flux, noNodes1D, s_noVariables);
 
  // TODO labels:DG use node variables instead
  // Define mean flow constants
  const MFloat c = 1.0;
  const MFloat rho = 1.0;
 
  // Normal direction
  const std::array<MInt, nDim> n = {{1 - dirId, dirId}};
 
  // Characteristic boundary conditions
  for(MInt i = 0; i < noNodes1D; i++) {
    // FIXME labels:DG,totest The following is just a hack and probably does not work when the
    // Roe Riemann solver is used for internal fluxes
    // Determine mean flow
    const std::array<MFloat, nDim> meanVelocity = {
        {0.5 * (nL(i, CV::U0) + nR(i, CV::U0)), 0.5 * (nL(i, CV::V0) + nR(i, CV::V0))}};
 
    // Mean flow matrix with simplified transformation matrix
    const std::array<MFloat, nDim> um = {
        {meanVelocity[0] * n[0] + meanVelocity[1] * n[1], -meanVelocity[0] * n[1] + meanVelocity[1] * n[0]}};
 
    // stateL(vecL) and stateR(vecR) variables multiplied by transformation
    // matrix (see reference)
    const std::array<MFloat, nDim> vecL = {
        {uL(i, CV::UU[0]) * n[0] + uL(i, CV::UU[1]) * n[1], -uL(i, CV::UU[0]) * n[1] + uL(i, CV::UU[1]) * n[0]}};
    const std::array<MFloat, nDim> vecR = {
        {uR(i, CV::UU[0]) * n[0] + uR(i, CV::UU[1]) * n[1], -uR(i, CV::UU[0]) * n[1] + uR(i, CV::UU[1]) * n[0]}};
 
    // Calculate characteristic waves
    // L1 and L2 composite from vecL and vecR
    const MFloat L1 = 0.5 * (uR(i, CV::P) + rho * c * (-vecR[0] + um[1] / (um[0] - c) * (-vecR[1])));
    const MFloat L2 = 0.5 * (uL(i, CV::P) + rho * c * (vecL[0] + um[1] / (um[0] + c) * vecL[1]));
 
    // Calculate Riemann flux
    f(i, CV::P) = (um[0] - c) * (L1) + (um[0] + c) * (L2);
    f(i, CV::UU[0]) = 1 / rho / c * (-1 * ((um[0] - c) * L1) + ((um[0] + c) * L2)) * n[0];
    f(i, CV::UU[1]) = 1 / rho / c * (-1 * ((um[0] - c) * L1) + ((um[0] + c) * L2)) * n[1];
  }
}

◆ calcRiemannRoe() [3/3]

void DgSysEqnAcousticPerturb< 3 >::calcRiemannRoe	(	const MFloat *	nodeVarsL,
		const MFloat *	nodeVarsR,
		const MFloat *	stateL,
		const MFloat *	stateR,
		const MInt	noNodes1D,
		const MInt	dirId,
		MFloat *	flux
	)		const

inline

Definition at line 1606 of file dgcartesiansyseqnacousticperturb.h.

                                                                           {
  // TRACE();
 
  const MInt nDim = 3;
  const MFloatTensor uL(const_cast<MFloat*>(stateL), noNodes1D, noNodes1D, noVars());
  const MFloatTensor uR(const_cast<MFloat*>(stateR), noNodes1D, noNodes1D, noVars());
  const MFloatTensor nL(const_cast<MFloat*>(nodeVarsL), noNodes1D, noNodes1D, noNodeVars());
  const MFloatTensor nR(const_cast<MFloat*>(nodeVarsR), noNodes1D, noNodes1D, noNodeVars());
  MFloatTensor f(flux, noNodes1D, noNodes1D, s_noVariables);
 
  // TODO labels:DG use node variables instead
  // Define mean flow constants
  const MFloat c = 1.0;
  const MFloat rho = 1.0;
 
  // Normal direction
  const std::array<MInt, nDim> n = {{dirId == 0 ? 1 : 0, dirId == 1 ? 1 : 0, dirId == 2 ? 1 : 0}};
 
  // Characteristic boundary conditions
  // stateL and stateR multiplied by transformation matrix (see reference)
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D; j++) {
      // FIXME labels:DG,totest The following is just a hack and probably does not work when the
      // Roe Riemann solver is used for internal fluxes
      // Determine mean flow
      const std::array<MFloat, nDim> meanVelocity = {{0.5 * (nL(i, j, CV::U0) + nR(i, j, CV::U0)),
                                                      0.5 * (nL(i, j, CV::V0) + nR(i, j, CV::V0)),
                                                      0.5 * (nL(i, j, CV::W0) + nR(i, j, CV::W0))}};
 
      // Mean flow matrix with transformation matrix (see ref.)
      const std::array<MFloat, nDim> um = {{meanVelocity[0] * n[0] + meanVelocity[1] * n[1] + meanVelocity[2] * n[2],
                                            meanVelocity[0] * n[1] + meanVelocity[1] * n[2] + meanVelocity[2] * n[0],
                                            meanVelocity[0] * n[2] + meanVelocity[1] * n[0] + meanVelocity[2] * n[1]}};
 
      const std::array<MFloat, nDim> vecL = {
          {uL(i, j, CV::UU[0]) * n[0] + uL(i, j, CV::UU[1]) * n[1] + uL(i, j, CV::UU[2]) * n[2],
           uL(i, j, CV::UU[0]) * n[1] + uL(i, j, CV::UU[1]) * n[2] + uL(i, j, CV::UU[2]) * n[0],
           uL(i, j, CV::UU[0]) * n[2] + uL(i, j, CV::UU[1]) * n[0] + uL(i, j, CV::UU[2]) * n[1]}};
 
      const std::array<MFloat, nDim> vecR = {
          {uR(i, j, CV::UU[0]) * n[0] + uR(i, j, CV::UU[1]) * n[1] + uR(i, j, CV::UU[2]) * n[2],
           uR(i, j, CV::UU[0]) * n[1] + uR(i, j, CV::UU[1]) * n[2] + uR(i, j, CV::UU[2]) * n[0],
           uR(i, j, CV::UU[0]) * n[2] + uR(i, j, CV::UU[1]) * n[0] + uR(i, j, CV::UU[2]) * n[1]}};
 
      // Calculate characteristic waves
      // L1 and L2 composite from vecL and vecR
      const MFloat L1 =
          0.5
          * (uR(i, j, CV::P)
             + rho * c * (-vecL[0] + um[1] / (um[0] - c) * (-vecL[1]) + um[2] / (um[0] - c) * (-vecL[2])));
      const MFloat L2 =
          0.5
          * (uL(i, j, CV::P) + rho * c * (vecR[0] + um[1] / (um[0] + c) * vecR[1] + um[2] / (um[0] + c) * (vecR[2])));
      // Calculate Riemann flux
      f(i, j, CV::P) = (um[0] - c) * (L1) + (um[0] + c) * (L2);
      f(i, j, CV::UU[0]) = 1 / rho / c * (-1 * ((um[0] - c) * L1) + ((um[0] + c) * L2)) * n[0];
      f(i, j, CV::UU[1]) = 1 / rho / c * (-1 * ((um[0] - c) * L1) + ((um[0] + c) * L2)) * n[1];
      f(i, j, CV::UU[2]) = 1 / rho / c * (-1 * ((um[0] - c) * L1) + ((um[0] + c) * L2)) * n[2];
    }
  }
}

◆ calcSource()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcSource	(	const MFloat *const	nodeVars,
		const MFloat *const	u,
		const MInt	noNodes1D,
		const MFloat	t,
		const MFloat *const	x,
		MFloat *const	src
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	u	Pointer to solution variables (size: noNodes1D^nDim*noVars)).
[in]	noNodes1D	in the cell.
[in]	t	Current time.
[in]	x	Coordinates of the integration points of the current cell.
[out]	src	Pointer to where source terms are stored (size: noNodes1D^nDim*noVars)).

Definition at line 978 of file dgcartesiansyseqnacousticperturb.h.

◆ calcSpongeSource()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::calcSpongeSource	(	const MFloat *	nodeVars,
		const MFloat *	u,
		const MInt	noNodes1D,
		const MFloat *	eta,
		MFloat *	src
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	nodeVars	Pointer to node variables (size: noNodes1D^nDim*noNodeVars)
[in]	u	Pointer to solution variables (size: noNodes1D^nDim*noVars)).
[in]	noNodes1D	in the cell.
[in]	eta	Pointer to eta used for the sponge calculations.
[out]	src	Pointer to where source terms are stored (size: noNodes1D^nDim*noVars)).

Definition at line 1328 of file dgcartesiansyseqnacousticperturb.h.

                                                                                                                 {
  const MInt noNodes1D3 = (nDim == 3) ? noNodes1D : 1;
 
  // gamma = heat capacity ratio for monoatomic gas
  MFloat gammaMinusOne = 5.0 / 3.0 - 1.0;
  MFloat FgammaMinusOne = 1.0 / gammaMinusOne;
 
  const MFloatTensor spongeEta(const_cast<MFloat*>(eta), noNodes1D, noNodes1D, noNodes1D3);
  MFloatTensor uState(const_cast<MFloat*>(u), noNodes1D, noNodes1D, noNodes1D3, s_noVariables);
  MFloatTensor spongeSource(src, noNodes1D, noNodes1D, noNodes1D3, s_noVariables);
 
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D; j++) {
      for(MInt k = 0; k < noNodes1D3; k++) {
        MFloat deltaP = (m_spongePressureInfy - uState(i, j, k, CV::P)) * FgammaMinusOne;
 
        spongeSource(i, j, k, CV::P) = m_spongeSigma * spongeEta(i, j, k) * deltaP;
        for(MInt l = 0; l < nDim; l++) {
          spongeSource(i, j, k, CV::UU[l]) = 0.0;
        }
      }
    }
  }
}

◆ cflFactor()

template<MInt nDim>

MFloat DgSysEqnAcousticPerturb< nDim >::cflFactor ( const MInt polyDeg ) const

private

Yield an empirically-derived maximum stable CFL number. The user should be able to specify a CFL number of "1.0" in virtually all cases.

Author: Rodrigo Miguez (rodrigo) rodri.nosp@m.go.m.nosp@m.iguez.nosp@m.@rwt.nosp@m.h-aac.nosp@m.hen..nosp@m.de

Date: 2018-02-28

Parameters

[in]	polyDeg	Polynomial degree in the cell.
[out]	cfl	factor obtained from m_cflFactor.

Definition at line 942 of file dgcartesiansyseqnacousticperturb.h.

                                                                        {
  if(polyDeg > s_maxPolyDeg) {
    TERMM(1, "Polynomial degree exceeds maximum supported");
  }
  // Note: The factor "0.95" is used to avoid numerical instabilities when using the full
  //       cfl factor value.
 
  const MFloat factor = 0.95 * m_cflFactor[nDim - 2][polyDeg];
 
  return factor;
}

◆ consToPrim()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::consToPrim	(	const MFloat *	cons,
		MFloat *	prim
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	cons	Conservative variable storage (size: noVars).
[out]	prim	Primitive variable storage (size: noVars).

Definition at line 1686 of file dgcartesiansyseqnacousticperturb.h.

                                                                                     {
  // Very easy...
  std::copy(cons, cons + s_noVariables, prim);
}

◆ extendNodeVar()

template<MInt nDim>

MBool DgSysEqnAcousticPerturb< nDim >::extendNodeVar ( const MInt varId ) const

Definition at line 1744 of file dgcartesiansyseqnacousticperturb.h.

                                                                         {
  // Do not extend...
  if(varId == CV::C0) return false; // mean speed of sound ...
  for(MInt i = 0; i < nDim; i++) {
    if(varId == CV::DC0[i]) { // and its derivatives
      return false;
    }
  }
  return true;
}

◆ getDefaultNodeVars()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::getDefaultNodeVars ( MFloat *const nodeVars ) const

Author: Ansgar Niemoeller (ansgar) a.nie.nosp@m.moel.nosp@m.ler@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2016-07-30

Parameters

[out] nodeVars Contains the default node variables ordered by CV.

Definition at line 1716 of file dgcartesiansyseqnacousticperturb.h.

                                                                                   {
  // Mean velocities
  std::copy_n(std::begin(m_meanVelocity), nDim, &nodeVars[CV::UU0[0]]);
  // Mean density
  nodeVars[CV::RHO0] = m_meanDensity;
  // Mean speed of sound
  nodeVars[CV::C0] = m_meanSpeedOfSound;
  // Derivatives of mean speed of sound, always zero
  std::fill_n(&nodeVars[CV::DC0[0]], nDim, 0.0);
}

◆ getDefaultNodeVarsBody()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::getDefaultNodeVarsBody ( MFloat *const nodeVars ) const

Definition at line 1730 of file dgcartesiansyseqnacousticperturb.h.

                                                                                       {
  // Mean velocities
  std::fill_n(&nodeVars[CV::UU0[0]], nDim, 0.0);
  // Mean density
  nodeVars[CV::RHO0] = m_meanDensity;
  // Mean speed of sound
  nodeVars[CV::C0] = m_meanSpeedOfSound;
  // Derivatives of mean speed of sound, always zero
  std::fill_n(&nodeVars[CV::DC0[0]], nDim, 0.0);
}

◆ getTimeStep()

template<MInt nDim>

MFloat DgSysEqnAcousticPerturb< nDim >::getTimeStep	(	const MFloat *	nodeVars,
		const MFloat *const	u,
		const MInt	noNodes1D,
		const MFloat	invJacobian,
		const MInt	sbpMode
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	u	Pointer to element variables.
[in]	noNodes1D	Current number of nodes 1D.
[in]	invJacobian	Inverse Jacobian of element.

Returns: Calculated time step for the element.

Definition at line 1368 of file dgcartesiansyseqnacousticperturb.h.

                                                                            {
  // TRACE();
 
  const MInt noNodes1D3 = (nDim == 3) ? noNodes1D : 1;
  const MFloat inf = std::numeric_limits<MFloat>::infinity();
  const MFloatTensor U0(const_cast<MFloat*>(nodeVars), noNodes1D, noNodes1D, noNodes1D3, s_noNodeVars);
 
  MFloat maxLambda[] = {-inf, -inf, -inf};
  for(MInt i = 0; i < noNodes1D; i++) {
    for(MInt j = 0; j < noNodes1D; j++) {
      for(MInt k = 0; k < noNodes1D3; k++) {
        for(MInt d = 0; d < nDim; d++) {
          maxLambda[d] = std::max(maxLambda[d], fabs(U0(i, j, k, CV::UU0[d])) + U0(i, j, k, CV::C0));
        }
      }
    }
  }
 
  MFloat dt;
  if(sbpMode) {
    dt = this->cfl() * F2 / (invJacobian * std::accumulate(&maxLambda[0], &maxLambda[0] + nDim, F0) * (noNodes1D - 1));
  } else {
    const MInt polyDeg = noNodes1D - 1;
    dt = this->cflScaled(polyDeg) * cflFactor(polyDeg) * F2
         / (invJacobian * std::accumulate(&maxLambda[0], &maxLambda[0] + nDim, F0));
  }
  return dt;
}

◆ primToCons()

template<MInt nDim>

void DgSysEqnAcousticPerturb< nDim >::primToCons	(	const MFloat *	prim,
		MFloat *	cons
	)		const

Author: Michael Schlottke (mic) mic@a.nosp@m.ia.r.nosp@m.wth-a.nosp@m.ache.nosp@m.n.de

Date: 2013-04-18

Parameters

[in]	prim	Primitive variable storage (size: noVars).
[out]	cons	Conservative variable storage (size: noVars).

Definition at line 1703 of file dgcartesiansyseqnacousticperturb.h.

                                                                                     {
  // Very easy...
  std::copy(prim, prim + s_noVariables, cons);
}