18 #include <deal.II/base/utilities.h>
19 #include <deal.II/base/memory_consumption.h>
20 #include <deal.II/base/polynomial.h>
21 #include <deal.II/base/tensor_product_polynomials.h>
22 #include <deal.II/base/polynomials_piecewise.h>
23 #include <deal.II/fe/fe_poly.h>
25 #include <deal.II/matrix_free/shape_info.h>
33 namespace MatrixFreeFunctions
38 template <
typename Number>
47 template <
typename Number>
54 ExcMessage(
"FEEvaluation only works for scalar finite elements."));
56 const unsigned int n_dofs_1d = fe.
degree+1,
57 n_q_points_1d = quad.
size();
72 lexicographic = fe_poly != 0 ?
84 n_q_points = Utilities::fixed_power<dim>(n_q_points_1d);
85 dofs_per_cell = Utilities::fixed_power<dim>(n_dofs_1d);
89 const unsigned int array_size = n_dofs_1d*n_q_points_1d;
90 this->shape_gradients.resize_fast (array_size);
91 this->shape_values.resize_fast (array_size);
92 this->shape_hessians.resize_fast (array_size);
94 this->face_value[0].resize(n_dofs_1d);
95 this->face_gradient[0].resize(n_dofs_1d);
96 this->subface_value[0].resize(array_size);
97 this->face_value[1].resize(n_dofs_1d);
98 this->face_gradient[1].resize(n_dofs_1d);
99 this->subface_value[1].resize(array_size);
100 this->shape_values_number.resize (array_size);
101 this->shape_gradient_number.resize (array_size);
103 for (
unsigned int i=0; i<n_dofs_1d; ++i)
107 const unsigned int my_i = lexicographic[i];
108 for (
unsigned int q=0; q<n_q_points_1d; ++q)
116 shape_values_number[i*n_q_points_1d+q] = fe.
shape_value(my_i,q_point);
117 shape_gradient_number[i*n_q_points_1d+q] = fe.
shape_grad (my_i,q_point)[0];
118 shape_values [i*n_q_points_1d+q] =
119 shape_values_number [i*n_q_points_1d+q];
120 shape_gradients[i*n_q_points_1d+q] =
121 shape_gradient_number[i*n_q_points_1d+q];
122 shape_hessians[i*n_q_points_1d+q] =
125 subface_value[0][i*n_q_points_1d+q] = fe.
shape_value(my_i,q_point);
127 subface_value[1][i*n_q_points_1d+q] = fe.
shape_value(my_i,q_point);
130 this->face_value[0][i] = fe.
shape_value(my_i,q_point);
131 this->face_gradient[0][i] = fe.
shape_grad(my_i,q_point)[0];
133 this->face_value[1][i] = fe.
shape_value(my_i,q_point);
134 this->face_gradient[1][i] = fe.
shape_grad(my_i,q_point)[0];
139 this->face_indices.reinit(n_faces, this->dofs_per_face);
144 for (
unsigned int i=0; i<this->dofs_per_face; i++)
146 const unsigned int jump_term =
147 this->dofs_per_face*((i*n_dofs_1d)/this->dofs_per_face);
148 this->face_indices(0,i) = i*n_dofs_1d;
149 this->face_indices(1,i) = i*n_dofs_1d + n_dofs_1d-1;
150 this->face_indices(2,i) = i%n_dofs_1d + jump_term;
151 this->face_indices(3,i) = (i%n_dofs_1d + jump_term +
152 (n_dofs_1d-1)*n_dofs_1d);
153 this->face_indices(4,i) = i;
154 this->face_indices(5,i) = (n_dofs_1d-1)*this->dofs_per_face+i;
160 for (
unsigned int i=0; i<n_dofs_1d; i++)
162 this->face_indices(0,i) = n_dofs_1d*i;
163 this->face_indices(1,i) = n_dofs_1d*i + n_dofs_1d-1;
164 this->face_indices(2,i) = i;
165 this->face_indices(3,i) = (n_dofs_1d-1)*n_dofs_1d+i;
171 this->face_indices(0,0) = 0;
172 this->face_indices(1,0) = n_dofs_1d-1;
182 template <
typename Number>
186 std::size_t memory =
sizeof(*this);
190 memory += face_indices.memory_consumption();
191 for (
unsigned int i=0; i<2; ++i)
207 DEAL_II_NAMESPACE_CLOSE
#define AssertDimension(dim1, dim2)
::ExceptionBase & ExcMessage(std::string arg1)
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
unsigned int n_components() const
#define Assert(cond, exc)
std::size_t memory_consumption(const T &t)
const unsigned int degree
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
const std::vector< Point< dim > > & get_points() const
void reinit(const Quadrature< 1 > &quad, const FiniteElement< dim > &fe_dim)
unsigned int size() const
std::size_t memory_consumption() const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
std::vector< unsigned int > get_poly_space_numbering_inverse() const
::ExceptionBase & ExcNotImplemented()
::ExceptionBase & ExcInternalError()
const unsigned int dofs_per_cell