Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-2846-g6fb608615f 2025-03-15 04:10:00+00:00
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
fe_q_iso_q1.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2013 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
18
19#include <deal.II/fe/fe_dgq.h>
22#include <deal.II/fe/fe_tools.h>
23
24#include <deal.II/lac/vector.h>
25
26#include <memory>
27#include <sstream>
28#include <vector>
29
31
32
33namespace internal
34{
35 template <int dim>
37 make_local_sparsity_pattern(unsigned int n_points_1d)
38 {
39 const unsigned int n_per_cell = Utilities::pow(n_points_1d, dim);
40
41 const std::vector<unsigned int> R =
42 FETools::lexicographic_to_hierarchic_numbering<dim>(n_points_1d - 1);
43
44 Table<2, bool> table;
45 table.reinit(n_per_cell, n_per_cell);
46 table.fill(false);
47
48 const int N1d = n_points_1d;
49 for (unsigned int i = 0; i < n_per_cell; ++i)
50 for (unsigned int j = 0; j < n_per_cell; ++j)
51 {
52 // compute l1 distance:
53 int distance = 0;
54
55 int xi = i;
56 int xj = j;
57 for (unsigned int d = 0; d < dim; ++d)
58 {
59 int current_distance = std::abs((xi % N1d) - (xj % N1d));
60 xi /= N1d;
61 xj /= N1d;
62 distance = std::max(distance, current_distance);
63 if (distance > 1)
64 break;
65 }
66
67 if (distance <= 1)
68 table(R[i], R[j]) = true;
69 }
70
71 return table;
72 }
73} // namespace internal
74
75
76
77template <int dim, int spacedim>
78FE_Q_iso_Q1<dim, spacedim>::FE_Q_iso_Q1(const unsigned int subdivisions)
79 : FE_Q_Base<dim, spacedim>(
80 TensorProductPolynomials<dim, Polynomials::PiecewisePolynomial<double>>(
81 Polynomials::generate_complete_Lagrange_basis_on_subdivisions(
82 subdivisions,
83 1)),
84 FiniteElementData<dim>(this->get_dpo_vector(subdivisions),
85 1,
86 subdivisions,
87 FiniteElementData<dim>::H1),
88 std::vector<bool>(1, false))
89{
90 Assert(subdivisions > 0,
91 ExcMessage("This element can only be used with a positive number of "
92 "subelements"));
93
94 const QTrapezoid<1> trapez;
95 const QIterated<1> points(trapez, subdivisions);
96
97 this->initialize(points.get_points());
99 internal::make_local_sparsity_pattern<dim>(subdivisions + 1);
100}
101
102
103
104template <int dim, int spacedim>
106 const std::vector<Point<1>> &support_points)
107 : FE_Q_Base<dim, spacedim>(
108 TensorProductPolynomials<dim, Polynomials::PiecewisePolynomial<double>>(
109 Polynomials::generate_complete_linear_basis_on_subdivisions(
110 support_points)),
111 FiniteElementData<dim>(this->get_dpo_vector(support_points.size() - 1),
112 1,
113 support_points.size() - 1,
114 FiniteElementData<dim>::H1),
115 std::vector<bool>(1, false))
116{
117 Assert(support_points.size() > 1,
118 ExcMessage("This element can only be used with a positive number of "
119 "subelements"));
120
121 this->initialize(support_points);
123 internal::make_local_sparsity_pattern<dim>(support_points.size());
124}
125
126
127
128template <int dim, int spacedim>
129std::string
131{
132 // note that the FETools::get_fe_by_name function depends on the
133 // particular format of the string this function returns, so they have to be
134 // kept in sync
135
136 std::ostringstream namebuf;
137 namebuf << "FE_Q_iso_Q1<" << Utilities::dim_string(dim, spacedim) << ">("
138 << this->degree << ")";
139 return namebuf.str();
140}
141
142
143
144template <int dim, int spacedim>
145void
148 const std::vector<Vector<double>> &support_point_values,
149 std::vector<double> &nodal_values) const
150{
151 AssertDimension(support_point_values.size(),
152 this->get_unit_support_points().size());
153 AssertDimension(support_point_values.size(), nodal_values.size());
154 AssertDimension(this->n_dofs_per_cell(), nodal_values.size());
155
156 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
157 {
158 AssertDimension(support_point_values[i].size(), 1);
159
160 nodal_values[i] = support_point_values[i](0);
161 }
162}
163
164
165
166template <int dim, int spacedim>
167std::unique_ptr<FiniteElement<dim, spacedim>>
169{
170 return std::make_unique<FE_Q_iso_Q1<dim, spacedim>>(*this);
171}
172
173
174
175template <int dim, int spacedim>
178 const FiniteElement<dim, spacedim> &fe_other,
179 const unsigned int codim) const
180{
181 Assert(codim <= dim, ExcImpossibleInDim(dim));
182 (void)codim;
183
184 // vertex/line/face domination
185 // (if fe_other is derived from FE_DGQ)
186 // ------------------------------------
187 if (codim > 0)
188 if (dynamic_cast<const FE_DGQ<dim, spacedim> *>(&fe_other) != nullptr)
189 // there are no requirements between continuous and discontinuous elements
191
192 // vertex/line/face domination
193 // (if fe_other is not derived from FE_DGQ)
194 // & cell domination
195 // ----------------------------------------
196 if (const FE_Q_iso_Q1<dim, spacedim> *fe_q_iso_q1_other =
197 dynamic_cast<const FE_Q_iso_Q1<dim, spacedim> *>(&fe_other))
198 {
199 // different behavior as in FE_Q: as FE_Q_iso_Q1(2) is not a subspace of
200 // FE_Q_iso_Q1(3), need that the element degrees are multiples of each
201 // other
202 if (this->degree < fe_q_iso_q1_other->degree &&
203 fe_q_iso_q1_other->degree % this->degree == 0)
205 else if (this->degree == fe_q_iso_q1_other->degree)
207 else if (this->degree > fe_q_iso_q1_other->degree &&
208 this->degree % fe_q_iso_q1_other->degree == 0)
210 else
212 }
213 else if (const FE_Nothing<dim> *fe_nothing =
214 dynamic_cast<const FE_Nothing<dim> *>(&fe_other))
215 {
216 if (fe_nothing->is_dominating())
218 else
219 // the FE_Nothing has no degrees of freedom and it is typically used
220 // in a context where we don't require any continuity along the
221 // interface
223 }
224
227}
228
229
230
231// explicit instantiations
232#include "fe/fe_q_iso_q1.inst"
233
void initialize(const std::vector< Point< 1 > > &support_points_1d)
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const override final
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
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const override
FE_Q_iso_Q1(const unsigned int n_subdivisions)
Table< 2, bool > local_dof_sparsity_pattern
Definition fe.h:2722
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
#define DEAL_II_NOT_IMPLEMENTED()
#define Assert(cond, exc)
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcMessage(std::string arg1)
std::size_t size
Definition mpi.cc:745
std::string dim_string(const int dim, const int spacedim)
Definition utilities.cc:551
constexpr T pow(const T base, const int iexp)
Definition utilities.h:967
Table< 2, bool > make_local_sparsity_pattern(unsigned int n_points_1d)
STL namespace.
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > abs(const ::VectorizedArray< Number, width > &)