Reference documentation for deal.II version Git 7439f2b0b4 2021-09-15 12:22:58 +0200
\(\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_raviart_thomas.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2021 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 
22 #include <deal.II/base/table.h>
23 #include <deal.II/base/utilities.h>
24 
26 
27 #include <deal.II/fe/fe.h>
29 #include <deal.II/fe/fe_tools.h>
30 #include <deal.II/fe/fe_values.h>
31 #include <deal.II/fe/mapping.h>
32 
33 #include <deal.II/grid/tria.h>
35 
36 #include <iostream>
37 #include <memory>
38 #include <sstream>
39 
40 
42 
43 
44 template <int dim>
46  : FE_PolyTensor<dim>(
47  PolynomialsRaviartThomas<dim>(deg),
49  dim,
50  deg + 1,
51  FiniteElementData<dim>::Hdiv),
52  std::vector<bool>(PolynomialsRaviartThomas<dim>::n_polynomials(deg),
53  true),
54  std::vector<ComponentMask>(PolynomialsRaviartThomas<dim>::n_polynomials(
55  deg),
56  std::vector<bool>(dim, true)))
57 {
58  Assert(dim >= 2, ExcImpossibleInDim(dim));
59  const unsigned int n_dofs = this->n_dofs_per_cell();
60 
62  // First, initialize the
63  // generalized support points and
64  // quadrature weights, since they
65  // are required for interpolation.
67 
68  // Now compute the inverse node matrix, generating the correct
69  // basis functions from the raw ones. For a discussion of what
70  // exactly happens here, see FETools::compute_node_matrix.
72  this->inverse_node_matrix.reinit(n_dofs, n_dofs);
73  this->inverse_node_matrix.invert(M);
74  // From now on, the shape functions provided by FiniteElement::shape_value
75  // and similar functions will be the correct ones, not
76  // the raw shape functions from the polynomial space anymore.
77 
78  // Reinit the vectors of
79  // restriction and prolongation
80  // matrices to the right sizes.
81  // Restriction only for isotropic
82  // refinement
84  // Fill prolongation matrices with embedding operators
87 
88  // TODO: the implementation makes the assumption that all faces have the
89  // same number of dofs
90  AssertDimension(this->n_unique_faces(), 1);
91  const unsigned int face_no = 0;
92 
93  // TODO[TL]: for anisotropic refinement we will probably need a table of
94  // submatrices with an array for each refine case
96  for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_face; ++i)
97  face_embeddings[i].reinit(this->n_dofs_per_face(face_no),
98  this->n_dofs_per_face(face_no));
99  FETools::compute_face_embedding_matrices<dim, double>(*this,
100  face_embeddings,
101  0,
102  0);
103  this->interface_constraints.reinit((1 << (dim - 1)) *
104  this->n_dofs_per_face(face_no),
105  this->n_dofs_per_face(face_no));
106  unsigned int target_row = 0;
107  for (unsigned int d = 0; d < GeometryInfo<dim>::max_children_per_face; ++d)
108  for (unsigned int i = 0; i < face_embeddings[d].m(); ++i)
109  {
110  for (unsigned int j = 0; j < face_embeddings[d].n(); ++j)
111  this->interface_constraints(target_row, j) = face_embeddings[d](i, j);
112  ++target_row;
113  }
114 
115  // We need to initialize the dof permutation table and the one for the sign
116  // change.
118 }
119 
120 
121 
122 template <int dim>
123 std::string
125 {
126  // note that the
127  // FETools::get_fe_by_name
128  // function depends on the
129  // particular format of the string
130  // this function returns, so they
131  // have to be kept in synch
132 
133  // note that this->degree is the maximal
134  // polynomial degree and is thus one higher
135  // than the argument given to the
136  // constructor
137  std::ostringstream namebuf;
138  namebuf << "FE_RaviartThomas<" << dim << ">(" << this->degree - 1 << ")";
139 
140  return namebuf.str();
141 }
142 
143 
144 template <int dim>
145 std::unique_ptr<FiniteElement<dim, dim>>
147 {
148  return std::make_unique<FE_RaviartThomas<dim>>(*this);
149 }
150 
151 
152 //---------------------------------------------------------------------------
153 // Auxiliary and internal functions
154 //---------------------------------------------------------------------------
155 
156 
157 template <int dim>
158 void
160 {
161  QGauss<dim> cell_quadrature(deg + 1);
162  const unsigned int n_interior_points = (deg > 0) ? cell_quadrature.size() : 0;
163 
164  // TODO: the implementation makes the assumption that all faces have the
165  // same number of dofs
166  AssertDimension(this->n_unique_faces(), 1);
167  const unsigned int face_no = 0;
168 
169  unsigned int n_face_points = (dim > 1) ? 1 : 0;
170  // compute (deg+1)^(dim-1)
171  for (unsigned int d = 1; d < dim; ++d)
172  n_face_points *= deg + 1;
173 
174 
175  this->generalized_support_points.resize(
176  GeometryInfo<dim>::faces_per_cell * n_face_points + n_interior_points);
177  this->generalized_face_support_points[face_no].resize(n_face_points);
178 
179  // Number of the point being entered
180  unsigned int current = 0;
181 
182  if (dim > 1)
183  {
184  QGauss<dim - 1> face_points(deg + 1);
187 
188  boundary_weights.reinit(n_face_points, legendre.n());
189 
190  for (unsigned int k = 0; k < n_face_points; ++k)
191  {
192  this->generalized_face_support_points[face_no][k] =
193  face_points.point(k);
194  // Compute its quadrature
195  // contribution for each
196  // moment.
197  for (unsigned int i = 0; i < legendre.n(); ++i)
198  {
199  boundary_weights(k, i) =
200  face_points.weight(k) *
201  legendre.compute_value(i, face_points.point(k));
202  }
203  }
204 
205  Quadrature<dim> faces =
207  face_points);
208  for (; current < GeometryInfo<dim>::faces_per_cell * n_face_points;
209  ++current)
210  {
211  // Enter the support point
212  // into the vector
213  this->generalized_support_points[current] = faces.point(
214  current +
216  this->reference_cell(), 0, true, false, false, n_face_points));
217  }
218  }
219 
220  if (deg == 0)
221  return;
222 
223  // Create Legendre basis for the space D_xi Q_k
224  std::unique_ptr<AnisotropicPolynomials<dim>> polynomials[dim];
225  for (unsigned int dd = 0; dd < dim; ++dd)
226  {
227  std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
228  for (unsigned int d = 0; d < dim; ++d)
231 
232  polynomials[dd] = std::make_unique<AnisotropicPolynomials<dim>>(poly);
233  }
234 
235  interior_weights.reinit(
236  TableIndices<3>(n_interior_points, polynomials[0]->n(), dim));
237 
238  for (unsigned int k = 0; k < cell_quadrature.size(); ++k)
239  {
240  this->generalized_support_points[current++] = cell_quadrature.point(k);
241  for (unsigned int i = 0; i < polynomials[0]->n(); ++i)
242  for (unsigned int d = 0; d < dim; ++d)
243  interior_weights(k, i, d) =
244  cell_quadrature.weight(k) *
245  polynomials[d]->compute_value(i, cell_quadrature.point(k));
246  }
247 
248  Assert(current == this->generalized_support_points.size(),
249  ExcInternalError());
250 }
251 
252 
253 template <int dim>
254 void
256 {
257  // For 1D do nothing.
258  //
259  // TODO: For 2D we simply keep the legacy behavior for now. This should be
260  // changed in the future and can be taken care of by similar means as the 3D
261  // case below. The legacy behavior can be found in fe_poly_tensor.cc in the
262  // function internal::FE_PolyTensor::get_dof_sign_change_h_div(...)
263  if (dim < 3)
264  return;
265 
266  // TODO: the implementation makes the assumption that all faces have the
267  // same number of dofs
268  AssertDimension(this->n_unique_faces(), 1);
269  const unsigned int face_no = 0;
270 
272  .n_elements() == 8 * this->n_dofs_per_quad(face_no),
273  ExcInternalError());
274 
275  // The 3D RaviartThomas space has tensor_degree*tensor_degree face dofs
276  const unsigned int n = this->tensor_degree();
277  Assert(n * n == this->n_dofs_per_quad(face_no), ExcInternalError());
278 
279  // The vector of tables adjust_quad_dof_index_for_face_orientation_table
280  // contains offsets for local face_dofs and is being filled with zeros in
281  // fe.cc. We need to fill it with the correct values in case of non-standard,
282  // flipped (rotated by +180 degrees) or rotated (rotated by +90 degrees) faces
283 
284  // The dofs on a face are connected to a n x n matrix. for example, for
285  // tensor_degree==3 we have the following dofs on a quad:
286 
287  // ___________
288  // | |
289  // | 6 7 8 |
290  // | |
291  // | 3 4 5 |
292  // | |
293  // | 0 1 2 |
294  // |___________|
295  //
296  // We have dof_index=i+n*j with index i in x-direction and index j in
297  // y-direction running from 0 to n-1. to extract i and j we can use
298  // i=dof_index%n and j=dof_index/n. The indices i and j can then be used to
299  // compute the offset.
300 
301  // Example: if the switches are (true | true | true) that means we rotate the
302  // face first by + 90 degree(counterclockwise) then by another +180
303  // degrees but we do not flip it since the face has standard
304  // orientation. The flip axis is the diagonal from the lower left to the upper
305  // right corner of the face. With these flags the configuration above becomes:
306  // ___________
307  // | |
308  // | 2 5 8 |
309  // | |
310  // | 1 4 7 |
311  // | |
312  // | 0 3 6 |
313  // |___________|
314  //
315  // This is exactly what is realized by the formulas implemented below. Note
316  // that the necessity of a permuattion depends on the three flags.
317 
318  // There is also a pattern for the sign change of the mapped face_dofs
319  // depending on the switches. In the above example it would be
320  // ___________
321  // | |
322  // | + - + |
323  // | |
324  // | + - + |
325  // | |
326  // | + - + |
327  // |___________|
328  //
329 
330  for (unsigned int local_face_dof = 0;
331  local_face_dof < this->n_dofs_per_quad(face_no);
332  ++local_face_dof)
333  {
334  // Row and column
335  unsigned int i = local_face_dof % n;
336  unsigned int j = local_face_dof / n;
337 
338  // We have 8 cases that are all treated the same way. Note that the
339  // corresponding case to case_no is just its binary representation.
340  // The above example of (false | true | true) would be case_no=3
341  for (unsigned int case_no = 0; case_no < 8; ++case_no)
342  {
343  // Get the binary representation of the case
344  const bool face_orientation = Utilities::get_bit(case_no, 2);
345  const bool face_flip = Utilities::get_bit(case_no, 1);
346  const bool face_rotation = Utilities::get_bit(case_no, 0);
347 
348  if (((!face_orientation) && (!face_rotation)) ||
349  ((face_orientation) && (face_rotation)))
350  {
351  // We flip across the diagonal
352  // This is the local face dof offset
354  local_face_dof, case_no) = j + i * n - local_face_dof;
355  }
356  else
357  {
358  // Offset is zero
360  local_face_dof, case_no) = 0;
361  } // if face needs dof permutation
362 
363  // Get new local face_dof by adding offset
364  const unsigned int new_local_face_dof =
365  local_face_dof +
367  local_face_dof, case_no);
368  // compute new row and column index
369  i = new_local_face_dof % n;
370  j = new_local_face_dof / n;
371 
372  /*
373  * Now compute if a sign change is necessary. This is done for the
374  * case of face_orientation==true
375  */
376  // flip = false, rotation=true
377  if (!face_flip && face_rotation)
378  {
380  local_face_dof, case_no) = ((j % 2) == 1);
381  }
382  // flip = true, rotation=false
383  else if (face_flip && !face_rotation)
384  {
385  // This case is symmetric (although row and column may be
386  // switched)
388  local_face_dof, case_no) = ((j % 2) == 1) != ((i % 2) == 1);
389  }
390  // flip = true, rotation=true
391  else if (face_flip && face_rotation)
392  {
394  local_face_dof, case_no) = ((i % 2) == 1);
395  }
396  /*
397  * flip = false, rotation=false => nothing to do
398  */
399 
400  /*
401  * If face_orientation==false the sign flip is exactly the opposite.
402  */
403  if (!face_orientation)
405  local_face_dof, case_no) =
407  local_face_dof, case_no);
408  } // case_no
409  } // local_face_dof
410 }
411 
412 
413 template <>
414 void
416 {
417  // there is only one refinement case in 1d,
418  // which is the isotropic one (first index of
419  // the matrix array has to be 0)
420  for (unsigned int i = 0; i < GeometryInfo<1>::max_children_per_cell; ++i)
421  this->restriction[0][i].reinit(0, 0);
422 }
423 
424 
425 
426 // This function is the same Raviart-Thomas interpolation performed by
427 // interpolate. Still, we cannot use interpolate, since it was written
428 // for smooth functions. The functions interpolated here are not
429 // smooth, maybe even not continuous. Therefore, we must double the
430 // number of quadrature points in each direction in order to integrate
431 // only smooth functions.
432 
433 // Then again, the interpolated function is chosen such that the
434 // moments coincide with the function to be interpolated.
435 
436 template <int dim>
437 void
439 {
440  const unsigned int iso = RefinementCase<dim>::isotropic_refinement - 1;
441 
442  QGauss<dim - 1> q_base(this->degree);
443  const unsigned int n_face_points = q_base.size();
444  // First, compute interpolation on
445  // subfaces
446  for (unsigned int face : GeometryInfo<dim>::face_indices())
447  {
448  // The shape functions of the
449  // child cell are evaluated
450  // in the quadrature points
451  // of a full face.
452  Quadrature<dim> q_face =
453  QProjector<dim>::project_to_face(this->reference_cell(), q_base, face);
454  // Store shape values, since the
455  // evaluation suffers if not
456  // ordered by point
457  Table<2, double> cached_values_on_face(this->n_dofs_per_cell(),
458  q_face.size());
459  for (unsigned int k = 0; k < q_face.size(); ++k)
460  for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
461  cached_values_on_face(i, k) = this->shape_value_component(
462  i, q_face.point(k), GeometryInfo<dim>::unit_normal_direction[face]);
463 
464  for (unsigned int sub = 0; sub < GeometryInfo<dim>::max_children_per_face;
465  ++sub)
466  {
467  // The weight functions for
468  // the coarse face are
469  // evaluated on the subface
470  // only.
472  this->reference_cell(), q_base, face, sub);
473  const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
475 
476  // On a certain face, we must
477  // compute the moments of ALL
478  // fine level functions with
479  // the coarse level weight
480  // functions belonging to
481  // that face. Due to the
482  // orthogonalization process
483  // when building the shape
484  // functions, these weights
485  // are equal to the
486  // corresponding shape
487  // functions.
488  for (unsigned int k = 0; k < n_face_points; ++k)
489  for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
490  ++i_child)
491  for (unsigned int i_face = 0;
492  i_face < this->n_dofs_per_face(face);
493  ++i_face)
494  {
495  // The quadrature
496  // weights on the
497  // subcell are NOT
498  // transformed, so we
499  // have to do it here.
500  this->restriction[iso][child](
501  face * this->n_dofs_per_face(face) + i_face, i_child) +=
502  Utilities::fixed_power<dim - 1>(.5) * q_sub.weight(k) *
503  cached_values_on_face(i_child, k) *
504  this->shape_value_component(
505  face * this->n_dofs_per_face(face) + i_face,
506  q_sub.point(k),
508  }
509  }
510  }
511 
512  if (this->degree == 1)
513  return;
514 
515  // Create Legendre basis for the space D_xi Q_k. Here, we cannot
516  // use the shape functions
517  std::unique_ptr<AnisotropicPolynomials<dim>> polynomials[dim];
518  for (unsigned int dd = 0; dd < dim; ++dd)
519  {
520  std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
521  for (unsigned int d = 0; d < dim; ++d)
522  poly[d] =
524  poly[dd] =
526 
527  polynomials[dd] = std::make_unique<AnisotropicPolynomials<dim>>(poly);
528  }
529 
530  // TODO: the implementation makes the assumption that all faces have the
531  // same number of dofs
532  AssertDimension(this->n_unique_faces(), 1);
533  const unsigned int face_no = 0;
534 
535  QGauss<dim> q_cell(this->degree);
536  const unsigned int start_cell_dofs =
538 
539  // Store shape values, since the
540  // evaluation suffers if not
541  // ordered by point
542  Table<3, double> cached_values_on_cell(this->n_dofs_per_cell(),
543  q_cell.size(),
544  dim);
545  for (unsigned int k = 0; k < q_cell.size(); ++k)
546  for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
547  for (unsigned int d = 0; d < dim; ++d)
548  cached_values_on_cell(i, k, d) =
549  this->shape_value_component(i, q_cell.point(k), d);
550 
551  for (unsigned int child = 0; child < GeometryInfo<dim>::max_children_per_cell;
552  ++child)
553  {
554  Quadrature<dim> q_sub =
556  q_cell,
557  child);
558 
559  for (unsigned int k = 0; k < q_sub.size(); ++k)
560  for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
561  ++i_child)
562  for (unsigned int d = 0; d < dim; ++d)
563  for (unsigned int i_weight = 0; i_weight < polynomials[d]->n();
564  ++i_weight)
565  {
566  this->restriction[iso][child](start_cell_dofs + i_weight * dim +
567  d,
568  i_child) +=
569  q_sub.weight(k) * cached_values_on_cell(i_child, k, d) *
570  polynomials[d]->compute_value(i_weight, q_sub.point(k));
571  }
572  }
573 }
574 
575 
576 
577 template <int dim>
578 std::vector<unsigned int>
580 {
581  // the element is face-based and we have
582  // (deg+1)^(dim-1) DoFs per face
583  unsigned int dofs_per_face = 1;
584  for (unsigned int d = 1; d < dim; ++d)
585  dofs_per_face *= deg + 1;
586 
587  // and then there are interior dofs
588  const unsigned int interior_dofs = dim * deg * dofs_per_face;
589 
590  std::vector<unsigned int> dpo(dim + 1);
591  dpo[dim - 1] = dofs_per_face;
592  dpo[dim] = interior_dofs;
593 
594  return dpo;
595 }
596 
597 
598 
599 template <int dim>
600 std::pair<Table<2, bool>, std::vector<unsigned int>>
602 {
603  Table<2, bool> constant_modes(dim, this->n_dofs_per_cell());
604  for (unsigned int d = 0; d < dim; ++d)
605  for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
606  constant_modes(d, i) = true;
607  std::vector<unsigned int> components;
608  for (unsigned int d = 0; d < dim; ++d)
609  components.push_back(d);
610  return std::pair<Table<2, bool>, std::vector<unsigned int>>(constant_modes,
611  components);
612 }
613 
614 
615 
616 //---------------------------------------------------------------------------
617 // Data field initialization
618 //---------------------------------------------------------------------------
619 
620 
621 template <int dim>
622 bool
623 FE_RaviartThomas<dim>::has_support_on_face(const unsigned int shape_index,
624  const unsigned int face_index) const
625 {
626  AssertIndexRange(shape_index, this->n_dofs_per_cell());
628 
629  // Return computed values if we
630  // know them easily. Otherwise, it
631  // is always safe to return true.
632  switch (this->degree)
633  {
634  case 1:
635  {
636  switch (dim)
637  {
638  case 2:
639  {
640  // only on the one
641  // non-adjacent face
642  // are the values
643  // actually zero. list
644  // these in a table
645  return (face_index !=
646  GeometryInfo<dim>::opposite_face[shape_index]);
647  }
648 
649  default:
650  return true;
651  }
652  }
653 
654  default: // other rt_order
655  return true;
656  }
657 
658  return true;
659 }
660 
661 
662 
663 template <int dim>
664 void
666  const std::vector<Vector<double>> &support_point_values,
667  std::vector<double> & nodal_values) const
668 {
669  Assert(support_point_values.size() == this->generalized_support_points.size(),
670  ExcDimensionMismatch(support_point_values.size(),
671  this->generalized_support_points.size()));
672  Assert(nodal_values.size() == this->n_dofs_per_cell(),
673  ExcDimensionMismatch(nodal_values.size(), this->n_dofs_per_cell()));
674  Assert(support_point_values[0].size() == this->n_components(),
675  ExcDimensionMismatch(support_point_values[0].size(),
676  this->n_components()));
677 
678  std::fill(nodal_values.begin(), nodal_values.end(), 0.);
679 
680  const unsigned int n_face_points = boundary_weights.size(0);
681  for (unsigned int face : GeometryInfo<dim>::face_indices())
682  for (unsigned int k = 0; k < n_face_points; ++k)
683  for (unsigned int i = 0; i < boundary_weights.size(1); ++i)
684  {
685  nodal_values[i + face * this->n_dofs_per_face(face)] +=
686  boundary_weights(k, i) *
687  support_point_values[face * n_face_points + k](
689  }
690 
691  // TODO: the implementation makes the assumption that all faces have the
692  // same number of dofs
693  AssertDimension(this->n_unique_faces(), 1);
694  const unsigned int face_no = 0;
695 
696  const unsigned int start_cell_dofs =
698  const unsigned int start_cell_points =
699  GeometryInfo<dim>::faces_per_cell * n_face_points;
700 
701  for (unsigned int k = 0; k < interior_weights.size(0); ++k)
702  for (unsigned int i = 0; i < interior_weights.size(1); ++i)
703  for (unsigned int d = 0; d < dim; ++d)
704  nodal_values[start_cell_dofs + i * dim + d] +=
705  interior_weights(k, i, d) *
706  support_point_values[k + start_cell_points](d);
707 }
708 
709 
710 
711 template <int dim>
712 std::size_t
714 {
715  Assert(false, ExcNotImplemented());
716  return 0;
717 }
718 
719 
720 
721 // explicit instantiations
722 #include "fe_raviart_thomas.inst"
723 
724 
void reinit(MatrixBlock< MatrixType > &v, const BlockSparsityPattern &p)
Definition: matrix_block.h:618
static Quadrature< dim > project_to_child(const Quadrature< dim > &quadrature, const unsigned int child_no)
Definition: qprojector.cc:1238
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2404
const unsigned int components
Definition: fe_base.h:430
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1653
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2455
void compute_embedding_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number >>> &matrices, const bool isotropic_only=false, const double threshold=1.e-12)
FullMatrix< double > interface_constraints
Definition: fe.h:2430
static void project_to_subface(const SubQuadrature &quadrature, const unsigned int face_no, const unsigned int subface_no, std::vector< Point< dim >> &q_points, const RefinementCase< dim - 1 > &ref_case=RefinementCase< dim - 1 >::isotropic_refinement)
static Quadrature< dim > project_to_all_faces(const Quadrature< dim - 1 > &quadrature)
Definition: qprojector.h:579
static std::vector< Polynomial< double > > generate_complete_basis(const unsigned int degree)
Definition: polynomial.cc:744
#define AssertIndexRange(index, range)
Definition: exceptions.h:1718
void initialize_quad_dof_index_permutation_and_sign_change()
const unsigned int degree
Definition: fe_base.h:436
const Point< dim > & point(const unsigned int i) const
void invert(const FullMatrix< number2 > &M)
STL namespace.
unsigned int n_dofs_per_quad(unsigned int face_no=0) const
FullMatrix< double > inverse_node_matrix
FullMatrix< double > compute_node_matrix(const FiniteElement< dim, spacedim > &fe)
virtual std::size_t memory_consumption() const override
virtual std::string get_name() const override
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
friend class FE_RaviartThomas
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2418
#define Assert(cond, exc)
Definition: exceptions.h:1461
void initialize_support_points(const unsigned int rt_degree)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition: fe.h:2461
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_unique_faces() const
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
std::vector< Table< 2, bool > > adjust_quad_dof_sign_for_face_orientation_table
unsigned int size() const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const override
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
static void project_to_face(const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim >> &q_points)
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 bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement)
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
const unsigned int dofs_per_face
Definition: fe_base.h:406
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:303
static ::ExceptionBase & ExcNotImplemented()
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const override
std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
double legendre(unsigned int l, double x)
Definition: cmath.h:75
Table< 3, double > interior_weights
bool get_bit(const unsigned char number, const unsigned int n)
Definition: utilities.h:1380
double weight(const unsigned int i) const
ReferenceCell reference_cell() const
unsigned int tensor_degree() const
Table< 2, double > boundary_weights
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2478
static ::ExceptionBase & ExcInternalError()
std::vector< MappingKind > mapping_kind