17 #ifndef __deal2__integrators_advection_h
18 #define __deal2__integrators_advection_h
21 #include <deal.II/base/config.h>
23 #include <deal.II/base/quadrature.h>
24 #include <deal.II/lac/full_matrix.h>
25 #include <deal.II/fe/mapping.h>
26 #include <deal.II/fe/fe_values.h>
27 #include <deal.II/meshworker/dof_info.h>
81 const VectorSlice<
const std::vector<std::vector<double> > > &velocity,
82 const double factor = 1.)
86 const unsigned int n_components = fe.
get_fe().n_components();
93 const unsigned int v_increment = (velocity[0].size() == 1) ? 0 : 1;
105 const double dx = factor * fe.
JxW(k);
106 const unsigned int vindex = k * v_increment;
108 for (
unsigned j=0; j<n_dofs; ++j)
109 for (
unsigned i=0; i<t_dofs; ++i)
110 for (
unsigned int c=0; c<n_components; ++c)
112 double wgradv = velocity[0][vindex]
114 for (
unsigned int d=1; d<dim; ++d)
115 wgradv += velocity[d][vindex]
135 const VectorSlice<
const std::vector<std::vector<double> > > &velocity,
144 const unsigned int v_increment = (velocity[0].size() == 1) ? 0 : 1;
145 if (v_increment == 1)
150 for (
unsigned k=0; k<nq; ++k)
152 const double dx = factor * fe.
JxW(k);
153 for (
unsigned i=0; i<n_dofs; ++i)
154 for (
unsigned int d=0; d<dim; ++d)
155 result(i) += dx * input[k][d]
156 * fe.
shape_value(i,k) * velocity[d][k * v_increment];
175 const VectorSlice<
const std::vector<std::vector<double> > > &velocity,
180 const unsigned int n_comp = fe.
get_fe().n_components();
186 const unsigned int v_increment = (velocity[0].size() == 1) ? 0 : 1;
187 if (v_increment == 1)
192 for (
unsigned k=0; k<nq; ++k)
194 const double dx = factor * fe.
JxW(k);
195 for (
unsigned i=0; i<n_dofs; ++i)
196 for (
unsigned int c=0; c<n_comp; ++c)
197 for (
unsigned int d=0; d<dim; ++d)
198 result(i) += dx * input[c][k][d]
237 const VectorSlice<
const std::vector<std::vector<double> > > &velocity,
242 unsigned int n_components = fe.
get_fe().n_components();
247 const unsigned int v_increment = (velocity[0].size() == 1) ? 0 : 1;
248 if (v_increment == 1)
255 const double dx = factor * fe.
JxW(k);
258 for (
unsigned int d=0; d<dim; ++d)
263 for (
unsigned i=0; i<t_dofs; ++i)
264 for (
unsigned j=0; j<n_dofs; ++j)
266 if (fe.
get_fe().is_primitive())
269 for (
unsigned int c=0; c<n_components; ++c)
319 const VectorSlice<
const std::vector<std::vector<double> > > &velocity,
320 const double factor = 1.)
328 const unsigned int v_increment = (velocity[0].size() == 1) ? 0 : 1;
329 if (v_increment == 1)
336 double nbeta = fe1.
normal_vector(k)[0] * velocity[0][k * v_increment];
337 for (
unsigned int d=1; d<dim; ++d)
338 nbeta += fe1.
normal_vector(k)[d] * velocity[d][k * v_increment];
339 const double dx_nbeta = factor * nbeta * fe1.
JxW(k);
341 for (
unsigned i=0; i<n1; ++i)
342 for (
unsigned j=0; j<n1; ++j)
343 if (fe1.
get_fe().is_primitive())
366 for (
unsigned int d=0; d<fe1.
get_fe().n_components(); ++d)
392 DEAL_II_NAMESPACE_CLOSE
double shape_value_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
#define AssertDimension(dim1, dim2)
Tensor< 1, spacedim > shape_grad_component(const unsigned int function_no, const unsigned int point_no, const unsigned int component) const
const unsigned int dofs_per_cell
const Point< spacedim > & normal_vector(const unsigned int i) const
#define AssertVectorVectorDimension(vec, dim1, dim2)
void cell_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< Tensor< 1, dim > > &input, const VectorSlice< const std::vector< std::vector< double > > > &velocity, double factor=1.)
Library of integrals over cells and faces.
#define Assert(cond, exc)
const unsigned int n_quadrature_points
const double & shape_value(const unsigned int function_no, const unsigned int point_no) const
void upwind_value_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const VectorSlice< const std::vector< std::vector< double > > > &velocity, double factor=1.)
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const VectorSlice< const std::vector< std::vector< double > > > &velocity, const double factor=1.)
const FiniteElement< dim, spacedim > & get_fe() const
::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
double JxW(const unsigned int quadrature_point) const