Reference documentation for deal.II version 9.2.0
\(\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\}}\)
fe_bernardi_raugel.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2004 - 2020 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 
21 
23 
24 #include <deal.II/fe/fe.h>
26 #include <deal.II/fe/fe_tools.h>
27 #include <deal.II/fe/fe_values.h>
28 #include <deal.II/fe/mapping.h>
29 
30 #include <deal.II/grid/tria.h>
32 
33 #include <iostream>
34 #include <sstream>
35 
36 
38 
39 template <int dim>
41  : FE_PolyTensor<dim>(
43  FiniteElementData<dim>(get_dpo_vector(),
44  dim,
45  2,
46  FiniteElementData<dim>::Hdiv),
47  std::vector<bool>(PolynomialsBernardiRaugel<dim>::n_polynomials(p), true),
48  std::vector<ComponentMask>(PolynomialsBernardiRaugel<dim>::n_polynomials(
49  p),
50  std::vector<bool>(dim, true)))
51 {
52  Assert(dim == 2 || dim == 3, ExcImpossibleInDim(dim));
53  Assert(p == 1, ExcMessage("Only BR1 elements are available"));
54 
55  // const unsigned int n_dofs = this->dofs_per_cell;
56 
57  this->mapping_kind = {mapping_none};
58  // These must be done first, since
59  // they change the evaluation of
60  // basis functions
61 
62  // Set up the generalized support
63  // points
65 }
66 
67 
68 
69 template <int dim>
70 std::string
72 {
73  std::ostringstream namebuf;
74  namebuf << "FE_BR<" << dim << ">(" << 1 << ")";
75 
76  return namebuf.str();
77 }
78 
79 
80 
81 template <int dim>
82 std::unique_ptr<FiniteElement<dim, dim>>
84 {
85  return std_cxx14::make_unique<FE_BernardiRaugel<dim>>(*this);
86 }
87 
88 
89 
90 template <int dim>
91 void
93  const std::vector<Vector<double>> &support_point_values,
94  std::vector<double> & nodal_values) const
95 {
96  Assert(support_point_values.size() == this->generalized_support_points.size(),
97  ExcDimensionMismatch(support_point_values.size(),
98  this->generalized_support_points.size()));
99  AssertDimension(support_point_values[0].size(), dim);
100  Assert(nodal_values.size() == this->dofs_per_cell,
101  ExcDimensionMismatch(nodal_values.size(), this->dofs_per_cell));
102 
103  std::vector<Tensor<1, dim>> normals;
104  for (unsigned int i : GeometryInfo<dim>::face_indices())
105  {
106  Tensor<1, dim> normal;
107  normal[i / 2] = 1;
108  normals.push_back(normal);
109  }
110 
111  for (unsigned int i = 0; i < dim * GeometryInfo<dim>::vertices_per_cell; ++i)
112  nodal_values[i] = support_point_values[i][i % dim];
113 
114  // Compute the support points values for the bubble functions
115  for (unsigned int i = dim * GeometryInfo<dim>::vertices_per_cell;
116  i < dim * GeometryInfo<dim>::vertices_per_cell +
118  ++i)
119  {
120  nodal_values[i] = 0;
121  for (unsigned int j = 0; j < dim; ++j)
122  nodal_values[i] +=
123  support_point_values[i][j] *
124  normals[i - dim * GeometryInfo<dim>::vertices_per_cell][j];
125  }
126 }
127 
128 
129 
130 template <int dim>
131 std::vector<unsigned int>
133 {
134  // compute the number of unknowns per cell interior/face/edge
135  //
136  // there are <tt>dim</tt> degrees of freedom per vertex and there
137  // is 1 degree of freedom per edge in 2D (face in 3D)
138  std::vector<unsigned int> dpo(dim + 1, 0u);
139  dpo[0] = dim;
140  dpo[dim - 1] = 1u;
141 
142  return dpo;
143 }
144 
145 
146 
147 template <int dim>
148 void
150 {
151  // The support points for our shape functions are the vertices and
152  // the face midpoints, for a total of #vertices + #faces points
153  this->generalized_support_points.resize(this->dofs_per_cell);
154 
155  // We need dim copies of each vertex for the first dim*vertices_per_cell
156  // generalized support points
157  for (unsigned int i = 0; i < dim * GeometryInfo<dim>::vertices_per_cell; ++i)
158  this->generalized_support_points[i] =
160 
161  // The remaining 2*dim points are the edge midpoints
162  for (unsigned int i = 0; i < dim; ++i)
163  {
164  for (unsigned int j = 0; j < 2; ++j)
165  {
166  Point<dim> p;
167  p[0] = 0.5;
168  p[1] = 0.5;
169  if (dim == 3)
170  p[2] = 0.5;
171  p[i] = j;
172 
173  const unsigned int k =
174  dim * GeometryInfo<dim>::vertices_per_cell + i * 2 + j;
175  this->generalized_support_points[k] = p;
176  }
177  }
178 }
179 
180 template class FE_BernardiRaugel<2>;
181 template class FE_BernardiRaugel<3>;
182 
fe_values.h
FE_BernardiRaugel::clone
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
Definition: fe_bernardi_raugel.cc:83
tria.h
FiniteElementData
Definition: fe_base.h:147
tria_iterator.h
GeometryInfo
Definition: geometry_info.h:1224
fe_bernardi_raugel.h
FE_PolyTensor< dim >::mapping_kind
std::vector< MappingKind > mapping_kind
Definition: fe_poly_tensor.h:220
ComponentMask
Definition: component_mask.h:83
FE_BernardiRaugel::get_name
virtual std::string get_name() const override
Definition: fe_bernardi_raugel.cc:71
quadrature_lib.h
bool
FE_BernardiRaugel::convert_generalized_support_point_values_to_dof_values
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
Definition: fe_bernardi_raugel.cc:92
FE_BernardiRaugel::get_dpo_vector
static std::vector< unsigned int > get_dpo_vector()
Definition: fe_bernardi_raugel.cc:132
StandardExceptions::ExcMessage
static ::ExceptionBase & ExcMessage(std::string arg1)
FE_BernardiRaugel::FE_BernardiRaugel
FE_BernardiRaugel(const unsigned int p=1)
Definition: fe_bernardi_raugel.cc:40
Tensor< 1, dim >
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
fe.h
AssertDimension
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1579
FE_BernardiRaugel::initialize_support_points
void initialize_support_points()
Definition: fe_bernardi_raugel.cc:149
fe_tools.h
mapping.h
Assert
#define Assert(cond, exc)
Definition: exceptions.h:1419
StandardExceptions::ExcImpossibleInDim
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
dof_accessor.h
StandardExceptions::ExcDimensionMismatch
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Point< dim >
FE_PolyTensor
Definition: fe_poly_tensor.h:143
FE_BernardiRaugel
Definition: fe_bernardi_raugel.h:87
memory.h
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
PolynomialsBernardiRaugel
Definition: polynomials_bernardi_raugel.h:87
Vector< double >
polynomials_p.h
mapping_none
@ mapping_none
Definition: mapping.h:68
qprojector.h