Reference documentation for deal.II version GIT relicensing-437-g81ec864850 2024-04-19 07:30:02+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\}}\)
Loading...
Searching...
No Matches
fe_data.cc
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2001 - 2023 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
16
17#include <deal.II/fe/fe.h>
18
20
21namespace internal
22{
24 expand(const unsigned int dim,
25 const std::vector<unsigned int> &dofs_per_object,
26 const ReferenceCell cell_type)
27 {
29
30 const unsigned int face_no = 0;
31
32 result.dofs_per_object_exclusive.resize(4, std::vector<unsigned int>(1));
33 result.dofs_per_object_inclusive.resize(4, std::vector<unsigned int>(1));
34 result.object_index.resize(4, std::vector<unsigned int>(1));
35 result.first_object_index_on_face.resize(3, std::vector<unsigned int>(1));
36
37 // dofs_per_vertex
38 const unsigned int dofs_per_vertex = dofs_per_object[0];
39 result.dofs_per_object_exclusive[0][0] = dofs_per_vertex;
40
41 // dofs_per_line
42 const unsigned int dofs_per_line = dofs_per_object[1];
43 result.dofs_per_object_exclusive[1][0] = dofs_per_line;
44
45 // dofs_per_quad
46 const unsigned int dofs_per_quad = dim > 1 ? dofs_per_object[2] : 0;
47 result.dofs_per_object_exclusive[2][0] = dofs_per_quad;
48
49 // dofs_per_hex
50 const unsigned int dofs_per_hex = dim > 2 ? dofs_per_object[3] : 0;
51 result.dofs_per_object_exclusive[3][0] = dofs_per_hex;
52
53
54 // first_line_index
55 const unsigned int first_line_index =
56 (cell_type.n_vertices() * dofs_per_vertex);
57 result.object_index[1][0] = first_line_index;
58
59 // first_quad_index
60 const unsigned int first_quad_index =
61 (first_line_index + cell_type.n_lines() * dofs_per_line);
62 result.object_index[2][0] = first_quad_index;
63
64 // first_hex_index
65 result.object_index[3][0] =
66 (first_quad_index +
67 (dim == 2 ? 1 : (dim == 3 ? cell_type.n_faces() : 0)) * dofs_per_quad);
68
69 // first_face_line_index
70 result.first_object_index_on_face[1][0] =
71 (cell_type.face_reference_cell(face_no).n_vertices() * dofs_per_vertex);
72
73 // first_face_quad_index
74 result.first_object_index_on_face[2][0] =
75 ((dim == 3 ? cell_type.face_reference_cell(face_no).n_vertices() *
76 dofs_per_vertex :
77 cell_type.n_vertices() * dofs_per_vertex) +
78 cell_type.face_reference_cell(face_no).n_lines() * dofs_per_line);
79
80 // dofs_per_face
81 result.dofs_per_object_inclusive[dim - 1][0] =
82 (cell_type.face_reference_cell(face_no).n_vertices() * dofs_per_vertex +
83 cell_type.face_reference_cell(face_no).n_lines() * dofs_per_line +
84 (dim == 3 ? 1 : 0) * dofs_per_quad);
85
86
87 // dofs_per_cell
88 result.dofs_per_object_inclusive[dim][0] =
89 (cell_type.n_vertices() * dofs_per_vertex +
90 cell_type.n_lines() * dofs_per_line +
91 (dim == 2 ? 1 : (dim == 3 ? cell_type.n_faces() : 0)) * dofs_per_quad +
92 (dim == 3 ? 1 : 0) * dofs_per_hex);
93
94 return result;
95 }
96
97 unsigned int
98 number_unique_entries(const std::vector<unsigned int> &vector)
99 {
100 if (std::all_of(vector.begin(), vector.end(), [&](const auto &e) {
101 return e == vector.front();
102 }))
103 {
104 return 1;
105 }
106 else
107 return vector.size();
108 }
109} // namespace internal
110
111
112
113template <int dim>
115 const std::vector<unsigned int> &dofs_per_object,
116 const unsigned int n_components,
117 const unsigned int degree,
118 const Conformity conformity,
119 const BlockIndices &block_indices)
120 : FiniteElementData(dofs_per_object,
121 dim == 0 ?
122 ReferenceCells::Vertex :
123 (dim == 1 ? ReferenceCells::Line :
124 (dim == 2 ? ReferenceCells::Quadrilateral :
125 ReferenceCells::Hexahedron)),
126 n_components,
127 degree,
128 conformity,
129 block_indices)
130{}
131
132template <int dim>
134 const std::vector<unsigned int> &dofs_per_object,
135 const ReferenceCell cell_type,
136 const unsigned int n_components,
137 const unsigned int degree,
138 const Conformity conformity,
139 const BlockIndices &block_indices)
140 : FiniteElementData(internal::expand(dim, dofs_per_object, cell_type),
141 cell_type,
142 n_components,
143 degree,
144 conformity,
145 block_indices)
146{}
147
148
149
150template <int dim>
153 const ReferenceCell reference_cell,
154 const unsigned int n_components,
155 const unsigned int degree,
156 const Conformity conformity,
157 const BlockIndices &block_indices)
158 : reference_cell_kind(reference_cell)
159 , number_of_unique_2d_subobjects(
160 internal::number_unique_entries(data.dofs_per_object_inclusive[2]))
161 , number_unique_faces(
162 internal::number_unique_entries(data.dofs_per_object_inclusive[dim - 1]))
163 , dofs_per_vertex(data.dofs_per_object_exclusive[0][0])
164 , dofs_per_line(data.dofs_per_object_exclusive[1][0])
165 , n_dofs_on_quad(data.dofs_per_object_exclusive[2])
166 , dofs_per_quad(n_dofs_on_quad[0])
167 , dofs_per_quad_max(
168 *max_element(n_dofs_on_quad.begin(), n_dofs_on_quad.end()))
169 , dofs_per_hex(data.dofs_per_object_exclusive[3][0])
170 , first_line_index(data.object_index[1][0])
171 , first_index_of_quads(data.object_index[2])
172 , first_quad_index(first_index_of_quads[0])
173 , first_hex_index(data.object_index[3][0])
174 , first_line_index_of_faces(data.first_object_index_on_face[1])
175 , first_face_line_index(first_line_index_of_faces[0])
176 , first_quad_index_of_faces(data.first_object_index_on_face[2])
177 , first_face_quad_index(first_quad_index_of_faces[0])
178 , n_dofs_on_face(data.dofs_per_object_inclusive[dim - 1])
179 , dofs_per_face(n_dofs_on_face[0])
180 , dofs_per_face_max(
181 *max_element(n_dofs_on_face.begin(), n_dofs_on_face.end()))
182 , dofs_per_cell(data.dofs_per_object_inclusive[dim][0])
183 , components(n_components)
184 , degree(degree)
185 , conforming_space(conformity)
186 , block_indices_data(block_indices.size() == 0 ?
187 BlockIndices(1, dofs_per_cell) :
188 block_indices)
189{}
190
191
192
193template <int dim>
194bool
196{
197 return ((dofs_per_vertex == f.dofs_per_vertex) &&
198 (dofs_per_line == f.dofs_per_line) &&
199 (dofs_per_quad == f.dofs_per_quad) &&
200 (dofs_per_hex == f.dofs_per_hex) && (components == f.components) &&
201 (degree == f.degree) && (conforming_space == f.conforming_space));
202}
203
204
205template class FiniteElementData<1>;
206template class FiniteElementData<2>;
207template class FiniteElementData<3>;
208
const unsigned int dofs_per_quad
Definition fe_data.h:346
const unsigned int components
Definition fe_data.h:446
const unsigned int degree
Definition fe_data.h:452
bool operator==(const FiniteElementData &) const
Definition fe_data.cc:195
const unsigned int dofs_per_line
Definition fe_data.h:332
const unsigned int dofs_per_hex
Definition fe_data.h:359
const Conformity conforming_space
Definition fe_data.h:457
const unsigned int dofs_per_vertex
Definition fe_data.h:326
FiniteElementData(const std::vector< unsigned int > &dofs_per_object, const unsigned int n_components, const unsigned int degree, const Conformity conformity=unknown, const BlockIndices &block_indices=BlockIndices())
Definition fe_data.cc:114
unsigned int n_vertices() const
unsigned int n_faces() const
unsigned int n_lines() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
unsigned int number_unique_entries(const std::vector< unsigned int > &vector)
Definition fe_data.cc:98
internal::GenericDoFsPerObject expand(const unsigned int dim, const std::vector< unsigned int > &dofs_per_object, const ReferenceCell reference_cell)
Definition fe_data.cc:24
std::vector< std::vector< unsigned int > > object_index
Definition fe_data.h:197
std::vector< std::vector< unsigned int > > first_object_index_on_face
Definition fe_data.h:202
std::vector< std::vector< unsigned int > > dofs_per_object_inclusive
Definition fe_data.h:192
std::vector< std::vector< unsigned int > > dofs_per_object_exclusive
Definition fe_data.h:187