17 #include <deal.II/base/qprojector.h> 18 #include <deal.II/base/quadrature.h> 19 #include <deal.II/base/quadrature_lib.h> 20 #include <deal.II/base/std_cxx14/memory.h> 21 #include <deal.II/base/template_constraints.h> 23 #include <deal.II/dofs/dof_accessor.h> 25 #include <deal.II/fe/fe_nothing.h> 26 #include <deal.II/fe/fe_q_dg0.h> 31 DEAL_II_NAMESPACE_OPEN
34 template <
int dim,
int spacedim>
44 get_riaf_vector(degree))
47 ExcMessage(
"This element can only be used for polynomial degrees " 48 "greater than zero"));
54 for (
unsigned int d = 0; d < dim; ++d)
62 template <
int dim,
int spacedim>
66 Polynomials::generate_complete_Lagrange_basis(points.get_points())),
77 ExcMessage(
"This element can only be used for polynomial degrees " 84 for (
unsigned int d = 0; d < dim; ++d)
92 template <
int dim,
int spacedim>
100 std::ostringstream namebuf;
102 const unsigned int n_points = this->
degree + 1;
103 std::vector<double> points(n_points);
107 unsigned int index = 0;
112 if ((dim > 1) ? (unit_support_points[j](1) == 0 &&
113 ((dim > 2) ? unit_support_points[j](2) == 0 :
true)) :
117 points[index] = unit_support_points[j](0);
119 points[n_points - 1] = unit_support_points[j](0);
121 points[index - 1] = unit_support_points[j](0);
127 Assert(index == n_points || (dim == 1 && index == n_points + 1),
129 "Could not decode support points in one coordinate direction."));
132 for (
unsigned int j = 0; j < n_points; j++)
133 if (std::fabs(points[j] - (
double)j / this->
degree) > 1e-15)
143 <<
">(QIterated(QTrapez()," << this->
degree <<
"))";
146 << this->degree <<
")";
153 for (
unsigned int j = 0; j < n_points; j++)
154 if (points[j] != points_gl.
point(j)(0))
164 <<
">(QUnknownNodes(" << this->degree <<
"))";
166 return namebuf.str();
171 template <
int dim,
int spacedim>
172 std::unique_ptr<FiniteElement<dim, spacedim>>
175 return std_cxx14::make_unique<FE_Q_DG0<dim, spacedim>>(*this);
180 template <
int dim,
int spacedim>
184 std::vector<double> & nodal_dofs)
const 197 const std::pair<unsigned int, unsigned int> index =
199 nodal_dofs[i] = support_point_values[i](index.first);
203 nodal_dofs[nodal_dofs.size() - 1] = 0.;
208 template <
int dim,
int spacedim>
218 (x_source_fe.
get_name().find(
"FE_Q_DG0<") == 0) ||
219 (dynamic_cast<const FEQDG0 *>(&x_source_fe) !=
nullptr),
234 template <
int dim,
int spacedim>
238 std::vector<bool> riaf(Utilities::fixed_power<dim>(deg + 1) + 1,
false);
239 riaf[riaf.size() - 1] =
true;
245 template <
int dim,
int spacedim>
246 std::vector<unsigned int>
249 std::vector<unsigned int> dpo(dim + 1, 1U);
250 for (
unsigned int i = 1; i < dpo.size(); ++i)
251 dpo[i] = dpo[i - 1] * (deg - 1);
259 template <
int dim,
int spacedim>
262 const unsigned int shape_index,
263 const unsigned int face_index)
const 275 template <
int dim,
int spacedim>
276 std::pair<Table<2, bool>, std::vector<unsigned int>>
283 constant_modes(0, i) =
true;
286 constant_modes(1, this->dofs_per_cell - 1) =
true;
288 return std::pair<Table<2, bool>, std::vector<unsigned int>>(
289 constant_modes, std::vector<unsigned int>(2, 0));
295 #include "fe_q_dg0.inst" 297 DEAL_II_NAMESPACE_CLOSE
#define AssertDimension(dim1, dim2)
const unsigned int degree
#define AssertThrow(cond, exc)
const Point< dim > & point(const unsigned int i) const
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
std::vector< Point< dim > > unit_support_points
const std::vector< Point< dim > > & get_points() const
unsigned int size() const
static::ExceptionBase & ExcMessage(std::string arg1)
#define Assert(cond, exc)
static::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
virtual std::string get_name() const =0
unsigned int n_components() const
std::string dim_string(const int dim, const int spacedim)
const unsigned int dofs_per_cell
static std::vector< bool > get_riaf_vector(const unsigned int degree)
void initialize(const std::vector< Point< 1 >> &support_points_1d)
FE_Q_DG0(const unsigned int p)
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const override
virtual std::string get_name() const override
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)