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