Reference documentation for deal.II version GIT relicensing-426-g7976cfd195 2024-04-18 21:10:01+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\}}\)
Loading...
Searching...
No Matches
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
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
122template <int dim>
123std::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
144template <int dim>
145std::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
157template <int dim>
158void
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);
185 TensorProductPolynomials<dim - 1> legendre =
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 =
206 QProjector<dim>::project_to_all_faces(this->reference_cell(),
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(
215 this->reference_cell(),
216 0,
218 n_face_points));
219 }
220 }
221
222 if (deg == 0)
223 return;
224
225 // Create Legendre basis for the space D_xi Q_k
226 std::unique_ptr<AnisotropicPolynomials<dim>> polynomials[dim];
227 for (unsigned int dd = 0; dd < dim; ++dd)
228 {
229 std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
230 for (unsigned int d = 0; d < dim; ++d)
233
234 polynomials[dd] = std::make_unique<AnisotropicPolynomials<dim>>(poly);
235 }
236
237 interior_weights.reinit(
238 TableIndices<3>(n_interior_points, polynomials[0]->n(), dim));
239
240 for (unsigned int k = 0; k < cell_quadrature.size(); ++k)
241 {
242 this->generalized_support_points[current++] = cell_quadrature.point(k);
243 for (unsigned int i = 0; i < polynomials[0]->n(); ++i)
244 for (unsigned int d = 0; d < dim; ++d)
245 interior_weights(k, i, d) =
246 cell_quadrature.weight(k) *
247 polynomials[d]->compute_value(i, cell_quadrature.point(k));
248 }
249
250 Assert(current == this->generalized_support_points.size(),
252}
253
254
255template <int dim>
256void
258{
259 // For 1d do nothing.
260 //
261 // TODO: For 2d we simply keep the legacy behavior for now. This should be
262 // changed in the future and can be taken care of by similar means as the 3d
263 // case below. The legacy behavior can be found in fe_poly_tensor.cc in the
264 // function internal::FE_PolyTensor::get_dof_sign_change_h_div(...)
265 if (dim < 3)
266 return;
267
268 // TODO: the implementation makes the assumption that all faces have the
269 // same number of dofs
270 AssertDimension(this->n_unique_faces(), 1);
271 const unsigned int face_no = 0;
272
273 Assert(
274 this->adjust_quad_dof_index_for_face_orientation_table[0].n_elements() ==
275 this->reference_cell().n_face_orientations(face_no) *
276 this->n_dofs_per_quad(face_no),
278
279 // The 3d RaviartThomas space has tensor_degree*tensor_degree face dofs
280 const unsigned int n = this->tensor_degree();
281 Assert(n * n == this->n_dofs_per_quad(face_no), ExcInternalError());
282
283 // The vector of tables adjust_quad_dof_index_for_face_orientation_table
284 // contains offsets for local face_dofs and is being filled with zeros in
285 // fe.cc. We need to fill it with the correct values in case of non-standard,
286 // flipped (rotated by +180 degrees) or rotated (rotated by +90 degrees) faces
287
288 // The dofs on a face are connected to a n x n matrix. for example, for
289 // tensor_degree==3 we have the following dofs on a quad:
290
291 // ___________
292 // | |
293 // | 6 7 8 |
294 // | |
295 // | 3 4 5 |
296 // | |
297 // | 0 1 2 |
298 // |___________|
299 //
300 // We have dof_index=i+n*j with index i in x-direction and index j in
301 // y-direction running from 0 to n-1. to extract i and j we can use
302 // i=dof_index%n and j=dof_index/n. The indices i and j can then be used to
303 // compute the offset.
304
305 // Example: if the switches are (true | true | true) that means we rotate the
306 // face first by + 90 degree(counterclockwise) then by another +180
307 // degrees but we do not flip it since the face has standard
308 // orientation. The flip axis is the diagonal from the lower left to the upper
309 // right corner of the face. With these flags the configuration above becomes:
310 // ___________
311 // | |
312 // | 2 5 8 |
313 // | |
314 // | 1 4 7 |
315 // | |
316 // | 0 3 6 |
317 // |___________|
318 //
319 // This is exactly what is realized by the formulas implemented below. Note
320 // that the necessity of a permuattion depends on the three flags.
321
322 // There is also a pattern for the sign change of the mapped face_dofs
323 // depending on the switches. In the above example it would be
324 // ___________
325 // | |
326 // | + - + |
327 // | |
328 // | + - + |
329 // | |
330 // | + - + |
331 // |___________|
332 //
333
334 for (unsigned int local_face_dof = 0;
335 local_face_dof < this->n_dofs_per_quad(face_no);
336 ++local_face_dof)
337 {
338 // Row and column
339 unsigned int i = local_face_dof % n;
340 unsigned int j = local_face_dof / n;
341
342 // We have 8 cases that are all treated the same way. Note that the
343 // corresponding case to case_no is just its binary representation.
344 // The above example of (false | true | true) would be case_no=3
345 for (const bool face_orientation : {false, true})
346 for (const bool face_flip : {false, true})
347 for (const bool face_rotation : {false, true})
348 {
349 const auto case_no =
351 face_rotation,
352 face_flip);
353
354 if (((!face_orientation) && (!face_rotation)) ||
355 ((face_orientation) && (face_rotation)))
356 {
357 // We flip across the diagonal
358 // This is the local face dof offset
359 this
360 ->adjust_quad_dof_index_for_face_orientation_table[face_no](
361 local_face_dof, case_no) = j + i * n - local_face_dof;
362 }
363 else
364 {
365 // Offset is zero
366 this
367 ->adjust_quad_dof_index_for_face_orientation_table[face_no](
368 local_face_dof, case_no) = 0;
369 } // if face needs dof permutation
370
371 // Get new local face_dof by adding offset
372 const unsigned int new_local_face_dof =
373 local_face_dof +
374 this->adjust_quad_dof_index_for_face_orientation_table[face_no](
375 local_face_dof, case_no);
376 // compute new row and column index
377 i = new_local_face_dof % n;
378 j = new_local_face_dof / n;
379
380 /*
381 * Now compute if a sign change is necessary. This is done for the
382 * case of face_orientation==true
383 */
384 // flip = false, rotation=true
385 if (!face_flip && face_rotation)
386 {
387 this
388 ->adjust_quad_dof_sign_for_face_orientation_table[face_no](
389 local_face_dof, case_no) = ((j % 2) == 1);
390 }
391 // flip = true, rotation=false
392 else if (face_flip && !face_rotation)
393 {
394 // This case is symmetric (although row and column may be
395 // switched)
396 this
397 ->adjust_quad_dof_sign_for_face_orientation_table[face_no](
398 local_face_dof, case_no) =
399 ((j % 2) == 1) != ((i % 2) == 1);
400 }
401 // flip = true, rotation=true
402 else if (face_flip && face_rotation)
403 {
404 this
405 ->adjust_quad_dof_sign_for_face_orientation_table[face_no](
406 local_face_dof, case_no) = ((i % 2) == 1);
407 }
408 /*
409 * flip = false, rotation=false => nothing to do
410 */
411
412 /*
413 * If face_orientation==false the sign flip is exactly the
414 * opposite.
415 */
416 if (!face_orientation)
417 this->adjust_quad_dof_sign_for_face_orientation_table[face_no](
418 local_face_dof, case_no) =
419 !this
420 ->adjust_quad_dof_sign_for_face_orientation_table[face_no](
421 local_face_dof, case_no);
422 } // case_no
423 } // local_face_dof
424}
425
426
427template <>
428void
430{
431 // there is only one refinement case in 1d,
432 // which is the isotropic one (first index of
433 // the matrix array has to be 0)
434 for (unsigned int i = 0; i < GeometryInfo<1>::max_children_per_cell; ++i)
435 this->restriction[0][i].reinit(0, 0);
436}
437
438
439
440// This function is the same Raviart-Thomas interpolation performed by
441// interpolate. Still, we cannot use interpolate, since it was written
442// for smooth functions. The functions interpolated here are not
443// smooth, maybe even not continuous. Therefore, we must double the
444// number of quadrature points in each direction in order to integrate
445// only smooth functions.
446
447// Then again, the interpolated function is chosen such that the
448// moments coincide with the function to be interpolated.
449
450template <int dim>
451void
453{
454 const unsigned int iso = RefinementCase<dim>::isotropic_refinement - 1;
455
456 QGauss<dim - 1> q_base(this->degree);
457 const unsigned int n_face_points = q_base.size();
458 // First, compute interpolation on
459 // subfaces
460 for (const unsigned int face : GeometryInfo<dim>::face_indices())
461 {
462 // The shape functions of the
463 // child cell are evaluated
464 // in the quadrature points
465 // of a full face.
466 Quadrature<dim> q_face =
467 QProjector<dim>::project_to_face(this->reference_cell(), q_base, face);
468 // Store shape values, since the
469 // evaluation suffers if not
470 // ordered by point
471 Table<2, double> cached_values_on_face(this->n_dofs_per_cell(),
472 q_face.size());
473 for (unsigned int k = 0; k < q_face.size(); ++k)
474 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
475 cached_values_on_face(i, k) = this->shape_value_component(
477
478 for (unsigned int sub = 0; sub < GeometryInfo<dim>::max_children_per_face;
479 ++sub)
480 {
481 // The weight functions for
482 // the coarse face are
483 // evaluated on the subface
484 // only.
486 this->reference_cell(), q_base, face, sub);
487 const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
489
490 // On a certain face, we must
491 // compute the moments of ALL
492 // fine level functions with
493 // the coarse level weight
494 // functions belonging to
495 // that face. Due to the
496 // orthogonalization process
497 // when building the shape
498 // functions, these weights
499 // are equal to the
500 // corresponding shape
501 // functions.
502 for (unsigned int k = 0; k < n_face_points; ++k)
503 for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
504 ++i_child)
505 for (unsigned int i_face = 0;
506 i_face < this->n_dofs_per_face(face);
507 ++i_face)
508 {
509 // The quadrature
510 // weights on the
511 // subcell are NOT
512 // transformed, so we
513 // have to do it here.
514 this->restriction[iso][child](
515 face * this->n_dofs_per_face(face) + i_face, i_child) +=
516 Utilities::fixed_power<dim - 1>(.5) * q_sub.weight(k) *
517 cached_values_on_face(i_child, k) *
518 this->shape_value_component(
519 face * this->n_dofs_per_face(face) + i_face,
520 q_sub.point(k),
522 }
523 }
524 }
525
526 if (this->degree == 1)
527 return;
528
529 // Create Legendre basis for the space D_xi Q_k. Here, we cannot
530 // use the shape functions
531 std::unique_ptr<AnisotropicPolynomials<dim>> polynomials[dim];
532 for (unsigned int dd = 0; dd < dim; ++dd)
533 {
534 std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
535 for (unsigned int d = 0; d < dim; ++d)
536 poly[d] =
538 poly[dd] =
540
541 polynomials[dd] = std::make_unique<AnisotropicPolynomials<dim>>(poly);
542 }
543
544 // TODO: the implementation makes the assumption that all faces have the
545 // same number of dofs
546 AssertDimension(this->n_unique_faces(), 1);
547 const unsigned int face_no = 0;
548
549 QGauss<dim> q_cell(this->degree);
550 const unsigned int start_cell_dofs =
551 GeometryInfo<dim>::faces_per_cell * this->n_dofs_per_face(face_no);
552
553 // Store shape values, since the
554 // evaluation suffers if not
555 // ordered by point
556 Table<3, double> cached_values_on_cell(this->n_dofs_per_cell(),
557 q_cell.size(),
558 dim);
559 for (unsigned int k = 0; k < q_cell.size(); ++k)
560 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
561 for (unsigned int d = 0; d < dim; ++d)
562 cached_values_on_cell(i, k, d) =
563 this->shape_value_component(i, q_cell.point(k), d);
564
565 for (unsigned int child = 0; child < GeometryInfo<dim>::max_children_per_cell;
566 ++child)
567 {
568 Quadrature<dim> q_sub =
569 QProjector<dim>::project_to_child(this->reference_cell(),
570 q_cell,
571 child);
572
573 for (unsigned int k = 0; k < q_sub.size(); ++k)
574 for (unsigned int i_child = 0; i_child < this->n_dofs_per_cell();
575 ++i_child)
576 for (unsigned int d = 0; d < dim; ++d)
577 for (unsigned int i_weight = 0; i_weight < polynomials[d]->n();
578 ++i_weight)
579 {
580 this->restriction[iso][child](start_cell_dofs + i_weight * dim +
581 d,
582 i_child) +=
583 q_sub.weight(k) * cached_values_on_cell(i_child, k, d) *
584 polynomials[d]->compute_value(i_weight, q_sub.point(k));
585 }
586 }
587}
588
589
590
591template <int dim>
592std::vector<unsigned int>
594{
595 // the element is face-based and we have
596 // (deg+1)^(dim-1) DoFs per face
597 unsigned int dofs_per_face = 1;
598 for (unsigned int d = 1; d < dim; ++d)
599 dofs_per_face *= deg + 1;
600
601 // and then there are interior dofs
602 const unsigned int interior_dofs = dim * deg * dofs_per_face;
603
604 std::vector<unsigned int> dpo(dim + 1);
605 dpo[dim - 1] = dofs_per_face;
606 dpo[dim] = interior_dofs;
607
608 return dpo;
609}
610
611
612
613template <int dim>
614std::pair<Table<2, bool>, std::vector<unsigned int>>
616{
617 Table<2, bool> constant_modes(dim, this->n_dofs_per_cell());
618 for (unsigned int d = 0; d < dim; ++d)
619 for (unsigned int i = 0; i < this->n_dofs_per_cell(); ++i)
620 constant_modes(d, i) = true;
621 std::vector<unsigned int> components;
622 for (unsigned int d = 0; d < dim; ++d)
623 components.push_back(d);
624 return std::pair<Table<2, bool>, std::vector<unsigned int>>(constant_modes,
625 components);
626}
627
628
629
630//---------------------------------------------------------------------------
631// Data field initialization
632//---------------------------------------------------------------------------
633
634
635template <int dim>
636bool
637FE_RaviartThomas<dim>::has_support_on_face(const unsigned int shape_index,
638 const unsigned int face_index) const
639{
640 AssertIndexRange(shape_index, this->n_dofs_per_cell());
642
643 // Return computed values if we
644 // know them easily. Otherwise, it
645 // is always safe to return true.
646 switch (this->degree)
647 {
648 case 1:
649 {
650 switch (dim)
651 {
652 case 2:
653 {
654 // only on the one
655 // non-adjacent face
656 // are the values
657 // actually zero. list
658 // these in a table
659 return (face_index !=
661 }
662
663 default:
664 return true;
665 }
666 }
667
668 default: // other rt_order
669 return true;
670 }
671
672 return true;
673}
674
675
676
677template <int dim>
678void
680 const std::vector<Vector<double>> &support_point_values,
681 std::vector<double> &nodal_values) const
682{
683 Assert(support_point_values.size() == this->generalized_support_points.size(),
684 ExcDimensionMismatch(support_point_values.size(),
685 this->generalized_support_points.size()));
686 Assert(nodal_values.size() == this->n_dofs_per_cell(),
687 ExcDimensionMismatch(nodal_values.size(), this->n_dofs_per_cell()));
688 Assert(support_point_values[0].size() == this->n_components(),
689 ExcDimensionMismatch(support_point_values[0].size(),
690 this->n_components()));
691
692 std::fill(nodal_values.begin(), nodal_values.end(), 0.);
693
694 const unsigned int n_face_points = boundary_weights.size(0);
695 for (const unsigned int face : GeometryInfo<dim>::face_indices())
696 for (unsigned int k = 0; k < n_face_points; ++k)
697 for (unsigned int i = 0; i < boundary_weights.size(1); ++i)
698 {
699 nodal_values[i + face * this->n_dofs_per_face(face)] +=
700 boundary_weights(k, i) *
701 support_point_values[face * n_face_points + k](
703 }
704
705 // TODO: the implementation makes the assumption that all faces have the
706 // same number of dofs
707 AssertDimension(this->n_unique_faces(), 1);
708 const unsigned int face_no = 0;
709
710 const unsigned int start_cell_dofs =
711 GeometryInfo<dim>::faces_per_cell * this->n_dofs_per_face(face_no);
712 const unsigned int start_cell_points =
713 GeometryInfo<dim>::faces_per_cell * n_face_points;
714
715 for (unsigned int k = 0; k < interior_weights.size(0); ++k)
716 for (unsigned int i = 0; i < interior_weights.size(1); ++i)
717 for (unsigned int d = 0; d < dim; ++d)
718 nodal_values[start_cell_dofs + i * dim + d] +=
719 interior_weights(k, i, d) *
720 support_point_values[k + start_cell_points](d);
721}
722
723
724
725template <int dim>
726std::size_t
732
733
734
735// explicit instantiations
736#include "fe_raviart_thomas.inst"
737
738
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
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:2424
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition fe.h:2412
size_type n() const
void invert(const FullMatrix< number2 > &M)
size_type m() const
static std::vector< Polynomial< double > > generate_complete_basis(const unsigned int degree)
static Quadrature< dim > project_to_child(const ReferenceCell &reference_cell, const Quadrature< dim > &quadrature, const unsigned int child_no)
static Quadrature< dim > project_to_all_faces(const ReferenceCell &reference_cell, const hp::QCollection< dim - 1 > &quadrature)
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
static constexpr unsigned char default_combined_face_orientation()
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#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:134
#define DEAL_II_NOT_IMPLEMENTED()
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 > compute_node_matrix(const FiniteElement< dim, spacedim > &fe)
unsigned char combined_face_orientation(const bool face_orientation, const bool face_rotation, const bool face_flip)
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)
static std_cxx20::ranges::iota_view< unsigned int, unsigned int > face_indices()