17 #include <deal.II/base/quadrature_lib.h> 18 #include <deal.II/base/std_cxx14/memory.h> 20 #include <deal.II/fe/fe_q.h> 22 #include <deal.II/lac/vector.h> 27 DEAL_II_NAMESPACE_OPEN
37 get_QGaussLobatto_points(
const unsigned int degree)
46 return std::vector<Point<1>>();
54 template <
int dim,
int spacedim>
64 std::vector<bool>(1, false))
66 this->
initialize(internal::FE_Q::get_QGaussLobatto_points(degree));
71 template <
int dim,
int spacedim>
75 Polynomials::generate_complete_Lagrange_basis(points.get_points())),
80 std::vector<bool>(1, false))
87 template <
int dim,
int spacedim>
95 std::ostringstream namebuf;
96 bool equidistant =
true;
97 std::vector<double> points(this->
degree + 1);
100 std::vector<unsigned int> lexicographic =
102 for (
unsigned int j = 0; j <= this->
degree; j++)
106 for (
unsigned int j = 0; j <= this->
degree; j++)
107 if (std::fabs(points[j] - (
double)j / this->
degree) > 1e-15)
113 if (equidistant ==
true)
117 <<
">(QIterated(QTrapez()," << this->
degree <<
"))";
120 << this->degree <<
")";
126 bool gauss_lobatto =
true;
127 for (
unsigned int j = 0; j <= this->
degree; j++)
128 if (points[j] != points_gl.
point(j)(0))
130 gauss_lobatto =
false;
133 if (gauss_lobatto ==
true)
138 <<
">(QUnknownNodes(" << this->degree <<
"))";
140 return namebuf.str();
145 template <
int dim,
int spacedim>
149 std::vector<double> & nodal_values)
const 160 nodal_values[i] = support_point_values[i](0);
166 template <
int dim,
int spacedim>
167 std::unique_ptr<FiniteElement<dim, spacedim>>
170 return std_cxx14::make_unique<FE_Q<dim, spacedim>>(*this);
177 DEAL_II_NAMESPACE_CLOSE
const std::vector< Point< dim > > & get_unit_support_points() const
virtual std::string get_name() const override
#define AssertDimension(dim1, dim2)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
FE_Q(const unsigned int p)
const unsigned int degree
#define AssertThrow(cond, exc)
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
const std::vector< unsigned int > & get_numbering_inverse() const
const Point< dim > & point(const unsigned int i) const
std::vector< Point< dim > > unit_support_points
const std::vector< Point< dim > > & get_points() const
std::string dim_string(const int dim, const int spacedim)
const unsigned int dofs_per_cell
void initialize(const std::vector< Point< 1 >> &support_points_1d)
TensorProductPolynomials< dim > poly_space