15#ifndef dealii_integrators_advection_h
16#define dealii_integrators_advection_h
78 const double factor = 1.)
80 const unsigned int n_dofs = fe.dofs_per_cell;
82 const unsigned int n_components = fe.get_fe().n_components();
99 for (
unsigned k = 0;
k < fe.n_quadrature_points; ++
k)
101 const double dx = factor * fe.JxW(
k);
104 for (
unsigned j = 0;
j < n_dofs; ++
j)
105 for (
unsigned i = 0; i <
t_dofs; ++i)
106 for (
unsigned int c = 0; c < n_components; ++c)
110 for (
unsigned int d = 1; d < dim; ++d)
113 M(i,
j) -= dx *
wgradv * fe.shape_value_component(
j,
k, c);
136 const unsigned int nq = fe.n_quadrature_points;
137 const unsigned int n_dofs = fe.dofs_per_cell;
149 for (
unsigned k = 0;
k <
nq; ++
k)
151 const double dx = factor * fe.JxW(
k);
152 for (
unsigned i = 0; i < n_dofs; ++i)
153 for (
unsigned int d = 0; d < dim; ++d)
154 result(i) += dx * input[
k][d] * fe.shape_value(i,
k) *
179 const unsigned int nq = fe.n_quadrature_points;
180 const unsigned int n_dofs = fe.dofs_per_cell;
181 const unsigned int n_comp = fe.get_fe().n_components();
194 for (
unsigned k = 0;
k <
nq; ++
k)
196 const double dx = factor * fe.JxW(
k);
197 for (
unsigned i = 0; i < n_dofs; ++i)
198 for (
unsigned int c = 0; c <
n_comp; ++c)
199 for (
unsigned int d = 0; d < dim; ++d)
200 result(i) += dx * input[c][
k][d] *
201 fe.shape_value_component(i,
k, c) *
217 const std::vector<double> &input,
221 const unsigned int nq = fe.n_quadrature_points;
222 const unsigned int n_dofs = fe.dofs_per_cell;
234 for (
unsigned k = 0;
k <
nq; ++
k)
236 const double dx = factor * fe.JxW(
k);
237 for (
unsigned i = 0; i < n_dofs; ++i)
238 for (
unsigned int d = 0; d < dim; ++d)
239 result(i) -= dx * input[
k] * fe.shape_grad(i,
k)[d] *
257 const ArrayView<
const std::vector<double>> &input,
261 const unsigned int nq = fe.n_quadrature_points;
262 const unsigned int n_dofs = fe.dofs_per_cell;
263 const unsigned int n_comp = fe.get_fe().n_components();
276 for (
unsigned k = 0;
k <
nq; ++
k)
278 const double dx = factor * fe.JxW(
k);
279 for (
unsigned i = 0; i < n_dofs; ++i)
280 for (
unsigned int c = 0; c <
n_comp; ++c)
281 for (
unsigned int d = 0; d < dim; ++d)
282 result(i) -= dx * input[c][
k] *
283 fe.shape_grad_component(i,
k, c)[d] *
315 const unsigned int n_dofs = fe.dofs_per_cell;
317 unsigned int n_components = fe.get_fe().n_components();
328 for (
unsigned k = 0;
k < fe.n_quadrature_points; ++
k)
330 const double dx = factor * fe.JxW(
k);
333 for (
unsigned int d = 0; d < dim; ++d)
338 for (
unsigned i = 0; i <
t_dofs; ++i)
339 for (
unsigned j = 0;
j < n_dofs; ++
j)
341 if (fe.get_fe().is_primitive())
343 dx *
nv * fe.shape_value(i,
k) * fe.shape_value(
j,
k);
345 for (
unsigned int c = 0; c < n_components; ++c)
347 fetest.shape_value_component(i,
k, c) *
348 fe.shape_value_component(
j,
k, c);
384 const std::vector<double> &input,
385 const std::vector<double> &
data,
389 const unsigned int n_dofs = fe.dofs_per_cell;
402 for (
unsigned k = 0;
k < fe.n_quadrature_points; ++
k)
404 const double dx = factor * fe.JxW(
k);
407 for (
unsigned int d = 0; d < dim; ++d)
411 const double val = (
nv > 0.) ? input[
k] : -
data[
k];
413 for (
unsigned i = 0; i < n_dofs; ++i)
415 const double v = fe.shape_value(i,
k);
451 const ArrayView<
const std::vector<double>> &input,
456 const unsigned int n_dofs = fe.dofs_per_cell;
457 const unsigned int n_comp = fe.get_fe().n_components();
470 for (
unsigned k = 0;
k < fe.n_quadrature_points; ++
k)
472 const double dx = factor * fe.JxW(
k);
475 for (
unsigned int d = 0; d < dim; ++d)
478 std::vector<double> val(
n_comp);
480 for (
unsigned int d = 0; d <
n_comp; ++d)
482 val[d] = (
nv > 0.) ? input[d][
k] : -
data[d][
k];
483 for (
unsigned i = 0; i < n_dofs; ++i)
485 const double v = fe.shape_value_component(i,
k, d);
525 const double factor = 1.)
527 const unsigned int n1 =
fe1.dofs_per_cell;
539 for (
unsigned k = 0;
k <
fe1.n_quadrature_points; ++
k)
542 for (
unsigned int d = 1; d < dim; ++d)
550 for (
unsigned i = 0; i <
n1; ++i)
551 for (
unsigned j = 0;
j <
n1; ++
j)
553 if (
fe1.get_fe().is_primitive())
562 for (
unsigned int d = 0; d <
fe1.get_fe().n_components();
566 fe.shape_value_component(
j,
k, d) *
567 fetest.shape_value_component(i,
k, d);
569 fe.shape_value_component(
j,
k, d) *
570 fetestn.shape_value_component(i,
k, d);
605 const std::vector<double> &
input1,
606 const std::vector<double> &
input2,
608 const double factor = 1.)
615 const unsigned int n1 =
fe1.dofs_per_cell;
627 for (
unsigned k = 0;
k <
fe1.n_quadrature_points; ++
k)
630 for (
unsigned int d = 1; d < dim; ++d)
634 for (
unsigned i = 0; i <
n1; ++i)
636 const double v1 =
fe1.shape_value(i,
k);
637 const double v2 =
fe2.shape_value(i,
k);
685 const double factor = 1.)
687 const unsigned int n_comp =
fe1.get_fe().n_components();
688 const unsigned int n1 =
fe1.dofs_per_cell;
703 for (
unsigned k = 0;
k <
fe1.n_quadrature_points; ++
k)
706 for (
unsigned int d = 1; d < dim; ++d)
710 for (
unsigned i = 0; i <
n1; ++i)
711 for (
unsigned int d = 0; d <
n_comp; ++d)
713 const double v1 =
fe1.shape_value_component(i,
k, d);
714 const double v2 =
fe2.shape_value_component(i,
k, d);
#define DEAL_II_DEPRECATED
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
#define Assert(cond, exc)
#define AssertVectorVectorDimension(VEC, DIM1, DIM2)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
std::vector< index_type > data
void upwind_face_residual(Vector< double > &result1, Vector< double > &result2, const FEValuesBase< dim > &fe1, const FEValuesBase< dim > &fe2, const std::vector< double > &input1, const std::vector< double > &input2, const ArrayView< const std::vector< double > > &velocity, const double factor=1.)
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const ArrayView< const std::vector< double > > &velocity, const double factor=1.)
void cell_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< Tensor< 1, dim > > &input, const ArrayView< const std::vector< double > > &velocity, double factor=1.)
void upwind_value_residual(Vector< double > &result, const FEValuesBase< dim > &fe, const std::vector< double > &input, const std::vector< double > &data, const ArrayView< const std::vector< double > > &velocity, double factor=1.)
void upwind_value_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const ArrayView< const std::vector< double > > &velocity, double factor=1.)
Library of integrals over cells and faces.
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)