Reference documentation for deal.II version Git b43ef918fe 2022-01-29 10:49:23 +0100
\(\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 - 2021 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
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 {
88  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,
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  {
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 
111  for (unsigned int f = 0; f < this->n_unique_faces(); ++f)
112  {
114  for (unsigned int j = 0; j < this->n_dofs_per_face(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 
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)
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.
134 
135  // initialize the restriction and prolongation matrices. the default
136  // constructor of FullMatrix<dim> initializes them with size zero
139  for (unsigned int ref = RefinementCase<dim>::cut_x;
140  ref < RefinementCase<dim>::isotropic_refinement + 1;
141  ++ref)
142  {
144  RefinementCase<dim>(ref)),
147  RefinementCase<dim>(ref)),
149  }
150 
151 
152  if (dim == 3)
153  {
155  this->n_unique_quads());
156 
157  for (unsigned int f = 0; f < this->n_unique_quads(); ++f)
158  {
161  this->reference_cell().face_reference_cell(f) ==
163  8 :
164  6);
166  }
167  }
168 
171 }
172 
173 
174 
175 template <int dim, int spacedim>
176 std::pair<std::unique_ptr<FiniteElement<dim, spacedim>>, unsigned int>
177 FiniteElement<dim, spacedim>::operator^(const unsigned int multiplicity) const
178 {
179  return {this->clone(), multiplicity};
180 }
181 
182 
183 
184 template <int dim, int spacedim>
185 double
187  const Point<dim> &) const
188 {
190  return 0.;
191 }
192 
193 
194 
195 template <int dim, int spacedim>
196 double
198  const Point<dim> &,
199  const unsigned int) const
200 {
202  return 0.;
203 }
204 
205 
206 
207 template <int dim, int spacedim>
210  const Point<dim> &) const
211 {
213  return Tensor<1, dim>();
214 }
215 
216 
217 
218 template <int dim, int spacedim>
221  const Point<dim> &,
222  const unsigned int) const
223 {
225  return Tensor<1, dim>();
226 }
227 
228 
229 
230 template <int dim, int spacedim>
233  const Point<dim> &) const
234 {
236  return Tensor<2, dim>();
237 }
238 
239 
240 
241 template <int dim, int spacedim>
244  const unsigned int,
245  const Point<dim> &,
246  const unsigned int) const
247 {
249  return Tensor<2, dim>();
250 }
251 
252 
253 
254 template <int dim, int spacedim>
257  const Point<dim> &) const
258 {
260  return Tensor<3, dim>();
261 }
262 
263 
264 
265 template <int dim, int spacedim>
268  const unsigned int,
269  const Point<dim> &,
270  const unsigned int) const
271 {
273  return Tensor<3, dim>();
274 }
275 
276 
277 
278 template <int dim, int spacedim>
281  const Point<dim> &) const
282 {
284  return Tensor<4, dim>();
285 }
286 
287 
288 
289 template <int dim, int spacedim>
292  const unsigned int,
293  const Point<dim> &,
294  const unsigned int) const
295 {
297  return Tensor<4, dim>();
298 }
299 
300 
301 template <int dim, int spacedim>
302 void
304  const bool isotropic_restriction_only,
305  const bool isotropic_prolongation_only)
306 {
307  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
308  ref_case <= RefinementCase<dim>::isotropic_refinement;
309  ++ref_case)
310  {
311  const unsigned int nc =
313 
314  for (unsigned int i = 0; i < nc; ++i)
315  {
316  if (this->restriction[ref_case - 1][i].m() !=
317  this->n_dofs_per_cell() &&
318  (!isotropic_restriction_only ||
320  this->restriction[ref_case - 1][i].reinit(this->n_dofs_per_cell(),
321  this->n_dofs_per_cell());
322  if (this->prolongation[ref_case - 1][i].m() !=
323  this->n_dofs_per_cell() &&
324  (!isotropic_prolongation_only ||
326  this->prolongation[ref_case - 1][i].reinit(this->n_dofs_per_cell(),
327  this->n_dofs_per_cell());
328  }
329  }
330 }
331 
332 
333 template <int dim, int spacedim>
334 const FullMatrix<double> &
336  const unsigned int child,
337  const RefinementCase<dim> &refinement_case) const
338 {
339  AssertIndexRange(refinement_case,
341  Assert(refinement_case != RefinementCase<dim>::no_refinement,
342  ExcMessage(
343  "Restriction matrices are only available for refined cells!"));
345  child, GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)));
346  // we use refinement_case-1 here. the -1 takes care of the origin of the
347  // vector, as for RefinementCase<dim>::no_refinement (=0) there is no data
348  // available and so the vector indices are shifted
349  Assert(restriction[refinement_case - 1][child].n() == this->n_dofs_per_cell(),
351  return restriction[refinement_case - 1][child];
352 }
353 
354 
355 
356 template <int dim, int spacedim>
357 const FullMatrix<double> &
359  const unsigned int child,
360  const RefinementCase<dim> &refinement_case) const
361 {
362  AssertIndexRange(refinement_case,
364  Assert(refinement_case != RefinementCase<dim>::no_refinement,
365  ExcMessage(
366  "Prolongation matrices are only available for refined cells!"));
368  child, GeometryInfo<dim>::n_children(RefinementCase<dim>(refinement_case)));
369  // we use refinement_case-1 here. the -1 takes care
370  // of the origin of the vector, as for
371  // RefinementCase::no_refinement (=0) there is no
372  // data available and so the vector indices
373  // are shifted
374  Assert(prolongation[refinement_case - 1][child].n() ==
375  this->n_dofs_per_cell(),
376  ExcEmbeddingVoid());
377  return prolongation[refinement_case - 1][child];
378 }
379 
380 
381 // TODO:[GK] This is probably not the most efficient way of doing this.
382 template <int dim, int spacedim>
383 unsigned int
385  const unsigned int index) const
386 {
387  AssertIndexRange(index, this->n_components());
388 
389  return first_block_of_base(component_to_base_table[index].first.first) +
390  component_to_base_table[index].second;
391 }
392 
393 
394 template <int dim, int spacedim>
397  const FEValuesExtractors::Scalar &scalar) const
398 {
399  AssertIndexRange(scalar.component, this->n_components());
400 
401  // TODO: it would be nice to verify that it is indeed possible
402  // to select this scalar component, i.e., that it is not part
403  // of a non-primitive element. unfortunately, there is no simple
404  // way to write such a condition...
405 
406  std::vector<bool> mask(this->n_components(), false);
407  mask[scalar.component] = true;
408  return mask;
409 }
410 
411 
412 template <int dim, int spacedim>
415  const FEValuesExtractors::Vector &vector) const
416 {
417  AssertIndexRange(vector.first_vector_component + dim - 1,
418  this->n_components());
419 
420  // TODO: it would be nice to verify that it is indeed possible
421  // to select these vector components, i.e., that they don't span
422  // beyond the beginning or end of anon-primitive element.
423  // unfortunately, there is no simple way to write such a condition...
424 
425  std::vector<bool> mask(this->n_components(), false);
426  for (unsigned int c = vector.first_vector_component;
427  c < vector.first_vector_component + dim;
428  ++c)
429  mask[c] = true;
430  return mask;
431 }
432 
433 
434 template <int dim, int spacedim>
437  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
438 {
441  this->n_components());
442 
443  // TODO: it would be nice to verify that it is indeed possible
444  // to select these vector components, i.e., that they don't span
445  // beyond the beginning or end of anon-primitive element.
446  // unfortunately, there is no simple way to write such a condition...
447 
448  std::vector<bool> mask(this->n_components(), false);
449  for (unsigned int c = sym_tensor.first_tensor_component;
450  c < sym_tensor.first_tensor_component +
452  ++c)
453  mask[c] = true;
454  return mask;
455 }
456 
457 
458 
459 template <int dim, int spacedim>
462 {
463  // if we get a block mask that represents all blocks, then
464  // do the same for the returned component mask
465  if (block_mask.represents_the_all_selected_mask())
466  return {};
467 
468  AssertDimension(block_mask.size(), this->n_blocks());
469 
470  std::vector<bool> component_mask(this->n_components(), false);
471  for (unsigned int c = 0; c < this->n_components(); ++c)
472  if (block_mask[component_to_block_index(c)] == true)
473  component_mask[c] = true;
474 
475  return component_mask;
476 }
477 
478 
479 
480 template <int dim, int spacedim>
481 BlockMask
483  const FEValuesExtractors::Scalar &scalar) const
484 {
485  // simply create the corresponding component mask (a simpler
486  // process) and then convert it to a block mask
487  return block_mask(component_mask(scalar));
488 }
489 
490 
491 template <int dim, int spacedim>
492 BlockMask
494  const FEValuesExtractors::Vector &vector) const
495 {
496  // simply create the corresponding component mask (a simpler
497  // process) and then convert it to a block mask
498  return block_mask(component_mask(vector));
499 }
500 
501 
502 template <int dim, int spacedim>
503 BlockMask
505  const FEValuesExtractors::SymmetricTensor<2> &sym_tensor) const
506 {
507  // simply create the corresponding component mask (a simpler
508  // process) and then convert it to a block mask
509  return block_mask(component_mask(sym_tensor));
510 }
511 
512 
513 
514 template <int dim, int spacedim>
515 BlockMask
517  const ComponentMask &component_mask) const
518 {
519  // if we get a component mask that represents all component, then
520  // do the same for the returned block mask
521  if (component_mask.represents_the_all_selected_mask())
522  return {};
523 
524  AssertDimension(component_mask.size(), this->n_components());
525 
526  // walk over all of the components
527  // of this finite element and see
528  // if we need to set the
529  // corresponding block. inside the
530  // block, walk over all the
531  // components that correspond to
532  // this block and make sure the
533  // component mask is set for all of
534  // them
535  std::vector<bool> block_mask(this->n_blocks(), false);
536  for (unsigned int c = 0; c < this->n_components();)
537  {
538  const unsigned int block = component_to_block_index(c);
539  if (component_mask[c] == true)
540  block_mask[block] = true;
541 
542  // now check all of the other
543  // components that correspond
544  // to this block
545  ++c;
546  while ((c < this->n_components()) &&
547  (component_to_block_index(c) == block))
548  {
549  Assert(component_mask[c] == block_mask[block],
550  ExcMessage(
551  "The component mask argument given to this function "
552  "is not a mask where the individual components belonging "
553  "to one block of the finite element are either all "
554  "selected or not selected. You can't call this function "
555  "with a component mask that splits blocks."));
556  ++c;
557  }
558  }
559 
560 
561  return block_mask;
562 }
563 
564 
565 
566 template <int dim, int spacedim>
567 unsigned int
569  const unsigned int face,
570  const bool face_orientation,
571  const bool face_flip,
572  const bool face_rotation) const
573 {
574  AssertIndexRange(face_index, this->n_dofs_per_face(face));
575  AssertIndexRange(face, this->reference_cell().n_faces());
576 
577  // TODO: we could presumably solve the 3d case below using the
578  // adjust_quad_dof_index_for_face_orientation_table field. for the
579  // 2d case, we can't use adjust_line_dof_index_for_line_orientation_table
580  // since that array is empty (presumably because we thought that
581  // there are no flipped edges in 2d, but these can happen in
582  // DoFTools::make_periodicity_constraints, for example). so we
583  // would need to either fill this field, or rely on derived classes
584  // implementing this function, as we currently do
585 
586  // see the function's documentation for an explanation of this
587  // assertion -- in essence, derived classes have to implement
588  // an overloaded version of this function if we are to use any
589  // other than standard orientation
590  if ((face_orientation != true) || (face_flip != false) ||
591  (face_rotation != false))
592  Assert((this->n_dofs_per_line() <= 1) && (this->n_dofs_per_quad(face) <= 1),
593  ExcMessage(
594  "The function in this base class can not handle this case. "
595  "Rather, the derived class you are using must provide "
596  "an overloaded version but apparently hasn't done so. See "
597  "the documentation of this function for more information."));
598 
599  // we need to distinguish between DoFs on vertices, lines and in 3d quads.
600  // do so in a sequence of if-else statements
601  if (face_index < this->get_first_face_line_index(face))
602  // DoF is on a vertex
603  {
604  // get the number of the vertex on the face that corresponds to this DoF,
605  // along with the number of the DoF on this vertex
606  const unsigned int face_vertex = face_index / this->n_dofs_per_vertex();
607  const unsigned int dof_index_on_vertex =
608  face_index % this->n_dofs_per_vertex();
609 
610  // then get the number of this vertex on the cell and translate
611  // this to a DoF number on the cell
612  return (this->reference_cell().face_to_cell_vertices(
613  face,
614  face_vertex,
615  (face_orientation ? 1 : 0) + (face_rotation ? 2 : 0) +
616  (face_flip ? 4 : 0)) *
617  this->n_dofs_per_vertex() +
618  dof_index_on_vertex);
619  }
620  else if (face_index < this->get_first_face_quad_index(face))
621  // DoF is on a face
622  {
623  // do the same kind of translation as before. we need to only consider
624  // DoFs on the lines, i.e., ignoring those on the vertices
625  const unsigned int index =
626  face_index - this->get_first_face_line_index(face);
627 
628  const unsigned int face_line = index / this->n_dofs_per_line();
629  const unsigned int dof_index_on_line = index % this->n_dofs_per_line();
630 
631  return (
632  this->get_first_line_index() +
633  this->reference_cell().face_to_cell_lines(face,
634  face_line,
635  (face_orientation ? 1 : 0) +
636  (face_rotation ? 2 : 0) +
637  (face_flip ? 4 : 0)) *
638  this->n_dofs_per_line() +
639  dof_index_on_line);
640  }
641  else
642  // DoF is on a quad
643  {
644  Assert(dim >= 3, ExcInternalError());
645 
646  // ignore vertex and line dofs
647  const unsigned int index =
648  face_index - this->get_first_face_quad_index(face);
649 
650  return (this->get_first_quad_index(face) + index);
651  }
652 }
653 
654 
655 
656 template <int dim, int spacedim>
657 unsigned int
659  const unsigned int index,
660  const unsigned int face,
661  const bool face_orientation,
662  const bool face_flip,
663  const bool face_rotation) const
664 {
665  // general template for 1D and 2D: not
666  // implemented. in fact, the function
667  // shouldn't even be called unless we are
668  // in 3d, so throw an internal error
669  Assert(dim == 3, ExcInternalError());
670  if (dim < 3)
671  return index;
672 
673  // adjust dofs on 3d faces if the face is
674  // flipped. note that we query a table that
675  // derived elements need to have set up
676  // front. the exception are discontinuous
677  // elements for which there should be no
678  // face dofs anyway (i.e. dofs_per_quad==0
679  // in 3d), so we don't need the table, but
680  // the function should also not have been
681  // called
682  AssertIndexRange(index, this->n_dofs_per_quad(face));
684  [this->n_unique_quads() == 1 ? 0 : face]
685  .n_elements() == (this->reference_cell().face_reference_cell(
687  8 :
688  6) *
689  this->n_dofs_per_quad(face),
690  ExcInternalError());
691  return index +
693  [this->n_unique_quads() == 1 ? 0 : face](index,
694  (face_orientation ? 4 : 0) +
695  (face_flip ? 2 : 0) +
696  (face_rotation ? 1 : 0));
697 }
698 
699 
700 
701 template <int dim, int spacedim>
702 unsigned int
704  const unsigned int index,
705  const bool line_orientation) const
706 {
707  // general template for 1D and 2D: do
708  // nothing. Do not throw an Assertion,
709  // however, in order to allow to call this
710  // function in 2D as well
711  if (dim < 3)
712  return index;
713 
714  AssertIndexRange(index, this->n_dofs_per_line());
716  this->n_dofs_per_line(),
717  ExcInternalError());
718  if (line_orientation)
719  return index;
720  else
722 }
723 
724 
725 
726 template <int dim, int spacedim>
727 bool
729 {
730  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
731  ref_case < RefinementCase<dim>::isotropic_refinement + 1;
732  ++ref_case)
733  for (unsigned int c = 0;
734  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
735  ++c)
736  {
737  // make sure also the lazily initialized matrices are created
739  Assert((prolongation[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
740  (prolongation[ref_case - 1][c].m() == 0),
741  ExcInternalError());
742  Assert((prolongation[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
743  (prolongation[ref_case - 1][c].n() == 0),
744  ExcInternalError());
745  if ((prolongation[ref_case - 1][c].m() == 0) ||
746  (prolongation[ref_case - 1][c].n() == 0))
747  return false;
748  }
749  return true;
750 }
751 
752 
753 
754 template <int dim, int spacedim>
755 bool
757 {
758  for (unsigned int ref_case = RefinementCase<dim>::cut_x;
759  ref_case < RefinementCase<dim>::isotropic_refinement + 1;
760  ++ref_case)
761  for (unsigned int c = 0;
762  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
763  ++c)
764  {
765  // make sure also the lazily initialized matrices are created
767  Assert((restriction[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
768  (restriction[ref_case - 1][c].m() == 0),
769  ExcInternalError());
770  Assert((restriction[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
771  (restriction[ref_case - 1][c].n() == 0),
772  ExcInternalError());
773  if ((restriction[ref_case - 1][c].m() == 0) ||
774  (restriction[ref_case - 1][c].n() == 0))
775  return false;
776  }
777  return true;
778 }
779 
780 
781 
782 template <int dim, int spacedim>
783 bool
785 {
786  const RefinementCase<dim> ref_case =
788 
789  for (unsigned int c = 0;
790  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
791  ++c)
792  {
793  // make sure also the lazily initialized matrices are created
795  Assert((prolongation[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
796  (prolongation[ref_case - 1][c].m() == 0),
797  ExcInternalError());
798  Assert((prolongation[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
799  (prolongation[ref_case - 1][c].n() == 0),
800  ExcInternalError());
801  if ((prolongation[ref_case - 1][c].m() == 0) ||
802  (prolongation[ref_case - 1][c].n() == 0))
803  return false;
804  }
805  return true;
806 }
807 
808 
809 
810 template <int dim, int spacedim>
811 bool
813 {
814  const RefinementCase<dim> ref_case =
816 
817  for (unsigned int c = 0;
818  c < GeometryInfo<dim>::n_children(RefinementCase<dim>(ref_case));
819  ++c)
820  {
821  // make sure also the lazily initialized matrices are created
823  Assert((restriction[ref_case - 1][c].m() == this->n_dofs_per_cell()) ||
824  (restriction[ref_case - 1][c].m() == 0),
825  ExcInternalError());
826  Assert((restriction[ref_case - 1][c].n() == this->n_dofs_per_cell()) ||
827  (restriction[ref_case - 1][c].n() == 0),
828  ExcInternalError());
829  if ((restriction[ref_case - 1][c].m() == 0) ||
830  (restriction[ref_case - 1][c].n() == 0))
831  return false;
832  }
833  return true;
834 }
835 
836 
837 
838 template <int dim, int spacedim>
839 bool
841  const internal::SubfaceCase<dim> &subface_case) const
842 {
843  // TODO: the implementation makes the assumption that all faces have the
844  // same number of dofs
845  AssertDimension(this->n_unique_faces(), 1);
846  const unsigned int face_no = 0;
847 
848  if (subface_case == internal::SubfaceCase<dim>::case_isotropic)
849  return (this->n_dofs_per_face(face_no) == 0) ||
850  (interface_constraints.m() != 0);
851  else
852  return false;
853 }
854 
855 
856 
857 template <int dim, int spacedim>
858 bool
860 {
861  return false;
862 }
863 
864 
865 
866 template <int dim, int spacedim>
867 const FullMatrix<double> &
869  const internal::SubfaceCase<dim> &subface_case) const
870 {
871  // TODO: the implementation makes the assumption that all faces have the
872  // same number of dofs
873  AssertDimension(this->n_unique_faces(), 1);
874  const unsigned int face_no = 0;
875  (void)face_no;
876 
877  (void)subface_case;
879  ExcMessage("Constraints for this element are only implemented "
880  "for the case that faces are refined isotropically "
881  "(which is always the case in 2d, and in 3d requires "
882  "that the neighboring cell of a coarse cell presents "
883  "exactly four children on the common face)."));
884  Assert((this->n_dofs_per_face(face_no) == 0) ||
885  (interface_constraints.m() != 0),
886  ExcMessage("The finite element for which you try to obtain "
887  "hanging node constraints does not appear to "
888  "implement them."));
889 
890  if (dim == 1)
894 
895  return interface_constraints;
896 }
897 
898 
899 
900 template <int dim, int spacedim>
903 {
904  // TODO: the implementation makes the assumption that all faces have the
905  // same number of dofs
906  AssertDimension(this->n_unique_faces(), 1);
907  const unsigned int face_no = 0;
908 
909  switch (dim)
910  {
911  case 1:
912  return {0U, 0U};
913  case 2:
914  return {this->n_dofs_per_vertex() + 2 * this->n_dofs_per_line(),
915  this->n_dofs_per_face(face_no)};
916  case 3:
917  return {5 * this->n_dofs_per_vertex() + 12 * this->n_dofs_per_line() +
918  4 * this->n_dofs_per_quad(face_no),
919  this->n_dofs_per_face(face_no)};
920  default:
921  Assert(false, ExcNotImplemented());
922  }
924 }
925 
926 
927 
928 template <int dim, int spacedim>
929 void
932  FullMatrix<double> &) const
933 {
934  // by default, no interpolation
935  // implemented. so throw exception,
936  // as documentation says
937  AssertThrow(
938  false,
940 }
941 
942 
943 
944 template <int dim, int spacedim>
945 void
949  const unsigned int) const
950 {
951  // by default, no interpolation
952  // implemented. so throw exception,
953  // as documentation says
954  AssertThrow(
955  false,
957 }
958 
959 
960 
961 template <int dim, int spacedim>
962 void
965  const unsigned int,
967  const unsigned int) const
968 {
969  // by default, no interpolation
970  // implemented. so throw exception,
971  // as documentation says
972  AssertThrow(
973  false,
975 }
976 
977 
978 
979 template <int dim, int spacedim>
980 std::vector<std::pair<unsigned int, unsigned int>>
982  const FiniteElement<dim, spacedim> &) const
983 {
984  Assert(false, ExcNotImplemented());
985  return std::vector<std::pair<unsigned int, unsigned int>>();
986 }
987 
988 
989 
990 template <int dim, int spacedim>
991 std::vector<std::pair<unsigned int, unsigned int>>
993  const FiniteElement<dim, spacedim> &) const
994 {
995  Assert(false, ExcNotImplemented());
996  return std::vector<std::pair<unsigned int, unsigned int>>();
997 }
998 
999 
1000 
1001 template <int dim, int spacedim>
1002 std::vector<std::pair<unsigned int, unsigned int>>
1005  const unsigned int) const
1006 {
1007  Assert(false, ExcNotImplemented());
1008  return std::vector<std::pair<unsigned int, unsigned int>>();
1009 }
1010 
1011 
1012 
1013 template <int dim, int spacedim>
1017  const unsigned int) const
1018 {
1019  Assert(false, ExcNotImplemented());
1021 }
1022 
1023 
1024 
1025 template <int dim, int spacedim>
1026 bool
1028  const FiniteElement<dim, spacedim> &f) const
1029 {
1030  // Compare fields in roughly increasing order of how expensive the
1031  // comparison is
1032  return ((typeid(*this) == typeid(f)) && (this->get_name() == f.get_name()) &&
1033  (static_cast<const FiniteElementData<dim> &>(*this) ==
1034  static_cast<const FiniteElementData<dim> &>(f)) &&
1036 }
1037 
1038 
1039 
1040 template <int dim, int spacedim>
1041 bool
1043  const FiniteElement<dim, spacedim> &f) const
1044 {
1045  return !(*this == f);
1046 }
1047 
1048 
1049 
1050 template <int dim, int spacedim>
1051 const std::vector<Point<dim>> &
1053 {
1054  // a finite element may define
1055  // support points, but only if
1056  // there are as many as there are
1057  // degrees of freedom
1058  Assert((unit_support_points.size() == 0) ||
1059  (unit_support_points.size() == this->n_dofs_per_cell()),
1060  ExcInternalError());
1061  return unit_support_points;
1062 }
1063 
1064 
1065 
1066 template <int dim, int spacedim>
1067 bool
1069 {
1070  return (unit_support_points.size() != 0);
1071 }
1072 
1073 
1074 
1075 template <int dim, int spacedim>
1076 const std::vector<Point<dim>> &
1078 {
1079  // If the finite element implements generalized support points, return
1080  // those. Otherwise fall back to unit support points.
1081  return ((generalized_support_points.size() == 0) ?
1084 }
1085 
1086 
1087 
1088 template <int dim, int spacedim>
1089 bool
1091 {
1092  return (get_generalized_support_points().size() != 0);
1093 }
1094 
1095 
1096 
1097 template <int dim, int spacedim>
1098 Point<dim>
1100 {
1101  AssertIndexRange(index, this->n_dofs_per_cell());
1102  Assert(unit_support_points.size() == this->n_dofs_per_cell(),
1104  return unit_support_points[index];
1105 }
1106 
1107 
1108 
1109 template <int dim, int spacedim>
1110 const std::vector<Point<dim - 1>> &
1112  const unsigned int face_no) const
1113 {
1114  // a finite element may define
1115  // support points, but only if
1116  // there are as many as there are
1117  // degrees of freedom on a face
1118  Assert((unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1119  .size() == 0) ||
1120  (unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1121  .size() == this->n_dofs_per_face(face_no)),
1122  ExcInternalError());
1123  return unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no];
1124 }
1125 
1126 
1127 
1128 template <int dim, int spacedim>
1129 bool
1131  const unsigned int face_no) const
1132 {
1133  return (unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1134  .size() != 0);
1135 }
1136 
1137 
1138 
1139 template <int dim, int spacedim>
1140 Point<dim - 1>
1142  const unsigned int index,
1143  const unsigned int face_no) const
1144 {
1145  AssertIndexRange(index, this->n_dofs_per_face(face_no));
1146  Assert(unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1147  .size() == this->n_dofs_per_face(face_no),
1149  return unit_face_support_points[this->n_unique_faces() == 1 ? 0 : face_no]
1150  [index];
1151 }
1152 
1153 
1154 
1155 template <int dim, int spacedim>
1156 bool
1158  const unsigned int) const
1159 {
1160  return true;
1161 }
1162 
1163 
1164 
1165 template <int dim, int spacedim>
1168 {
1169  // Translate the ComponentMask into first_selected and n_components after
1170  // some error checking:
1171  const unsigned int n_total_components = this->n_components();
1172  Assert((n_total_components == mask.size()) || (mask.size() == 0),
1173  ExcMessage("The given ComponentMask has the wrong size."));
1174 
1175  const unsigned int n_selected =
1176  mask.n_selected_components(n_total_components);
1177  Assert(n_selected > 0,
1178  ExcMessage("You need at least one selected component."));
1179 
1180  const unsigned int first_selected =
1181  mask.first_selected_component(n_total_components);
1182 
1183 #ifdef DEBUG
1184  // check that it is contiguous:
1185  for (unsigned int c = 0; c < n_total_components; ++c)
1186  Assert((c < first_selected && (!mask[c])) ||
1187  (c >= first_selected && c < first_selected + n_selected &&
1188  mask[c]) ||
1189  (c >= first_selected + n_selected && !mask[c]),
1190  ExcMessage("Error: the given ComponentMask is not contiguous!"));
1191 #endif
1192 
1193  return get_sub_fe(first_selected, n_selected);
1194 }
1195 
1196 
1197 
1198 template <int dim, int spacedim>
1201  const unsigned int first_component,
1202  const unsigned int n_selected_components) const
1203 {
1204  // No complicated logic is needed here, because it is overridden in
1205  // FESystem<dim,spacedim>. Just make sure that what the user chose is valid:
1206  Assert(first_component == 0 && n_selected_components == this->n_components(),
1207  ExcMessage(
1208  "You can only select a whole FiniteElement, not a part of one."));
1209 
1210  (void)first_component;
1211  (void)n_selected_components;
1212  return *this;
1213 }
1214 
1215 
1216 
1217 template <int dim, int spacedim>
1218 std::pair<Table<2, bool>, std::vector<unsigned int>>
1220 {
1221  Assert(false, ExcNotImplemented());
1222  return std::pair<Table<2, bool>, std::vector<unsigned int>>(
1223  Table<2, bool>(this->n_components(), this->n_dofs_per_cell()),
1224  std::vector<unsigned int>(this->n_components()));
1225 }
1226 
1227 
1228 
1229 template <int dim, int spacedim>
1230 void
1233  const std::vector<Vector<double>> &,
1234  std::vector<double> &) const
1235 {
1237  ExcMessage("The element for which you are calling the current "
1238  "function does not have generalized support points (see "
1239  "the glossary for a definition of generalized support "
1240  "points). Consequently, the current function can not "
1241  "be defined and is not implemented by the element."));
1242  Assert(false, ExcNotImplemented());
1243 }
1244 
1245 
1246 
1247 template <int dim, int spacedim>
1248 std::size_t
1250 {
1251  return (
1252  sizeof(FiniteElementData<dim>) +
1264 }
1265 
1266 
1267 
1268 template <int dim, int spacedim>
1269 std::vector<unsigned int>
1271  const std::vector<ComponentMask> &nonzero_components)
1272 {
1273  std::vector<unsigned int> retval(nonzero_components.size());
1274  for (unsigned int i = 0; i < nonzero_components.size(); ++i)
1275  retval[i] = nonzero_components[i].n_selected_components();
1276  return retval;
1277 }
1278 
1279 
1280 
1281 /*------------------------------- FiniteElement ----------------------*/
1282 
1283 #ifndef DOXYGEN
1284 template <int dim, int spacedim>
1285 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1287  const UpdateFlags flags,
1288  const Mapping<dim, spacedim> & mapping,
1289  const hp::QCollection<dim - 1> &quadrature,
1291  spacedim>
1292  &output_data) const
1293 {
1294  return get_data(flags,
1295  mapping,
1297  quadrature),
1298  output_data);
1299 }
1300 
1301 
1302 
1303 template <int dim, int spacedim>
1304 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1306  const UpdateFlags flags,
1307  const Mapping<dim, spacedim> &mapping,
1308  const Quadrature<dim - 1> & quadrature,
1310  spacedim>
1311  &output_data) const
1312 {
1313  return get_data(flags,
1314  mapping,
1316  quadrature),
1317  output_data);
1318 }
1319 
1320 
1321 
1322 template <int dim, int spacedim>
1323 inline void
1325  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1326  const unsigned int face_no,
1327  const hp::QCollection<dim - 1> & quadrature,
1328  const Mapping<dim, spacedim> & mapping,
1329  const typename Mapping<dim, spacedim>::InternalDataBase & mapping_internal,
1330  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
1331  spacedim>
1332  & mapping_data,
1333  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
1335  spacedim>
1336  &output_data) const
1337 {
1338  // base class version, implement overridden function in derived classes
1339  AssertDimension(quadrature.size(), 1);
1340  fill_fe_face_values(cell,
1341  face_no,
1342  quadrature[0],
1343  mapping,
1344  mapping_internal,
1345  mapping_data,
1346  fe_internal,
1347  output_data);
1348 }
1349 
1350 
1351 
1352 template <int dim, int spacedim>
1353 inline void
1355  const typename Triangulation<dim, spacedim>::cell_iterator &cell,
1356  const unsigned int face_no,
1357  const Quadrature<dim - 1> & quadrature,
1358  const Mapping<dim, spacedim> & mapping,
1359  const typename Mapping<dim, spacedim>::InternalDataBase & mapping_internal,
1360  const ::internal::FEValuesImplementation::MappingRelatedData<dim,
1361  spacedim>
1362  & mapping_data,
1363  const typename FiniteElement<dim, spacedim>::InternalDataBase &fe_internal,
1365  spacedim>
1366  &output_data) const
1367 {
1368  Assert(false,
1369  ExcMessage("Use of a deprecated interface, please implement "
1370  "fill_fe_face_values taking a hp::QCollection argument"));
1371  (void)cell;
1372  (void)face_no;
1373  (void)quadrature;
1374  (void)mapping;
1375  (void)mapping_internal;
1376  (void)mapping_data;
1377  (void)fe_internal;
1378  (void)output_data;
1379 }
1380 #endif
1381 
1382 
1383 
1384 template <int dim, int spacedim>
1385 std::unique_ptr<typename FiniteElement<dim, spacedim>::InternalDataBase>
1387  const UpdateFlags flags,
1388  const Mapping<dim, spacedim> &mapping,
1389  const Quadrature<dim - 1> & quadrature,
1391  spacedim>
1392  &output_data) const
1393 {
1394  return get_data(flags,
1395  mapping,
1397  this->reference_cell(), quadrature),
1398  output_data);
1399 }
1400 
1401 
1402 
1403 template <int dim, int spacedim>
1405 FiniteElement<dim, spacedim>::base_element(const unsigned int index) const
1406 {
1407  (void)index;
1408  AssertIndexRange(index, 1);
1409  // This function should not be
1410  // called for a system element
1412  return *this;
1413 }
1414 
1415 
1416 
1417 /*------------------------------- Explicit Instantiations -------------*/
1418 #include "fe.inst"
1419 
1420 
static ::ExceptionBase & ExcFEHasNoSupportPoints()
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > face_system_to_component_table
Definition: fe.h:2518
bool represents_the_all_selected_mask() const
size_type m() const
static const unsigned int invalid_unsigned_int
Definition: types.h:196
BlockMask block_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:482
virtual Tensor< 2, dim > shape_grad_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:243
std::vector< std::vector< FullMatrix< double > > > restriction
Definition: fe.h:2412
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
Definition: fe.cc:658
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1655
virtual std::size_t memory_consumption() const
Definition: fe.cc:49
virtual bool operator==(const FiniteElement< dim, spacedim > &fe) const
Definition: fe.cc:1027
std::vector< Point< dim > > generalized_support_points
Definition: fe.h:2463
virtual Tensor< 3, dim > shape_3rd_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:256
unsigned int get_first_line_index() const
FullMatrix< double > interface_constraints
Definition: fe.h:2438
unsigned int size() const
Definition: collection.h:263
bool operator!=(const FiniteElement< dim, spacedim > &) const
Definition: fe.cc:1042
FiniteElement(const FiniteElementData< dim > &fe_data, const std::vector< bool > &restriction_is_additive_flags, const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:57
const FiniteElement< dim, spacedim > & get_sub_fe(const ComponentMask &mask) const
Definition: fe.cc:1167
virtual double shape_value_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:197
virtual void get_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix) const
Definition: fe.cc:930
static ::ExceptionBase & ExcUnitShapeValuesDoNotExist()
virtual Point< dim > unit_support_point(const unsigned int index) const
Definition: fe.cc:1099
bool restriction_is_implemented() const
Definition: fe.cc:756
const std::vector< Point< dim > > & get_generalized_support_points() const
Definition: fe.cc:1077
virtual Tensor< 1, dim > shape_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:209
virtual Tensor< 3, dim > shape_3rd_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:267
bool has_face_support_points(const unsigned int face_no=0) const
Definition: fe.cc:1130
#define AssertIndexRange(index, range)
Definition: exceptions.h:1720
bool isotropic_prolongation_is_implemented() const
Definition: fe.cc:784
unsigned int n_dofs_per_vertex() const
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > component_to_base_table
Definition: fe.h:2573
bool has_generalized_support_points() const
Definition: fe.cc:1090
bool represents_the_all_selected_mask() const
unsigned int n_unique_quads() const
STL namespace.
ComponentMask component_mask(const FEValuesExtractors::Scalar &scalar) const
Definition: fe.cc:396
static const char U
unsigned int get_first_face_line_index(const unsigned int face_no=0) const
virtual std::unique_ptr< FiniteElement< dim, spacedim > > clone() const =0
#define AssertThrow(cond, exc)
Definition: exceptions.h:1571
unsigned int n_dofs_per_quad(unsigned int face_no=0) const
static std::vector< unsigned int > compute_n_nonzero_components(const std::vector< ComponentMask > &nonzero_components)
Definition: fe.cc:1270
virtual const FullMatrix< double > & get_prolongation_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: fe.cc:358
bool is_primitive() const
Definition: fe.h:3311
const std::vector< unsigned int > n_nonzero_components_table
Definition: fe.h:2598
unsigned int get_first_face_quad_index(const unsigned int face_no=0) const
unsigned int n_blocks() const
std::vector< Point< dim > > unit_support_points
Definition: fe.h:2450
bool has_support_points() const
Definition: fe.cc:1068
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
Definition: fe.cc:1386
unsigned int component_to_block_index(const unsigned int component) const
Definition: fe.cc:384
virtual Point< dim - 1 > unit_face_support_point(const unsigned int index, const unsigned int face_no=0) const
Definition: fe.cc:1141
static ::ExceptionBase & ExcMessage(std::string arg1)
const FullMatrix< double > & constraints(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:868
virtual Tensor< 1, dim > shape_grad_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:220
virtual std::size_t memory_consumption() const
Definition: fe.cc:1249
size_type n() const
unsigned int n_dofs_per_line() const
No update.
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2426
virtual void get_face_interpolation_matrix(const FiniteElement< dim, spacedim > &source, FullMatrix< double > &matrix, const unsigned int face_no=0) const
Definition: fe.cc:946
#define Assert(cond, exc)
Definition: exceptions.h:1461
UpdateFlags
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
Definition: fe.cc:335
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
Abstract base class for mapping classes.
Definition: mapping.h:310
void reinit(const unsigned int n_blocks, const size_type n_elements_per_block)
const std::vector< Point< dim - 1 > > & get_unit_face_support_points(const unsigned int face_no=0) const
Definition: fe.cc:1111
std::vector< std::vector< Point< dim - 1 > > > unit_face_support_points
Definition: fe.h:2457
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:407
bool constraints_are_implemented(const ::internal::SubfaceCase< dim > &subface_case=::internal::SubfaceCase< dim >::case_isotropic) const
Definition: fe.cc:840
static ::ExceptionBase & ExcWrongInterfaceMatrixSize(int arg1, int arg2)
VectorType::value_type * end(VectorType &V)
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
unsigned int size() const
virtual const FiniteElement< dim, spacedim > & base_element(const unsigned int index) const
Definition: fe.cc:1405
virtual std::string get_name() const =0
static ::ExceptionBase & ExcProjectionVoid()
std::vector< std::vector< Point< dim - 1 > > > generalized_face_support_points
Definition: fe.h:2469
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_vertex_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:981
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
Definition: fe.cc:568
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_unique_faces() 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
Definition: fe.cc:1232
Point< 2 > first
Definition: grid_out.cc:4614
types::global_dof_index first_block_of_base(const unsigned int b) const
Definition: fe.h:3231
static unsigned int n_children(const RefinementCase< dim > &refinement_case)
unsigned int adjust_line_dof_index_for_line_orientation(const unsigned int index, const bool line_orientation) const
Definition: fe.cc:703
const std::vector< Point< dim > > & get_unit_support_points() const
Definition: fe.cc:1052
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
Definition: fe.cc:963
virtual std::vector< std::pair< unsigned int, unsigned int > > hp_line_dof_identities(const FiniteElement< dim, spacedim > &fe_other) const
Definition: fe.cc:992
virtual FiniteElementDomination::Domination compare_for_domination(const FiniteElement< dim, spacedim > &fe_other, const unsigned int codim=0) const
Definition: fe.cc:1015
static ::ExceptionBase & ExcEmbeddingVoid()
virtual double shape_value(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:186
unsigned int n_components() const
unsigned int n_dofs_per_cell() const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:406
VectorType::value_type * begin(VectorType &V)
constexpr const ReferenceCell Quadrilateral
std::vector< std::pair< unsigned int, unsigned int > > system_to_component_table
Definition: fe.h:2506
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
virtual bool hp_constraints_are_implemented() const
Definition: fe.cc:859
std::vector< std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > > face_system_to_base_table
Definition: fe.h:2544
unsigned int first_selected_component(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
unsigned int n_selected_components(const unsigned int overall_number_of_components=numbers::invalid_unsigned_int) const
bool prolongation_is_implemented() const
Definition: fe.cc:728
virtual Tensor< 4, dim > shape_4th_derivative_component(const unsigned int i, const Point< dim > &p, const unsigned int component) const
Definition: fe.cc:291
virtual Tensor< 2, dim > shape_grad_grad(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:232
const std::vector< ComponentMask > nonzero_components
Definition: fe.h:2589
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:303
std::vector< std::pair< std::pair< unsigned int, unsigned int >, unsigned int > > system_to_base_table
Definition: fe.h:2537
static ::ExceptionBase & ExcNotImplemented()
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
Definition: fe.cc:1003
bool isotropic_restriction_is_implemented() const
Definition: fe.cc:812
BlockIndices base_to_block_indices
Definition: fe.h:2550
virtual Tensor< 4, dim > shape_4th_derivative(const unsigned int i, const Point< dim > &p) const
Definition: fe.cc:280
unsigned int size() const
unsigned int size() const
ReferenceCell face_reference_cell(const unsigned int face_no) const
std::vector< int > adjust_line_dof_index_for_line_orientation_table
Definition: fe.h:2501
TableIndices< 2 > interface_constraints_size() const
Definition: fe.cc:902
virtual std::pair< Table< 2, bool >, std::vector< unsigned int > > get_constant_modes() const
Definition: fe.cc:1219
const std::vector< bool > restriction_is_additive_flags
Definition: fe.h:2580
std::pair< std::unique_ptr< FiniteElement< dim, spacedim > >, unsigned int > operator^(const unsigned int multiplicity) const
Definition: fe.cc:177
std::enable_if< std::is_fundamental< T >::value, std::size_t >::type memory_consumption(const T &t)
ReferenceCell reference_cell() const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const
Definition: fe.cc:1157
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
std::vector< Table< 2, int > > adjust_quad_dof_index_for_face_orientation_table
Definition: fe.h:2486
static ::ExceptionBase & ExcInternalError()