static
void compute_flux_matrix (const InputVector &W,
std::array <std::array <Number, dim>, EulerEquations<dim>::n_components > &flux)
-// Number (&flux)[n_components][dim])
{
// First compute the pressure that appears in the flux matrix, and then
// compute the first <code>dim</code> columns of the matrix that
// $\alpha$. It's form has also been given already in the introduction:
template <typename InputVector, typename Number>
static
- void numerical_normal_flux (const Point<dim> &normal,
- const InputVector &Wplus,
- const InputVector &Wminus,
- const double alpha,
+ void numerical_normal_flux (const Point<dim> &normal,
+ const InputVector &Wplus,
+ const InputVector &Wminus,
+ const double alpha,
std::array < Number, n_components> &normal_flux)
-// Number (&normal_flux)[n_components])
{
-// Number iflux[n_components][dim];
-// Number oflux[n_components][dim];
std::array <std::array <Number, dim>, EulerEquations<dim>::n_components > iflux, oflux;
compute_flux_matrix (Wplus, iflux);
template <typename InputVector, typename Number>
static
void compute_forcing_vector (const InputVector &W,
- std::array < Number, n_components> &forcing)
-// Number (&forcing)[n_components])
+ std::array < Number, n_components> &forcing)
{
const double gravity = -1.0;
// that we compute the flux matrices and right hand sides in terms of
// autodifferentiation variables, so that the Jacobian contributions can
// later easily be computed from it:
-// typedef Sacado::Fad::DFad<double> FluxMatrix[EulerEquations<dim>::n_components][dim];
-// FluxMatrix *flux = new FluxMatrix[n_q_points];
std::vector <
- std::array <std::array <Sacado::Fad::DFad<double>, dim>, EulerEquations<dim>::n_components >
-// std::array <std::array <Sacado::Fad::DFad<double>, EulerEquations<dim>::n_components>, dim >
- > flux(n_q_points);
-
-
-// typedef double FluxMatrixOld[EulerEquations<dim>::n_components][dim];
-// FluxMatrixOld *flux_old = new FluxMatrixOld[n_q_points];
+ std::array <std::array <Sacado::Fad::DFad<double>, dim>, EulerEquations<dim>::n_components >
+ > flux(n_q_points);
std::vector <
- std::array <std::array <double, dim>, EulerEquations<dim>::n_components >
- > flux_old(n_q_points);
-
-// typedef Sacado::Fad::DFad<double> ForcingVector[EulerEquations<dim>::n_components];
-// ForcingVector *forcing = new ForcingVector[n_q_points];
+ std::array <std::array <double, dim>, EulerEquations<dim>::n_components >
+ > flux_old(n_q_points);
std::vector < std::array< Sacado::Fad::DFad<double>, EulerEquations<dim>::n_components> > forcing(n_q_points);
-// typedef double ForcingVectorOld[EulerEquations<dim>::n_components];
-// ForcingVectorOld *forcing_old = new ForcingVectorOld[n_q_points];
std::vector < std::array< double, EulerEquations<dim>::n_components> > forcing_old(n_q_points);
for (unsigned int q=0; q<n_q_points; ++q)
system_matrix.add(dof_indices[i], dof_indices, residual_derivatives);
right_hand_side(dof_indices[i]) -= R_i.val();
}
-
-// delete[] forcing;
-// delete[] flux;
-// delete[] forcing_old;
-// delete[] flux_old;
-
}
// w^-, \mathbf n)$ for each quadrature point. Before calling the function
// that does so, we also need to determine the Lax-Friedrich's stability
// parameter:
-// typedef Sacado::Fad::DFad<double> NormalFlux[EulerEquations<dim>::n_components];
-// NormalFlux *normal_fluxes; = new NormalFlux[n_q_points];
std::vector< std::array < Sacado::Fad::DFad<double>, EulerEquations<dim>::n_components> > normal_fluxes(n_q_points);
-
-// typedef double NormalFluxOld[EulerEquations<dim>::n_components];
-// NormalFluxOld *normal_fluxes_old = new NormalFluxOld[n_q_points];
-
std::vector< std::array < double, EulerEquations<dim>::n_components> > normal_fluxes_old(n_q_points);
double alpha;
right_hand_side(dof_indices[i]) -= R_i.val();
}
-
-// delete[] normal_fluxes;
-// delete[] normal_fluxes_old;
}