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