17 #include <deal.II/base/quadrature.h> 18 #include <deal.II/base/qprojector.h> 19 #include <deal.II/base/template_constraints.h> 20 #include <deal.II/fe/fe_q_dg0.h> 21 #include <deal.II/fe/fe_nothing.h> 22 #include <deal.II/base/quadrature_lib.h> 23 #include <deal.II/dofs/dof_accessor.h> 29 DEAL_II_NAMESPACE_OPEN
32 template <
int dim,
int spacedim>
40 get_riaf_vector(degree))
43 ExcMessage (
"This element can only be used for polynomial degrees " 44 "greater than zero"));
50 for (
unsigned int d=0; d<dim; ++d)
58 template <
int dim,
int spacedim>
72 ExcMessage (
"This element can only be used for polynomial degrees " 79 for (
unsigned int d=0; d<dim; ++d)
87 template <
int dim,
int spacedim>
95 std::ostringstream namebuf;
97 const unsigned int n_points = this->
degree +1;
98 std::vector<double> points(n_points);
101 unsigned int index = 0;
106 if ((dim>1) ? (unit_support_points[j](1)==0 &&
107 ((dim>2) ? unit_support_points[j](2)==0:
true)) :
true)
110 points[index] = unit_support_points[j](0);
112 points[n_points-1] = unit_support_points[j](0);
114 points[index-1] = unit_support_points[j](0);
120 Assert (index == n_points || (dim==1 && index == n_points+1),
121 ExcMessage (
"Could not decode support points in one coordinate direction."));
124 for (
unsigned int j=0; j<n_points; j++)
125 if (std::fabs(points[j] - (
double)j/this->
degree) > 1e-15)
134 namebuf <<
"FE_Q_DG0<" 136 <<
">(QIterated(QTrapez()," << this->
degree <<
"))";
138 namebuf <<
"FE_Q_DG0<" 140 <<
">(" << this->degree <<
")";
148 for (
unsigned int j=0; j<n_points; j++)
149 if (points[j] != points_gl.
point(j)(0))
155 namebuf <<
"FE_Q_DG0<" 157 <<
">(" << this->
degree <<
")";
159 namebuf <<
"FE_Q_DG0<" 161 <<
">(QUnknownNodes(" << this->degree <<
"))";
163 return namebuf.str();
168 template <
int dim,
int spacedim>
177 template <
int dim,
int spacedim>
181 std::vector<double> &nodal_dofs)
const 193 const std::pair<unsigned int, unsigned int> index
195 nodal_dofs[i] = support_point_values[i](index.first);
199 nodal_dofs[nodal_dofs.size()-1] = 0.;
204 template <
int dim,
int spacedim>
207 const std::vector<double> &values)
const 217 std::copy(values.begin(), values.end(), local_dofs.begin());
220 local_dofs[local_dofs.size()-1] = 0.;
225 template <
int dim,
int spacedim>
229 const unsigned int offset)
const 241 const std::pair<unsigned int, unsigned int> index
243 local_dofs[i] = values[i](offset+index.first);
247 local_dofs[local_dofs.size()-1] = 0.;
252 template <
int dim,
int spacedim>
255 std::vector<double> &local_dofs,
256 const VectorSlice<
const std::vector<std::vector<double> > > &values)
const 268 const std::pair<unsigned int, unsigned int> index
270 local_dofs[i] = values[index.first][i];
274 local_dofs[local_dofs.size()-1] = 0.;
279 template <
int dim,
int spacedim>
291 (dynamic_cast<const FEQDG0 *>(&x_source_fe) != 0),
308 template <
int dim,
int spacedim>
312 std::vector<bool> riaf(Utilities::fixed_power<dim>(deg+1)+1,
false);
313 riaf[riaf.size()-1]=
true;
319 template <
int dim,
int spacedim>
320 std::vector<unsigned int>
323 std::vector<unsigned int> dpo(dim+1, 1U);
324 for (
unsigned int i=1; i<dpo.size(); ++i)
325 dpo[i]=dpo[i-1]*(deg-1);
333 template <
int dim,
int spacedim>
336 const unsigned int face_index)
const 347 template <
int dim,
int spacedim>
348 std::pair<Table<2,bool>, std::vector<unsigned int> >
355 constant_modes(0, i) =
true;
358 constant_modes(1, this->dofs_per_cell-1) =
true;
360 return std::pair<Table<2,bool>, std::vector<unsigned int> >
361 (constant_modes, std::vector<unsigned int> (2, 0));
367 #include "fe_q_dg0.inst" 369 DEAL_II_NAMESPACE_CLOSE
virtual void convert_generalized_support_point_values_to_nodal_values(const std::vector< Vector< double > > &support_point_values, std::vector< double > &nodal_values) const
#define AssertDimension(dim1, dim2)
virtual FiniteElement< dim, spacedim > * clone() const
const unsigned int degree
#define AssertThrow(cond, exc)
static::ExceptionBase & ExcInterpolationNotImplemented()
const Point< dim > & point(const unsigned int i) const
std::vector< Point< dim > > unit_support_points
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual void interpolate(std::vector< double > &local_dofs, const std::vector< double > &values) const 1
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
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 void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const