Reference documentation for deal.II version GIT c9976103bc 2022-12-09 17:30:02+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\}}\)
fe.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2022 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 
19 
21 
22 #include <deal.II/fe/fe.h>
23 #include <deal.II/fe/fe_values.h>
24 #include <deal.II/fe/mapping.h>
25 
26 #include <deal.II/grid/tria.h>
28 
29 #include <algorithm>
30 #include <functional>
31 #include <numeric>
32 #include <typeinfo>
33 
35 
36 
37 /*------------------------------- FiniteElement ----------------------*/
38 
39 
40 template <int dim, int spacedim>
42  : update_each(update_default)
43 {}
44 
45 
46 
47 template <int dim, int spacedim>
48 std::size_t
50 {
51  return sizeof(*this);
52 }
53 
54 
55 
56 template <int dim, int spacedim>
58  const FiniteElementData<dim> & fe_data,
59  const std::vector<bool> & r_i_a_f,
60  const std::vector<ComponentMask> &nonzero_c)
61  : FiniteElementData<dim>(fe_data)
62  , adjust_line_dof_index_for_line_orientation_table(
63  dim == 3 ? this->n_dofs_per_line() : 0)
64  , system_to_base_table(this->n_dofs_per_cell())
65  , component_to_base_table(this->components,
66  std::make_pair(std::make_pair(0U, 0U), 0U))
67  ,
68 
69  // Special handling of vectors of length one: in this case, we
70  // assume that all entries were supposed to be equal
71  restriction_is_additive_flags(
72  r_i_a_f.size() == 1 ?
73  std::vector<bool>(fe_data.n_dofs_per_cell(), r_i_a_f[0]) :
74  r_i_a_f)
75  , nonzero_components(
76  nonzero_c.size() == 1 ?
77  std::vector<ComponentMask>(fe_data.n_dofs_per_cell(), nonzero_c[0]) :
78  nonzero_c)
79  , n_nonzero_components_table(compute_n_nonzero_components(nonzero_components))
80  , cached_primitivity(std::find_if(n_nonzero_components_table.begin(),
81  n_nonzero_components_table.end(),
82  [](const unsigned int n_components) {
83  return n_components != 1U;
84  }) == n_nonzero_components_table.end())
85 {
86  Assert(restriction_is_additive_flags.size() == this->n_dofs_per_cell(),
87  ExcDimensionMismatch(restriction_is_additive_flags.size(),
88  this->n_dofs_per_cell()));
89  AssertDimension(nonzero_components.size(), this->n_dofs_per_cell());
90  for (unsigned int i = 0; i < nonzero_components.size(); ++i)
91  {
92  Assert(nonzero_components[i].size() == this->n_components(),
94  Assert(nonzero_components[i].n_selected_components() >= 1,
96  Assert(n_nonzero_components_table[i] >= 1, ExcInternalError());
97  Assert(n_nonzero_components_table[i] <= this->n_components(),
99  }
100 
101  // initialize some tables in the default way, i.e. if there is only one
102  // (vector-)component; if the element is not primitive, leave these tables
103  // empty.
104  if (this->is_primitive())
105  {
106  system_to_component_table.resize(this->n_dofs_per_cell());
107  for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
108  system_to_component_table[j] = std::pair<unsigned, unsigned>(0, j);
109 
110  face_system_to_component_table.resize(this->n_unique_faces());
111  for (unsigned int f = 0; f < this->n_unique_faces(); ++f)
112  {
113  face_system_to_component_table[f].resize(this->n_dofs_per_face(f));
114  for (unsigned int j = 0; j < this->n_dofs_per_face(f); ++j)
115  face_system_to_component_table[f][j] =
116  std::pair<unsigned, unsigned>(0, j);
117  }
118  }
119 
120  for (unsigned int j = 0; j < this->n_dofs_per_cell(); ++j)
121  system_to_base_table[j] = std::make_pair(std::make_pair(0U, 0U), j);
122 
123  face_system_to_base_table.resize(this->n_unique_faces());
124  for (unsigned int f = 0; f < this->n_unique_faces(); ++f)
125  {
126  face_system_to_base_table[f].resize(this->n_dofs_per_face(f));
127  for (unsigned int j = 0; j < this->n_dofs_per_face(f); ++j)
128  face_system_to_base_table[f][j] =
129  std::make_pair(std::make_pair(0U, 0U), j);
130  }
131 
132  // Fill with default value; may be changed by constructor of derived class.
133  base_to_block_indices.reinit(1, 1);
134 
135  // initialize the restriction and prolongation matrices. the default
136  // constructor of FullMatrix<dim> initializes them with size zero
137  prolongation.resize(RefinementCase<dim>::isotropic_refinement);
138  restriction.resize(RefinementCase<dim>::isotropic_refinement);
139  for (unsigned int ref = RefinementCase<dim>::cut_x;
140  ref < RefinementCase<dim>::isotropic_refinement + 1;
141  ++ref)
142  {
143  prolongation[ref - 1].resize(GeometryInfo<dim>::n_children(
144  RefinementCase<dim>(ref)),
146  restriction[ref - 1].resize(GeometryInfo<dim>::n_children(
147  RefinementCase<dim>(ref)),
149  }
150 
151 
152  if (dim == 3)
153  {
154  adjust_quad_dof_index_for_face_orientation_table.resize(
155  this->n_unique_quads());
156 
157  for (unsigned int f = 0; f < this->n_unique_quads(); ++f)
158  {
159  adjust_quad_dof_index_for_face_orientation_table[f] =
160  Table<2, int>(this->n_dofs_per_quad(f),
161  this->reference_cell().n_face_orientations(f));
162  adjust_quad_dof_index_for_face_orientation_table[f].fill(0);
163  }
164  }
165 
166  unit_face_support_points.resize(this->n_unique_faces());
167  generalized_face_support_points.resize(this->n_unique_faces());
168 }
169 
170 
171 
172 template <int dim, int spacedim>
173 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
174 FiniteElement<dim, spacedim>::operator^(const unsigned int multiplicity) const
175 {
176  return {this->clone(), multiplicity};
177 }
178 
179 
180 
181 template <int dim, int spacedim>
182 double
184  const Point<dim> &) const
185 {
187  return 0.;
188 }
189 
190 
191 
192 template <int dim, int spacedim>
193 double
195  const Point<dim> &,
196  const unsigned int) const
197 {
199  return 0.;
200 }
201 
202 
203 
204 template <int dim, int spacedim>
207  const Point<dim> &) const
208 {
210  return Tensor<1, dim>();
211 }
212 
213 
214 
215 template <int dim, int spacedim>
218  const Point<dim> &,
219  const unsigned int) const
220 {
222  return Tensor<1, dim>();
223 }
224 
225 
226 
227 template <int dim, int spacedim>
230  const Point<dim> &) const
231 {
233  return Tensor<2, dim>();
234 }
235 
236 
237 
238 template <int dim, int spacedim>
241  const unsigned int,
242  const Point<dim> &,
243  const unsigned int) const
244 {
246  return Tensor<2, dim>();
247 }
248 
249 
250 
251 template <int dim, int spacedim>
254  const Point<dim> &) const
255 {
257  return Tensor<3, dim>();
258 }
259 
260 
261 
262 template <int dim, int spacedim>
265  const unsigned int,
266  const Point<dim> &,
267  const unsigned int) const
268 {
270  return Tensor<3, dim>();
271 }
272 
273 
274 
275 template <int dim, int spacedim>
278  const Point<dim> &) const
279 {
281  return Tensor<4, dim>();
282 }
283 
284 
285 
286 template <int dim, int spacedim>
289  const unsigned int,
290  const Point<dim> &,
291  const unsigned int) const
292 {
294  return Tensor<4, dim>();
295 }
296 
297 
298 template <int dim, int spacedim>
299 void
301  const bool isotropic_restriction_only,
302  const bool isotropic_prolongation_only)
303 {
304  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
305  ref_case <= RefinementCase<dim>::isotropic_refinement;
306  ++ref_case)
307  {
308  const unsigned int nc =
310 
311  for (unsigned int i = 0; i < nc; ++i)
312  {
313  if (this->restriction[ref_case - 1][i].m() !=
314  this->n_dofs_per_cell() &&
315  (!isotropic_restriction_only ||
317  this->restriction[ref_case - 1][i].reinit(this->n_dofs_per_cell(),
318  this->n_dofs_per_cell());
319  if (this->prolongation[ref_case - 1][i].m() !=
320  this->n_dofs_per_cell() &&
321  (!isotropic_prolongation_only ||
323  this->prolongation[ref_case - 1][i].reinit(this->n_dofs_per_cell(),
324  this->n_dofs_per_cell());
325  }
326  }
327 }
328 
329 
330 template <int dim, int spacedim>
331 const FullMatrix<double> &
333  const unsigned int child,
334  const RefinementCase<dim> &refinement_case) const
335 {
336  AssertIndexRange(refinement_case,
338  Assert(refinement_case != RefinementCase<dim>::no_refinement,
339  ExcMessage(
340  "Restriction matrices are only available for refined cells!"));
342  child, GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)));
343  // we use refinement_case-1 here. the -1 takes care of the origin of the
344  // vector, as for RefinementCase<dim>::no_refinement (=0) there is no data
345  // available and so the vector indices are shifted
346  Assert(restriction[refinement_case - 1][child].n() == this->n_dofs_per_cell(),
348  return restriction[refinement_case - 1][child];
349 }
350 
351 
352 
353 template <int dim, int spacedim>
354 const FullMatrix<double> &
356  const unsigned int child,
357  const RefinementCase<dim> &refinement_case) const
358 {
359  AssertIndexRange(refinement_case,
361  Assert(refinement_case != RefinementCase<dim>::no_refinement,
362  ExcMessage(
363  "Prolongation matrices are only available for refined cells!"));
365  child, GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)));
366  // we use refinement_case-1 here. the -1 takes care
367  // of the origin of the vector, as for
368  // RefinementCase::no_refinement (=0) there is no
369  // data available and so the vector indices
370  // are shifted
371  Assert(prolongation[refinement_case - 1][child].n() ==
372  this->n_dofs_per_cell(),
373  ExcEmbeddingVoid());
374  return prolongation[refinement_case - 1][child];
375 }
376 
377 
378 // TODO:[GK] This is probably not the most efficient way of doing this.
379 template <int dim, int spacedim>
380 unsigned int
382  const unsigned int index) const
383 {
384  AssertIndexRange(index, this->n_components());
385 
386  return first_block_of_base(component_to_base_table[index].first.first) +
388 }
389 
390 
391 template <int dim, int spacedim>
394  const FEValuesExtractors::Scalar &scalar) const
395 {
396  AssertIndexRange(scalar.component, this->n_components());
397 
398  // TODO: it would be nice to verify that it is indeed possible
399  // to select this scalar component, i.e., that it is not part
400  // of a non-primitive element. unfortunately, there is no simple
401  // way to write such a condition...
402 
403  std::vector<bool> mask(this->n_components(), false);
404  mask[scalar.component] = true;
405  return mask;
406 }
407 
408 
409 template <int dim, int spacedim>
412  const FEValuesExtractors::Vector &vector) const
413 {
414  AssertIndexRange(vector.first_vector_component + dim - 1,
415  this->n_components());
416 
417  // TODO: it would be nice to verify that it is indeed possible
418  // to select these vector components, i.e., that they don't span
419  // beyond the beginning or end of anon-primitive element.
420  // unfortunately, there is no simple way to write such a condition...
421 
422  std::vector<bool> mask(this->n_components(), false);
423  for (unsigned int c = vector.first_vector_component;
424  c < vector.first_vector_component + dim;
425  ++c)
426  mask[c] = true;
427  return mask;
428 }
429 
430 
431 template <int dim, int spacedim>
434  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
435 {
438  this->n_components());
439 
440  // TODO: it would be nice to verify that it is indeed possible
441  // to select these vector components, i.e., that they don't span
442  // beyond the beginning or end of anon-primitive element.
443  // unfortunately, there is no simple way to write such a condition...
444 
445  std::vector<bool> mask(this->n_components(), false);
446  for (unsigned int c = sym_tensor.first_tensor_component;
447  c < sym_tensor.first_tensor_component +
449  ++c)
450  mask[c] = true;
451  return mask;
452 }
453 
454 
455 
456 template <int dim, int spacedim>
459 {
460  // if we get a block mask that represents all blocks, then
461  // do the same for the returned component mask
463  return {};
464 
465  AssertDimension(block_mask.size(), this->n_blocks());
466 
467  std::vector<bool> component_mask(this->n_components(), false);
468  for (unsigned int c = 0; c < this->n_components(); ++c)
469  if (block_mask[component_to_block_index(c)] == true)
470  component_mask[c] = true;
471 
472  return component_mask;
473 }
474 
475 
476 
477 template <int dim, int spacedim>
478 BlockMask
480  const FEValuesExtractors::Scalar &scalar) const
481 {
482  // simply create the corresponding component mask (a simpler
483  // process) and then convert it to a block mask
484  return block_mask(component_mask(scalar));
485 }
486 
487 
488 template <int dim, int spacedim>
489 BlockMask
491  const FEValuesExtractors::Vector &vector) const
492 {
493  // simply create the corresponding component mask (a simpler
494  // process) and then convert it to a block mask
495  return block_mask(component_mask(vector));
496 }
497 
498 
499 template <int dim, int spacedim>
500 BlockMask
502  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
503 {
504  // simply create the corresponding component mask (a simpler
505  // process) and then convert it to a block mask
506  return block_mask(component_mask(sym_tensor));
507 }
508 
509 
510 
511 template <int dim, int spacedim>
512 BlockMask
514  const ComponentMask &component_mask) const
515 {
516  // if we get a component mask that represents all component, then
517  // do the same for the returned block mask
519  return {};
520 
521  AssertDimension(component_mask.size(), this->n_components());
522 
523  // walk over all of the components
524  // of this finite element and see
525  // if we need to set the
526  // corresponding block. inside the
527  // block, walk over all the
528  // components that correspond to
529  // this block and make sure the
530  // component mask is set for all of
531  // them
532  std::vector<bool> block_mask(this->n_blocks(), false);
533  for (unsigned int c = 0; c < this->n_components();)
534  {
535  const unsigned int block = component_to_block_index(c);
536  if (component_mask[c] == true)
537  block_mask[block] = true;
538 
539  // now check all of the other
540  // components that correspond
541  // to this block
542  ++c;
543  while ((c < this->n_components()) &&
544  (component_to_block_index(c) == block))
545  {
546  Assert(component_mask[c] == block_mask[block],
547  ExcMessage(
548  "The component mask argument given to this function "
549  "is not a mask where the individual components belonging "
550  "to one block of the finite element are either all "
551  "selected or not selected. You can't call this function "
552  "with a component mask that splits blocks."));
553  ++c;
554  }
555  }
556 
557 
558  return block_mask;
559 }
560 
561 
562 
563 template <int dim, int spacedim>
564 unsigned int
565 FiniteElement<dim, spacedim>::face_to_cell_index(const unsigned int face_index,
566  const unsigned int face,
567  const bool face_orientation,
568  const bool face_flip,
569  const bool face_rotation) const
570 {
571  AssertIndexRange(face_index, this->n_dofs_per_face(face));
572  AssertIndexRange(face, this->reference_cell().n_faces());
573 
574  // TODO: we could presumably solve the 3d case below using the
575  // adjust_quad_dof_index_for_face_orientation_table field. for the
576  // 2d case, we can't use adjust_line_dof_index_for_line_orientation_table
577  // since that array is empty (presumably because we thought that
578  // there are no flipped edges in 2d, but these can happen in
579  // DoFTools::make_periodicity_constraints, for example). so we
580  // would need to either fill this field, or rely on derived classes
581  // implementing this function, as we currently do
582 
583  // see the function's documentation for an explanation of this
584  // assertion -- in essence, derived classes have to implement
585  // an overloaded version of this function if we are to use any
586  // other than standard orientation
587  if ((face_orientation != true) || (face_flip != false) ||
588  (face_rotation != false))
589  Assert((this->n_dofs_per_line() <= 1) && (this->n_dofs_per_quad(face) <= 1),
590  ExcMessage(
591  "The function in this base class can not handle this case. "
592  "Rather, the derived class you are using must provide "
593  "an overloaded version but apparently hasn't done so. See "
594  "the documentation of this function for more information."));
595 
596  // we need to distinguish between DoFs on vertices, lines and in 3d quads.
597  // do so in a sequence of if-else statements
598  if (face_index < this->get_first_face_line_index(face))
599  // DoF is on a vertex
600  {
601  // get the number of the vertex on the face that corresponds to this DoF,
602  // along with the number of the DoF on this vertex
603  const unsigned int face_vertex = face_index / this->n_dofs_per_vertex();
604  const unsigned int dof_index_on_vertex =
605  face_index % this->n_dofs_per_vertex();
606 
607  // then get the number of this vertex on the cell and translate
608  // this to a DoF number on the cell
609  return (this->reference_cell().face_to_cell_vertices(
610  face,
611  face_vertex,
612  (face_orientation ? 1 : 0) + (face_rotation ? 2 : 0) +
613  (face_flip ? 4 : 0)) *
614  this->n_dofs_per_vertex() +
615  dof_index_on_vertex);
616  }
617  else if (face_index < this->get_first_face_quad_index(face))
618  // DoF is on a face
619  {
620  // do the same kind of translation as before. we need to only consider
621  // DoFs on the lines, i.e., ignoring those on the vertices
622  const unsigned int index =
623  face_index - this->get_first_face_line_index(face);
624 
625  const unsigned int face_line = index / this->n_dofs_per_line();
626  const unsigned int dof_index_on_line = index % this->n_dofs_per_line();
627 
628  return (
629  this->get_first_line_index() +
630  this->reference_cell().face_to_cell_lines(face,
631  face_line,
632  (face_orientation ? 1 : 0) +
633  (face_rotation ? 2 : 0) +
634  (face_flip ? 4 : 0)) *
635  this->n_dofs_per_line() +
636  dof_index_on_line);
637  }
638  else
639  // DoF is on a quad
640  {
641  Assert(dim >= 3, ExcInternalError());
642 
643  // ignore vertex and line dofs
644  const unsigned int index =
645  face_index - this->get_first_face_quad_index(face);
646 
647  return (this->get_first_quad_index(face) + index);
648  }
649 }
650 
651 
652 
653 template <int dim, int spacedim>
654 unsigned int
656  const unsigned int index,
657  const unsigned int face,
658  const bool face_orientation,
659  const bool face_flip,
660  const bool face_rotation) const
661 {
662  // general template for 1D and 2D: not
663  // implemented. in fact, the function
664  // shouldn't even be called unless we are
665  // in 3d, so throw an internal error
666  Assert(dim == 3, ExcInternalError());
667  if (dim < 3)
668  return index;
669 
670  // adjust dofs on 3d faces if the face is
671  // flipped. note that we query a table that
672  // derived elements need to have set up
673  // front. the exception are discontinuous
674  // elements for which there should be no
675  // face dofs anyway (i.e. dofs_per_quad==0
676  // in 3d), so we don't need the table, but
677  // the function should also not have been
678  // called
679  AssertIndexRange(index, this->n_dofs_per_quad(face));
681  [this->n_unique_quads() == 1 ? 0 : face]
682  .n_elements() == (this->reference_cell().face_reference_cell(
684  8 :
685  6) *
686  this->n_dofs_per_quad(face),
687  ExcInternalError());
688  return index +
690  [this->n_unique_quads() == 1 ? 0 : face](index,
691  (face_orientation ? 4 : 0) +
692  (face_flip ? 2 : 0) +
693  (face_rotation ? 1 : 0));
694 }
695 
696 
697 
698 template <int dim, int spacedim>
699 unsigned int
701  const unsigned int index,
702  const bool line_orientation) const
703 {
704  // general template for 1D and 2D: do
705  // nothing. Do not throw an Assertion,
706  // however, in order to allow to call this
707  // function in 2D as well
708  if (dim < 3)
709  return index;
710 
711  AssertIndexRange(index, this->n_dofs_per_line());
713  this->n_dofs_per_line(),
714  ExcInternalError());
715  if (line_orientation)
716  return index;
717  else
719 }
720 
721 
722 
723 template <int dim, int spacedim>
724 bool
726 {
727  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
728  ref_case < RefinementCase<dim>::isotropic_refinement + 1;
729  ++ref_case)
730  for (unsigned int c = 0;
731  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
732  ++c)
733  {
734  // make sure also the lazily initialized matrices are created
736  Assert((prolongation[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
737  (prolongation[ref_case - 1][c].m() == 0),
738  ExcInternalError());
739  Assert((prolongation[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
740  (prolongation[ref_case - 1][c].n() == 0),
741  ExcInternalError());
742  if ((prolongation[ref_case - 1][c].m() == 0) ||
743  (prolongation[ref_case - 1][c].n() == 0))
744  return false;
745  }
746  return true;
747 }
748 
749 
750 
751 template <int dim, int spacedim>
752 bool
754 {
755  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
756  ref_case < RefinementCase<dim>::isotropic_refinement + 1;
757  ++ref_case)
758  for (unsigned int c = 0;
759  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
760  ++c)
761  {
762  // make sure also the lazily initialized matrices are created
764  Assert((restriction[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
765  (restriction[ref_case - 1][c].m() == 0),
766  ExcInternalError());
767  Assert((restriction[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
768  (restriction[ref_case - 1][c].n() == 0),
769  ExcInternalError());
770  if ((restriction[ref_case - 1][c].m() == 0) ||
771  (restriction[ref_case - 1][c].n() == 0))
772  return false;
773  }
774  return true;
775 }
776 
777 
778 
779 template <int dim, int spacedim>
780 bool
782 {
783  const RefinementCase<dim> ref_case =
785 
786  for (unsigned int c = 0;
787  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
788  ++c)
789  {
790  // make sure also the lazily initialized matrices are created
792  Assert((prolongation[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
793  (prolongation[ref_case - 1][c].m() == 0),
794  ExcInternalError());
795  Assert((prolongation[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
796  (prolongation[ref_case - 1][c].n() == 0),
797  ExcInternalError());
798  if ((prolongation[ref_case - 1][c].m() == 0) ||
799  (prolongation[ref_case - 1][c].n() == 0))
800  return false;
801  }
802  return true;
803 }
804 
805 
806 
807 template <int dim, int spacedim>
808 bool
810 {
811  const RefinementCase<dim> ref_case =
813 
814  for (unsigned int c = 0;
815  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
816  ++c)
817  {
818  // make sure also the lazily initialized matrices are created
820  Assert((restriction[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
821  (restriction[ref_case - 1][c].m() == 0),
822  ExcInternalError());
823  Assert((restriction[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
824  (restriction[ref_case - 1][c].n() == 0),
825  ExcInternalError());
826  if ((restriction[ref_case - 1][c].m() == 0) ||
827  (restriction[ref_case - 1][c].n() == 0))
828  return false;
829  }
830  return true;
831 }
832 
833 
834 
835 template <int dim, int spacedim>
836 bool
838  const internal::SubfaceCase<dim> &subface_case) const
839 {
840  if (subface_case == internal::SubfaceCase<dim>::case_isotropic)
841  {
842  unsigned int n_dofs_on_faces = 0;
843 
844  for (const auto face_no : this->reference_cell().face_indices())
845  n_dofs_on_faces += this->n_dofs_per_face(face_no);
846 
847  return (n_dofs_on_faces == 0) || (interface_constraints.m() != 0);
848  }
849  else
850  return false;
851 }
852 
853 
854 
855 template <int dim, int spacedim>
856 bool
858 {
859  return false;
860 }
861 
862 
863 
864 template <int dim, int spacedim>
865 const FullMatrix<double> &
867  const internal::SubfaceCase<dim> &subface_case) const
868 {
869  // TODO: the implementation makes the assumption that all faces have the
870  // same number of dofs
871  AssertDimension(this->n_unique_faces(), 1);
872  const unsigned int face_no = 0;
873  (void)face_no;
874 
875  (void)subface_case;
877  ExcMessage("Constraints for this element are only implemented "
878  "for the case that faces are refined isotropically "
879  "(which is always the case in 2d, and in 3d requires "
880  "that the neighboring cell of a coarse cell presents "
881  "exactly four children on the common face)."));
882  Assert((this->n_dofs_per_face(face_no) == 0) ||
883  (interface_constraints.m() != 0),
884  ExcMessage("The finite element for which you try to obtain "
885  "hanging node constraints does not appear to "
886  "implement them."));
887 
888  if (dim == 1)
892 
893  return interface_constraints;
894 }
895 
896 
897 
898 template <int dim, int spacedim>
901 {
902  // TODO: the implementation makes the assumption that all faces have the
903  // same number of dofs
904  AssertDimension(this->n_unique_faces(), 1);
905  const unsigned int face_no = 0;
906 
907  switch (dim)
908  {
909  case 1:
910  return {0U, 0U};
911  case 2:
912  return {this->n_dofs_per_vertex() + 2 * this->n_dofs_per_line(),
913  this->n_dofs_per_face(face_no)};
914  case 3:
915  return {5 * this->n_dofs_per_vertex() + 12 * this->n_dofs_per_line() +
916  4 * this->n_dofs_per_quad(face_no),
917  this->n_dofs_per_face(face_no)};
918  default:
919  Assert(false, ExcNotImplemented());
920  }
922 }
923 
924 
925 
926 template <int dim, int spacedim>
927 void
930  FullMatrix<double> &) const
931 {
932  // by default, no interpolation
933  // implemented. so throw exception,
934  // as documentation says
935  AssertThrow(
936  false,
938 }
939 
940 
941 
942 template <int dim, int spacedim>
943 void
947  const unsigned int) const
948 {
949  // by default, no interpolation
950  // implemented. so throw exception,
951  // as documentation says
952  AssertThrow(
953  false,
955 }
956 
957 
958 
959 template <int dim, int spacedim>
960 void
963  const unsigned int,
965  const unsigned int) const
966 {
967  // by default, no interpolation
968  // implemented. so throw exception,
969  // as documentation says
970  AssertThrow(
971  false,
973 }
974 
975 
976 
977 template <int dim, int spacedim>
978 std::vector<std::pair<unsigned int, unsigned int>>
980  const FiniteElement<dim, spacedim> &) const
981 {
982  Assert(false, ExcNotImplemented());
983  return std::vector<std::pair<unsigned int, unsigned int>>();
984 }
985 
986 
987 
988 template <int dim, int spacedim>
989 std::vector<std::pair<unsigned int, unsigned int>>
991  const FiniteElement<dim, spacedim> &) const
992 {
993  Assert(false, ExcNotImplemented());
994  return std::vector<std::pair<unsigned int, unsigned int>>();
995 }
996 
997 
998 
999 template <int dim, int spacedim>
1000 std::vector<std::pair<unsigned int, unsigned int>>
1003  const unsigned int) const
1004 {
1005  Assert(false, ExcNotImplemented());
1006  return std::vector<std::pair<unsigned int, unsigned int>>();
1007 }
1008 
1009 
1010 
1011 template <int dim, int spacedim>
1015  const unsigned int) const
1016 {
1017  Assert(false, ExcNotImplemented());
1019 }
1020 
1021 
1022 
1023 template <int dim, int spacedim>
1024 bool
1026  const FiniteElement<dim, spacedim> &f) const
1027 {
1028  // Compare fields in roughly increasing order of how expensive the
1029  // comparison is
1030  return ((typeid(*this) == typeid(f)) && (this->get_name() == f.get_name()) &&
1031  (static_cast<const FiniteElementData<dim> &>(*this) ==
1032  static_cast<const FiniteElementData<dim> &>(f)) &&
1034 }
1035 
1036 
1037 
1038 template <int dim, int spacedim>
1039 bool
1041  const FiniteElement<dim, spacedim> &f) const
1042 {
1043  return !(*this == f);
1044 }
1045 
1046 
1047 
1048 template <int dim, int spacedim>
1049 const std::vector<Point<dim>> &
1051 {
1052  // a finite element may define
1053  // support points, but only if
1054  // there are as many as there are
1055  // degrees of freedom
1056  Assert((unit_support_points.size() == 0) ||
1057  (unit_support_points.size() == this->n_dofs_per_cell()),
1058  ExcInternalError());
1059  return unit_support_points;
1060 }
1061 
1062 
1063 
1064 template <int dim, int spacedim>
1065 bool
1067 {
1068  return (unit_support_points.size() != 0);
1069 }
1070 
1071 
1072 
1073 template <int dim, int spacedim>
1074 const std::vector<Point<dim>> &
1076 {
1077  // If the finite element implements generalized support points, return
1078  // those. Otherwise fall back to unit support points.
1079  return ((generalized_support_points.size() == 0) ?
1082 }
1083 
1084 
1085 
1086 template <int dim, int spacedim>
1087 bool
1089 {
1090  return (get_generalized_support_points().size() != 0);
1091 }
1092 
1093 
1094 
1095 template <int dim, int spacedim>
1096 Point<dim>
1097 FiniteElement<dim, spacedim>::unit_support_point(const unsigned int index) const
1098 {
1099  AssertIndexRange(index, this->n_dofs_per_cell());
1100  Assert(unit_support_points.size() == this->n_dofs_per_cell(),
1102  return unit_support_points[index];
1103 }
1104 
1105 
1106 
1107 template <int dim, int spacedim>
1108 const std::vector<Point<dim - 1>> &
1110  const unsigned int face_no) const
1111 {
1112  // a finite element may define
1113  // support points, but only if
1114  // there are as many as there are
1115  // degrees of freedom on a face
1116  Assert((unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1117  .size() == 0) ||
1118  (unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1119  .size() == this->n_dofs_per_face(face_no)),
1120  ExcInternalError());
1121  return unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no];
1122 }
1123 
1124 
1125 
1126 template <int dim, int spacedim>
1127 bool
1129  const unsigned int face_no) const
1130 {
1131  return (unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1132  .size() != 0);
1133 }
1134 
1135 
1136 
1137 template <int dim, int spacedim>
1138 Point<dim - 1>
1140  const unsigned int index,
1141  const unsigned int face_no) const
1142 {
1143  AssertIndexRange(index, this->n_dofs_per_face(face_no));
1144  Assert(unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1145  .size() == this->n_dofs_per_face(face_no),
1147  return unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1148  [index];
1149 }
1150 
1151 
1152 
1153 template <int dim, int spacedim>
1154 bool
1156  const unsigned int) const
1157 {
1158  return true;
1159 }
1160 
1161 
1162 
1163 template <int dim, int spacedim>
1166 {
1167  // Translate the ComponentMask into first_selected and n_components after
1168  // some error checking:
1169  const unsigned int n_total_components = this->n_components();
1170  Assert((n_total_components == mask.size()) || (mask.size() == 0),
1171  ExcMessage("The given ComponentMask has the wrong size."));
1172 
1173  const unsigned int n_selected =
1174  mask.n_selected_components(n_total_components);
1175  Assert(n_selected > 0,
1176  ExcMessage("You need at least one selected component."));
1177 
1178  const unsigned int first_selected =
1179  mask.first_selected_component(n_total_components);
1180 
1181 #ifdef DEBUG
1182  // check that it is contiguous:
1183  for (unsigned int c = 0; c < n_total_components; ++c)
1184  Assert((c < first_selected && (!mask[c])) ||
1185  (c >= first_selected && c < first_selected + n_selected &&
1186  mask[c]) ||
1187  (c >= first_selected + n_selected && !mask[c]),
1188  ExcMessage("Error: the given ComponentMask is not contiguous!"));
1189 #endif
1190 
1191  return get_sub_fe(first_selected, n_selected);
1192 }
1193 
1194 
1195 
1196 template <int dim, int spacedim>
1199  const unsigned int first_component,
1200  const unsigned int n_selected_components) const
1201 {
1202  // No complicated logic is needed here, because it is overridden in
1203  // FESystem<dim,spacedim>. Just make sure that what the user chose is valid:
1204  Assert(first_component == 0 && n_selected_components == this->n_components(),
1205  ExcMessage(
1206  "You can only select a whole FiniteElement, not a part of one."));
1207 
1208  (void)first_component;
1209  (void)n_selected_components;
1210  return *this;
1211 }
1212 
1213 
1214 
1215 template <int dim, int spacedim>
1216 std::pair<Table<2, bool>, std::vector<unsigned int>>
1218 {
1219  Assert(false, ExcNotImplemented());
1220  return std::pair<Table<2, bool>, std::vector<unsigned int>>(
1221  Table<2, bool>(this->n_components(), this->n_dofs_per_cell()),
1222  std::vector<unsigned int>(this->n_components()));
1223 }
1224 
1225 
1226 
1227 template <int dim, int spacedim>
1228 void
1231  const std::vector<Vector<double>> &,
1232  std::vector<double> &) const
1233 {
1235  ExcMessage("The element for which you are calling the current "
1236  "function does not have generalized support points (see "
1237  "the glossary for a definition of generalized support "
1238  "points). Consequently, the current function can not "
1239  "be defined and is not implemented by the element."));
1240  Assert(false, ExcNotImplemented());
1241 }
1242 
1243 
1244 
1245 template <int dim, int spacedim>
1246 std::size_t
1248 {
1249  return (
1250  sizeof(FiniteElementData<dim>) +
1262 }
1263 
1264 
1265 
1266 template <int dim, int spacedim>
1267 std::vector<unsigned int>
1269  const std::vector<ComponentMask> &nonzero_components)
1270 {
1271  std::vector<unsigned int> retval(nonzero_components.size());
1272  for (unsigned int i = 0; i < nonzero_components.size(); ++i)
1273  retval[i] = nonzero_components[i].n_selected_components();
1274  return retval;
1275 }
1276 
1277 
1278 
1279 /*------------------------------- FiniteElement ----------------------*/
1280 
1281 #ifndef DOXYGEN
1282 template <int dim, int spacedim>
1283 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1285  const UpdateFlags flags,
1286  const Mapping<dim, spacedim> & mapping,
1287  const hp::QCollection<dim - 1> &quadrature,
1289  spacedim>
1290  &output_data) const
1291 {
1292  return get_data(flags,
1293  mapping,
1295  quadrature),
1296  output_data);
1297 }
1298 
1299 
1300 
1301 template <int dim, int spacedim>
1302 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1304  const UpdateFlags flags,
1305  const Mapping<dim, spacedim> &mapping,
1306  const Quadrature<dim - 1> & quadrature,
1308  spacedim>
1309  &output_data) const
1310 {
1311  return get_data(flags,
1312  mapping,
1314  quadrature),
1315  output_data);
1316 }
1317 
1318 
1319 
1320 template <int dim, int spacedim>
1321 inline void
1323  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1324  const unsigned int face_no,
1325  const hp::QCollection<dim - 1> & quadrature,
1326  const Mapping<dim, spacedim> & mapping,
1327  const typename Mapping<dim, spacedim>::InternalDataBase & mapping_internal,
1328  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
1329  spacedim>
1330  & mapping_data,
1331  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
1333  spacedim>
1334  &output_data) const
1335 {
1336  // base class version, implement overridden function in derived classes
1337  AssertDimension(quadrature.size(), 1);
1338  fill_fe_face_values(cell,
1339  face_no,
1340  quadrature[0],
1341  mapping,
1342  mapping_internal,
1343  mapping_data,
1344  fe_internal,
1345  output_data);
1346 }
1347 
1348 
1349 
1350 template <int dim, int spacedim>
1351 inline void
1353  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1354  const unsigned int face_no,
1355  const Quadrature<dim - 1> & quadrature,
1356  const Mapping<dim, spacedim> & mapping,
1357  const typename Mapping<dim, spacedim>::InternalDataBase & mapping_internal,
1358  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
1359  spacedim>
1360  & mapping_data,
1361  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
1363  spacedim>
1364  &output_data) const
1365 {
1366  Assert(false,
1367  ExcMessage("Use of a deprecated interface, please implement "
1368  "fill_fe_face_values taking a hp::QCollection argument"));
1369  (void)cell;
1370  (void)face_no;
1371  (void)quadrature;
1372  (void)mapping;
1373  (void)mapping_internal;
1374  (void)mapping_data;
1375  (void)fe_internal;
1376  (void)output_data;
1377 }
1378 #endif
1379 
1380 
1381 
1382 template <int dim, int spacedim>
1383 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1385  const UpdateFlags flags,
1386  const Mapping<dim, spacedim> &mapping,
1387  const Quadrature<dim - 1> & quadrature,
1389  spacedim>
1390  &output_data) const
1391 {
1392  return get_data(flags,
1393  mapping,
1395  this->reference_cell(), quadrature),
1396  output_data);
1397 }
1398 
1399 
1400 
1401 template <int dim, int spacedim>
1403 FiniteElement<dim, spacedim>::base_element(const unsigned int index) const
1404 {
1405  (void)index;
1406  AssertIndexRange(index, 1);
1407  // This function should not be
1408  // called for a system element
1410  return *this;
1411 }
1412 
1413 
1414 
1415 /*------------------------------- Explicit Instantiations -------------*/
1416 #include "fe.inst"
1417 
1418 
unsigned int size() const
bool represents_the_all_selected_mask() const
unsigned int size() const
bool represents_the_all_selected_mask() const
unsigned int size() const
unsigned int get_first_line_index() const
unsigned int n_dofs_per_vertex() const
unsigned int n_dofs_per_cell() const
unsigned int n_dofs_per_line() const
unsigned int get_first_quad_index(const unsigned int quad_no=0) const
unsigned int n_dofs_per_face(unsigned int face_no=0, unsigned int child=0) const
unsigned int n_blocks() const
unsigned int n_components() const
unsigned int n_unique_quads() const
ReferenceCell reference_cell() const
unsigned int get_first_face_line_index(const unsigned int face_no=0) const
unsigned int get_first_face_quad_index(const unsigned int face_no=0) const
unsigned int n_unique_faces() const
unsigned int n_dofs_per_quad(unsigned int face_no=0) const
virtual std::size_t memory_consumption() const
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
bool isotropic_prolongation_is_implemented() const
virtual std::string get_name() const =0
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2587
virtual std::unique_ptr< InternalDataBase > get_face_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const hp::QCollection< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
const std::vector< Point< dim > > & get_generalized_support_points() const
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition: fe.h:2446
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
bool prolongation_is_implemented() const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const std::vector< Point< dim > > & get_unit_support_points() const
bool has_face_support_points(const unsigned int face_no=0) const
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
bool operator!=(const FiniteElement< dim, spacedim > &) const
virtual std::unique_ptr< InternalDataBase > get_subface_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim - 1 > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
bool has_support_points() const
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2401
TableIndices< 2 > interface_constraints_size() const
bool has_generalized_support_points() const
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2475
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const
BlockIndices base_to_block_indices
Definition: fe.h:2539
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
bool isotropic_restriction_is_implemented() const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2490
virtual Point< dim > unit_support_point(const unsigned int index) const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
virtual unsigned int face_to_cell_index(const unsigned int face_dof_index, const unsigned int face, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false) const
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
virtual void get_subface_interpolation_matrix(const FiniteElement< dim, spacedim > &source, const unsigned int subface, FullMatrix< double > &matrix, const unsigned int face_no=0) const
unsigned int component_to_block_index(const unsigned int component) const
virtual void fill_fe_face_values(const typename Triangulation< dim, spacedim >::cell_iterator &cell, const unsigned int face_no, const hp::QCollection< dim - 1 > &quadrature, const Mapping< dim, spacedim > &mapping, const typename Mapping< dim, spacedim >::InternalDataBase &mapping_internal, const ::internal::FEValuesImplementation::MappingRelatedData< dim, spacedim > &mapping_data, const InternalDataBase &fe_internal, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
virtual std::unique_ptr< InternalDataBase > get_data(const UpdateFlags update_flags, const Mapping< dim, spacedim > &mapping, const Quadrature< dim > &quadrature, ::internal::FEValuesImplementation::FiniteElementRelatedData< dim, spacedim > &output_data) const =0
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2526
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2495
const FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2569
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2562
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition: fe.h:2533
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2439
bool restriction_is_implemented() const
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition: fe.h:2507
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_quad_dof_identities(const FiniteElement< dim, spacedim > &fe_other, const unsigned int face_no=0) const
const std::vector< Point< dim - 1 > > & get_unit_face_support_points(const unsigned int face_no=0) const
FullMatrix< double > interface_constraints
Definition: fe.h:2427
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
virtual std::size_t memory_consumption() const
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2452
unsigned int adjust_quad_dof_index_for_face_orientation(const unsigned int index, const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation) const
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2578
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
virtual bool hp_constraints_are_implemented() const
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2415
types::global_dof_index first_block_of_base(const unsigned int b) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
size_type n() const
size_type m() const
Definition: point.h:111
Definition: vector.h:109
unsigned int size() const
Definition: collection.h:264
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:458
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:459
Point< 2 > first
Definition: grid_out.cc:4605
static ::ExceptionBase & ExcFEHasNoSupportPoints()
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
static ::ExceptionBase & ExcInternalError()
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
static ::ExceptionBase & ExcEmbeddingVoid()
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1501
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1695
#define AssertIndexRange(index, range)
Definition: exceptions.h:1760
static ::ExceptionBase & ExcProjectionVoid()
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
Definition: exceptions.h:1611
UpdateFlags
@ update_default
No update.
void reference_cell(Triangulation< dim, spacedim > &tria, const ReferenceCell &reference_cell)
static const char U
std::enable_if_t< std::is_fundamental< T >::value, std::size_t > memory_consumption(const T &t)
constexpr const ReferenceCell Quadrilateral
VectorType::value_type * begin(VectorType &V)
VectorType::value_type * end(VectorType &V)
static const unsigned int invalid_unsigned_int
Definition: types.h:206
static unsigned int n_children(const RefinementCase< dim > &refinement_case)