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