Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_data.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2001 - 2018 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #include <deal.II/base/geometry_info.h>
17 
18 #include <deal.II/fe/fe.h>
19 
20 DEAL_II_NAMESPACE_OPEN
21 
22 template <int dim>
24  const std::vector<unsigned int> &dofs_per_object,
25  const unsigned int n_components,
26  const unsigned int degree,
27  const Conformity conformity,
28  const BlockIndices & block_indices)
29  : dofs_per_vertex(dofs_per_object[0])
30  , dofs_per_line(dofs_per_object[1])
31  , dofs_per_quad(dim > 1 ? dofs_per_object[2] : 0)
32  , dofs_per_hex(dim > 2 ? dofs_per_object[3] : 0)
33  , first_line_index(GeometryInfo<dim>::vertices_per_cell * dofs_per_vertex)
34  , first_quad_index(first_line_index +
35  GeometryInfo<dim>::lines_per_cell * dofs_per_line)
36  , first_hex_index(first_quad_index +
37  GeometryInfo<dim>::quads_per_cell * dofs_per_quad)
38  , first_face_line_index(GeometryInfo<dim - 1>::vertices_per_cell *
39  dofs_per_vertex)
40  , first_face_quad_index(
41  (dim == 3 ? GeometryInfo<dim - 1>::vertices_per_cell * dofs_per_vertex :
42  GeometryInfo<dim>::vertices_per_cell * dofs_per_vertex) +
43  GeometryInfo<dim - 1>::lines_per_cell * dofs_per_line)
44  , dofs_per_face(GeometryInfo<dim>::vertices_per_face * dofs_per_vertex +
45  GeometryInfo<dim>::lines_per_face * dofs_per_line +
46  GeometryInfo<dim>::quads_per_face * dofs_per_quad)
47  , dofs_per_cell(GeometryInfo<dim>::vertices_per_cell * dofs_per_vertex +
48  GeometryInfo<dim>::lines_per_cell * dofs_per_line +
49  GeometryInfo<dim>::quads_per_cell * dofs_per_quad +
50  GeometryInfo<dim>::hexes_per_cell * dofs_per_hex)
51  , components(n_components)
52  , degree(degree)
53  , conforming_space(conformity)
54  , block_indices_data(block_indices.size() == 0 ?
55  BlockIndices(1, dofs_per_cell) :
56  block_indices)
57 {
58  Assert(dofs_per_object.size() == dim + 1,
59  ExcDimensionMismatch(dofs_per_object.size() - 1, dim));
60 }
61 
62 
63 
64 template <int dim>
65 bool
67 {
68  return ((dofs_per_vertex == f.dofs_per_vertex) &&
69  (dofs_per_line == f.dofs_per_line) &&
70  (dofs_per_quad == f.dofs_per_quad) &&
71  (dofs_per_hex == f.dofs_per_hex) && (components == f.components) &&
72  (degree == f.degree) && (conforming_space == f.conforming_space));
73 }
74 
75 
76 template class FiniteElementData<1>;
77 template class FiniteElementData<2>;
78 template class FiniteElementData<3>;
79 
80 DEAL_II_NAMESPACE_CLOSE
const unsigned int components
Definition: fe_base.h:292
const unsigned int dofs_per_quad
Definition: fe_base.h:237
const unsigned int degree
Definition: fe_base.h:298
const unsigned int dofs_per_line
Definition: fe_base.h:231
bool operator==(const FiniteElementData &) const
Definition: fe_data.cc:66
const unsigned int dofs_per_hex
Definition: fe_base.h:243
#define Assert(cond, exc)
Definition: exceptions.h:1407
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
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:23
const Conformity conforming_space
Definition: fe_base.h:303
const unsigned int dofs_per_vertex
Definition: fe_base.h:225
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)