Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
fe_abf.cc
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2019 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 
16 
17 #include <deal.II/base/qprojector.h>
18 #include <deal.II/base/quadrature.h>
19 #include <deal.II/base/quadrature_lib.h>
20 #include <deal.II/base/std_cxx14/memory.h>
21 #include <deal.II/base/table.h>
22 #include <deal.II/base/utilities.h>
23 
24 #include <deal.II/dofs/dof_accessor.h>
25 
26 #include <deal.II/fe/fe.h>
27 #include <deal.II/fe/fe_abf.h>
28 #include <deal.II/fe/fe_tools.h>
29 #include <deal.II/fe/fe_values.h>
30 #include <deal.II/fe/mapping.h>
31 
32 #include <deal.II/grid/tria.h>
33 #include <deal.II/grid/tria_iterator.h>
34 
35 #include <iostream>
36 #include <sstream>
37 
38 
39 // TODO: implement the adjust_quad_dof_index_for_face_orientation_table and
40 // adjust_line_dof_index_for_line_orientation_table fields, and write tests
41 // similar to bits/face_orientation_and_fe_q_*
42 
43 
44 DEAL_II_NAMESPACE_OPEN
45 
46 
47 template <int dim>
48 FE_ABF<dim>::FE_ABF(const unsigned int deg)
49  : FE_PolyTensor<PolynomialsABF<dim>, dim>(
50  deg,
51  FiniteElementData<dim>(get_dpo_vector(deg),
52  dim,
53  deg + 2,
54  FiniteElementData<dim>::Hdiv),
55  std::vector<bool>(PolynomialsABF<dim>::compute_n_pols(deg), true),
56  std::vector<ComponentMask>(PolynomialsABF<dim>::compute_n_pols(deg),
57  std::vector<bool>(dim, true)))
58  , rt_order(deg)
59 {
60  Assert(dim >= 2, ExcImpossibleInDim(dim));
61  const unsigned int n_dofs = this->dofs_per_cell;
62 
64  // First, initialize the
65  // generalized support points and
66  // quadrature weights, since they
67  // are required for interpolation.
69 
70  // Now compute the inverse node matrix, generating the correct
71  // basis functions from the raw ones. For a discussion of what
72  // exactly happens here, see FETools::compute_node_matrix.
74  this->inverse_node_matrix.reinit(n_dofs, n_dofs);
75  this->inverse_node_matrix.invert(M);
76  // From now on, the shape functions provided by FiniteElement::shape_value
77  // and similar functions will be the correct ones, not
78  // the raw shape functions from the polynomial space anymore.
79 
80  // Reinit the vectors of
81  // restriction and prolongation
82  // matrices to the right sizes.
83  // Restriction only for isotropic
84  // refinement
86  // Fill prolongation matrices with embedding operators
87  FETools::compute_embedding_matrices(*this, this->prolongation, false, 1.e-10);
88 
90 
91  // TODO[TL]: for anisotropic refinement we will probably need a table of
92  // submatrices with an array for each refine case
93  std::vector<FullMatrix<double>> face_embeddings(
94  1 << (dim - 1),
96  // TODO: Something goes wrong there. The error of the least squares fit
97  // is to large ...
98  // FETools::compute_face_embedding_matrices(*this, face_embeddings.data(), 0,
99  // 0);
100  this->interface_constraints.reinit((1 << (dim - 1)) * this->dofs_per_face,
101  this->dofs_per_face);
102  unsigned int target_row = 0;
103  for (const auto &face_embedding : face_embeddings)
104  for (unsigned int i = 0; i < face_embedding.m(); ++i)
105  {
106  for (unsigned int j = 0; j < face_embedding.n(); ++j)
107  this->interface_constraints(target_row, j) = face_embedding(i, j);
108  ++target_row;
109  }
110 }
111 
112 
113 
114 template <int dim>
115 std::string
117 {
118  // note that the
119  // FETools::get_fe_by_name
120  // function depends on the
121  // particular format of the string
122  // this function returns, so they
123  // have to be kept in synch
124 
125  std::ostringstream namebuf;
126 
127  namebuf << "FE_ABF<" << dim << ">(" << rt_order << ")";
128 
129  return namebuf.str();
130 }
131 
132 
133 
134 template <int dim>
135 std::unique_ptr<FiniteElement<dim, dim>>
137 {
138  return std_cxx14::make_unique<FE_ABF<dim>>(rt_order);
139 }
140 
141 
142 //---------------------------------------------------------------------------
143 // Auxiliary and internal functions
144 //---------------------------------------------------------------------------
145 
146 
147 
148 // Version for 2d and higher. See above for 1d version
149 template <int dim>
150 void
152 {
153  QGauss<dim> cell_quadrature(deg + 2);
154  const unsigned int n_interior_points = cell_quadrature.size();
155 
156  unsigned int n_face_points = (dim > 1) ? 1 : 0;
157  // compute (deg+1)^(dim-1)
158  for (unsigned int d = 1; d < dim; ++d)
159  n_face_points *= deg + 1;
160 
161  this->generalized_support_points.resize(
162  GeometryInfo<dim>::faces_per_cell * n_face_points + n_interior_points);
163  this->generalized_face_support_points.resize(n_face_points);
164 
165 
166  // These might be required when the faces contribution is computed
167  // Therefore they will be initialized at this point.
168  std::vector<AnisotropicPolynomials<dim> *> polynomials_abf(dim);
169 
170  // Generate x_1^{i} x_2^{r+1} ...
171  for (unsigned int dd = 0; dd < dim; ++dd)
172  {
173  std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
174  for (unsigned int d = 0; d < dim; ++d)
175  poly[d].push_back(Polynomials::Monomial<double>(deg + 1));
177 
178  polynomials_abf[dd] = new AnisotropicPolynomials<dim>(poly);
179  }
180 
181  // Number of the point being entered
182  unsigned int current = 0;
183 
184  if (dim > 1)
185  {
186  QGauss<dim - 1> face_points(deg + 1);
187  TensorProductPolynomials<dim - 1> legendre =
189 
190  boundary_weights.reinit(n_face_points, legendre.n());
191 
192  // Assert (face_points.size() == this->dofs_per_face,
193  // ExcInternalError());
194 
195  for (unsigned int k = 0; k < n_face_points; ++k)
196  {
197  this->generalized_face_support_points[k] = face_points.point(k);
198  // Compute its quadrature
199  // contribution for each
200  // moment.
201  for (unsigned int i = 0; i < legendre.n(); ++i)
202  {
203  boundary_weights(k, i) =
204  face_points.weight(k) *
205  legendre.compute_value(i, face_points.point(k));
206  }
207  }
208 
209  Quadrature<dim> faces =
211  for (; current < GeometryInfo<dim>::faces_per_cell * n_face_points;
212  ++current)
213  {
214  // Enter the support point
215  // into the vector
216  this->generalized_support_points[current] = faces.point(current);
217  }
218 
219 
220  // Now initialize edge interior weights for the ABF elements.
221  // These are completely independent from the usual edge moments. They
222  // stem from applying the Gauss theorem to the nodal values, which
223  // was necessary to cast the ABF elements into the deal.II framework
224  // for vector valued elements.
225  boundary_weights_abf.reinit(faces.size(), polynomials_abf[0]->n() * dim);
226  for (unsigned int k = 0; k < faces.size(); ++k)
227  {
228  for (unsigned int i = 0; i < polynomials_abf[0]->n() * dim; ++i)
229  {
230  boundary_weights_abf(k, i) =
231  polynomials_abf[i % dim]->compute_value(i / dim,
232  faces.point(k)) *
233  faces.weight(k);
234  }
235  }
236  }
237 
238  // Create Legendre basis for the
239  // space D_xi Q_k
240  if (deg > 0)
241  {
242  std::vector<AnisotropicPolynomials<dim> *> polynomials(dim);
243 
244  for (unsigned int dd = 0; dd < dim; ++dd)
245  {
246  std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
247  for (unsigned int d = 0; d < dim; ++d)
250 
251  polynomials[dd] = new AnisotropicPolynomials<dim>(poly);
252  }
253 
254  interior_weights.reinit(
255  TableIndices<3>(n_interior_points, polynomials[0]->n(), dim));
256 
257  for (unsigned int k = 0; k < cell_quadrature.size(); ++k)
258  {
259  for (unsigned int i = 0; i < polynomials[0]->n(); ++i)
260  for (unsigned int d = 0; d < dim; ++d)
261  interior_weights(k, i, d) =
262  cell_quadrature.weight(k) *
263  polynomials[d]->compute_value(i, cell_quadrature.point(k));
264  }
265 
266  for (unsigned int d = 0; d < dim; ++d)
267  delete polynomials[d];
268  }
269 
270 
271  // Decouple the creation of the generalized support points
272  // from computation of interior weights.
273  for (unsigned int k = 0; k < cell_quadrature.size(); ++k)
274  this->generalized_support_points[current++] = cell_quadrature.point(k);
275 
276  // Additional functionality for the ABF elements
277  // TODO: Here the canonical extension of the principle
278  // behind the ABF elements is implemented. It is unclear,
279  // if this really leads to the ABF spaces in 3D!
280  interior_weights_abf.reinit(TableIndices<3>(cell_quadrature.size(),
281  polynomials_abf[0]->n() * dim,
282  dim));
283  Tensor<1, dim> poly_grad;
284 
285  for (unsigned int k = 0; k < cell_quadrature.size(); ++k)
286  {
287  for (unsigned int i = 0; i < polynomials_abf[0]->n() * dim; ++i)
288  {
289  poly_grad =
290  polynomials_abf[i % dim]->compute_grad(i / dim,
291  cell_quadrature.point(k)) *
292  cell_quadrature.weight(k);
293  // The minus sign comes from the use of the Gauss theorem to replace
294  // the divergence.
295  for (unsigned int d = 0; d < dim; ++d)
296  interior_weights_abf(k, i, d) = -poly_grad[d];
297  }
298  }
299 
300  for (unsigned int d = 0; d < dim; ++d)
301  delete polynomials_abf[d];
302 
303  Assert(current == this->generalized_support_points.size(),
304  ExcInternalError());
305 }
306 
307 
308 
309 // This function is the same Raviart-Thomas interpolation performed by
310 // interpolate. Still, we cannot use interpolate, since it was written
311 // for smooth functions. The functions interpolated here are not
312 // smooth, maybe even not continuous. Therefore, we must double the
313 // number of quadrature points in each direction in order to integrate
314 // only smooth functions.
315 
316 // Then again, the interpolated function is chosen such that the
317 // moments coincide with the function to be interpolated.
318 
319 template <int dim>
320 void
322 {
323  if (dim == 1)
324  {
325  unsigned int iso = RefinementCase<dim>::isotropic_refinement - 1;
326  for (unsigned int i = 0; i < GeometryInfo<dim>::max_children_per_cell;
327  ++i)
328  this->restriction[iso][i].reinit(0, 0);
329  return;
330  }
331  unsigned int iso = RefinementCase<dim>::isotropic_refinement - 1;
332  QGauss<dim - 1> q_base(rt_order + 1);
333  const unsigned int n_face_points = q_base.size();
334  // First, compute interpolation on
335  // subfaces
336  for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
337  {
338  // The shape functions of the
339  // child cell are evaluated
340  // in the quadrature points
341  // of a full face.
342  Quadrature<dim> q_face = QProjector<dim>::project_to_face(q_base, face);
343  // Store shape values, since the
344  // evaluation suffers if not
345  // ordered by point
346  Table<2, double> cached_values_face(this->dofs_per_cell, q_face.size());
347  for (unsigned int k = 0; k < q_face.size(); ++k)
348  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
349  cached_values_face(i, k) = this->shape_value_component(
350  i, q_face.point(k), GeometryInfo<dim>::unit_normal_direction[face]);
351 
352  for (unsigned int sub = 0; sub < GeometryInfo<dim>::max_children_per_face;
353  ++sub)
354  {
355  // The weight functions for
356  // the coarse face are
357  // evaluated on the subface
358  // only.
359  Quadrature<dim> q_sub =
360  QProjector<dim>::project_to_subface(q_base, face, sub);
361  const unsigned int child = GeometryInfo<dim>::child_cell_on_face(
363 
364  // On a certain face, we must
365  // compute the moments of ALL
366  // fine level functions with
367  // the coarse level weight
368  // functions belonging to
369  // that face. Due to the
370  // orthogonalization process
371  // when building the shape
372  // functions, these weights
373  // are equal to the
374  // corresponding shape
375  // functions.
376  for (unsigned int k = 0; k < n_face_points; ++k)
377  for (unsigned int i_child = 0; i_child < this->dofs_per_cell;
378  ++i_child)
379  for (unsigned int i_face = 0; i_face < this->dofs_per_face;
380  ++i_face)
381  {
382  // The quadrature
383  // weights on the
384  // subcell are NOT
385  // transformed, so we
386  // have to do it here.
387  this->restriction[iso][child](face * this->dofs_per_face +
388  i_face,
389  i_child) +=
390  Utilities::fixed_power<dim - 1>(.5) * q_sub.weight(k) *
391  cached_values_face(i_child, k) *
392  this->shape_value_component(
393  face * this->dofs_per_face + i_face,
394  q_sub.point(k),
396  }
397  }
398  }
399 
400  if (rt_order == 0)
401  return;
402 
403  // Create Legendre basis for the
404  // space D_xi Q_k. Here, we cannot
405  // use the shape functions
406  std::vector<AnisotropicPolynomials<dim> *> polynomials(dim);
407  for (unsigned int dd = 0; dd < dim; ++dd)
408  {
409  std::vector<std::vector<Polynomials::Polynomial<double>>> poly(dim);
410  for (unsigned int d = 0; d < dim; ++d)
412  poly[dd] = Polynomials::Legendre::generate_complete_basis(rt_order - 1);
413 
414  polynomials[dd] = new AnisotropicPolynomials<dim>(poly);
415  }
416 
417  QGauss<dim> q_cell(rt_order + 1);
418  const unsigned int start_cell_dofs =
419  GeometryInfo<dim>::faces_per_cell * this->dofs_per_face;
420 
421  // Store shape values, since the
422  // evaluation suffers if not
423  // ordered by point
424  Table<3, double> cached_values_cell(this->dofs_per_cell, q_cell.size(), dim);
425  for (unsigned int k = 0; k < q_cell.size(); ++k)
426  for (unsigned int i = 0; i < this->dofs_per_cell; ++i)
427  for (unsigned int d = 0; d < dim; ++d)
428  cached_values_cell(i, k, d) =
429  this->shape_value_component(i, q_cell.point(k), d);
430 
431  for (unsigned int child = 0; child < GeometryInfo<dim>::max_children_per_cell;
432  ++child)
433  {
434  Quadrature<dim> q_sub = QProjector<dim>::project_to_child(q_cell, child);
435 
436  for (unsigned int k = 0; k < q_sub.size(); ++k)
437  for (unsigned int i_child = 0; i_child < this->dofs_per_cell; ++i_child)
438  for (unsigned int d = 0; d < dim; ++d)
439  for (unsigned int i_weight = 0; i_weight < polynomials[d]->n();
440  ++i_weight)
441  {
442  this->restriction[iso][child](start_cell_dofs + i_weight * dim +
443  d,
444  i_child) +=
445  q_sub.weight(k) * cached_values_cell(i_child, k, d) *
446  polynomials[d]->compute_value(i_weight, q_sub.point(k));
447  }
448  }
449 
450  for (unsigned int d = 0; d < dim; ++d)
451  delete polynomials[d];
452 }
453 
454 
455 
456 template <int dim>
457 std::vector<unsigned int>
458 FE_ABF<dim>::get_dpo_vector(const unsigned int rt_order)
459 {
460  if (dim == 1)
461  {
462  Assert(false, ExcImpossibleInDim(1));
463  return std::vector<unsigned int>();
464  }
465 
466  // the element is face-based (not
467  // to be confused with George
468  // W. Bush's Faith Based
469  // Initiative...), and we have
470  // (rt_order+1)^(dim-1) DoFs per face
471  unsigned int dofs_per_face = 1;
472  for (int d = 0; d < dim - 1; ++d)
473  dofs_per_face *= rt_order + 1;
474 
475  // and then there are interior dofs
476  const unsigned int interior_dofs = dim * (rt_order + 1) * dofs_per_face;
477 
478  std::vector<unsigned int> dpo(dim + 1);
479  dpo[dim - 1] = dofs_per_face;
480  dpo[dim] = interior_dofs;
481 
482  return dpo;
483 }
484 
485 //---------------------------------------------------------------------------
486 // Data field initialization
487 //---------------------------------------------------------------------------
488 
489 template <int dim>
490 bool
491 FE_ABF<dim>::has_support_on_face(const unsigned int shape_index,
492  const unsigned int face_index) const
493 {
494  Assert(shape_index < this->dofs_per_cell,
495  ExcIndexRange(shape_index, 0, this->dofs_per_cell));
498 
499  // Return computed values if we
500  // know them easily. Otherwise, it
501  // is always safe to return true.
502  switch (rt_order)
503  {
504  case 0:
505  {
506  switch (dim)
507  {
508  case 2:
509  {
510  // only on the one
511  // non-adjacent face
512  // are the values
513  // actually zero. list
514  // these in a table
515  return (face_index !=
516  GeometryInfo<dim>::opposite_face[shape_index]);
517  }
518 
519  default:
520  return true;
521  }
522  }
523 
524  default: // other rt_order
525  return true;
526  }
527 
528  return true;
529 }
530 
531 
532 
533 template <int dim>
534 void
536  const std::vector<Vector<double>> &support_point_values,
537  std::vector<double> & nodal_values) const
538 {
539  Assert(support_point_values.size() == this->generalized_support_points.size(),
540  ExcDimensionMismatch(support_point_values.size(),
541  this->generalized_support_points.size()));
542  Assert(support_point_values[0].size() == this->n_components(),
543  ExcDimensionMismatch(support_point_values[0].size(),
544  this->n_components()));
545  Assert(nodal_values.size() == this->dofs_per_cell,
546  ExcDimensionMismatch(nodal_values.size(), this->dofs_per_cell));
547 
548  std::fill(nodal_values.begin(), nodal_values.end(), 0.);
549 
550  const unsigned int n_face_points = boundary_weights.size(0);
551  for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
552  for (unsigned int k = 0; k < n_face_points; ++k)
553  for (unsigned int i = 0; i < boundary_weights.size(1); ++i)
554  {
555  nodal_values[i + face * this->dofs_per_face] +=
556  boundary_weights(k, i) *
557  support_point_values[face * n_face_points + k][GeometryInfo<
558  dim>::unit_normal_direction[face]];
559  }
560 
561  const unsigned int start_cell_dofs =
562  GeometryInfo<dim>::faces_per_cell * this->dofs_per_face;
563  const unsigned int start_cell_points =
564  GeometryInfo<dim>::faces_per_cell * n_face_points;
565 
566  for (unsigned int k = 0; k < interior_weights.size(0); ++k)
567  for (unsigned int i = 0; i < interior_weights.size(1); ++i)
568  for (unsigned int d = 0; d < dim; ++d)
569  nodal_values[start_cell_dofs + i * dim + d] +=
570  interior_weights(k, i, d) *
571  support_point_values[k + start_cell_points][d];
572 
573  const unsigned int start_abf_dofs =
574  start_cell_dofs + interior_weights.size(1) * dim;
575 
576  // Cell integral of ABF terms
577  for (unsigned int k = 0; k < interior_weights_abf.size(0); ++k)
578  for (unsigned int i = 0; i < interior_weights_abf.size(1); ++i)
579  for (unsigned int d = 0; d < dim; ++d)
580  nodal_values[start_abf_dofs + i] +=
581  interior_weights_abf(k, i, d) *
582  support_point_values[k + start_cell_points][d];
583 
584  // Face integral of ABF terms
585  for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
586  {
587  const double n_orient = GeometryInfo<dim>::unit_normal_orientation[face];
588  for (unsigned int fp = 0; fp < n_face_points; ++fp)
589  {
590  // TODO: Check what the face_orientation, face_flip and face_rotation
591  // have to be in 3D
593  face, false, false, false, n_face_points);
594  for (unsigned int i = 0; i < boundary_weights_abf.size(1); ++i)
595  nodal_values[start_abf_dofs + i] +=
596  n_orient * boundary_weights_abf(k + fp, i) *
597  support_point_values[face * n_face_points + fp][GeometryInfo<
598  dim>::unit_normal_direction[face]];
599  }
600  }
601 
602  // TODO: Check if this "correction" can be removed.
603  for (unsigned int i = 0; i < boundary_weights_abf.size(1); ++i)
604  if (std::fabs(nodal_values[start_abf_dofs + i]) < 1.0e-16)
605  nodal_values[start_abf_dofs + i] = 0.0;
606 }
607 
608 
609 
610 template <int dim>
611 std::size_t
613 {
614  Assert(false, ExcNotImplemented());
615  return 0;
616 }
617 
618 
619 
620 /*-------------- Explicit Instantiations -------------------------------*/
621 #include "fe_abf.inst"
622 
623 DEAL_II_NAMESPACE_CLOSE
void initialize_restriction()
Definition: fe_abf.cc:321
static Quadrature< dim > project_to_child(const Quadrature< dim > &quadrature, const unsigned int child_no)
Definition: quadrature.cc:1063
void compute_embedding_matrices(const FiniteElement< dim, spacedim > &fe, std::vector< std::vector< FullMatrix< number >>> &matrices, const bool isotropic_only=false, const double threshold=1.e-12)
FullMatrix< double > interface_constraints
Definition: fe.h:2468
virtual void convert_generalized_support_point_values_to_dof_values(const std::vector< Vector< double >> &support_point_values, std::vector< double > &nodal_values) const override
Definition: fe_abf.cc:535
static void project_to_subface(const SubQuadrature &quadrature, const unsigned int face_no, const unsigned int subface_no, std::vector< Point< dim >> &q_points, const RefinementCase< dim - 1 > &ref_case=RefinementCase< dim - 1 >::isotropic_refinement)
virtual std::string get_name() const override
Definition: fe_abf.cc:116
static std::vector< Polynomial< double > > generate_complete_basis(const unsigned int degree)
Definition: polynomial.cc:866
const Point< dim > & point(const unsigned int i) const
friend class FE_ABF
Definition: fe_abf.h:251
STL namespace.
FullMatrix< double > compute_node_matrix(const FiniteElement< dim, spacedim > &fe)
static ::ExceptionBase & ExcIndexRange(int arg1, int arg2, int arg3)
virtual std::unique_ptr< FiniteElement< dim, dim > > clone() const override
Definition: fe_abf.cc:136
static ::ExceptionBase & ExcImpossibleInDim(int arg1)
std::vector< std::vector< FullMatrix< double > > > prolongation
Definition: fe.h:2456
void invert(const FullMatrix< number2 > &M)
void reinit(const TableIndices< N > &new_size, const bool omit_default_initialization=false)
#define Assert(cond, exc)
Definition: exceptions.h:1407
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static std::vector< Polynomial< number > > generate_complete_basis(const unsigned int degree)
Definition: polynomial.cc:686
unsigned int size() const
const unsigned int dofs_per_cell
Definition: fe_base.h:282
static Quadrature< dim > project_to_all_faces(const SubQuadrature &quadrature)
static void project_to_face(const SubQuadrature &quadrature, const unsigned int face_no, std::vector< Point< dim >> &q_points)
static unsigned int child_cell_on_face(const RefinementCase< dim > &ref_case, const unsigned int face, const unsigned int subface, const bool face_orientation=true, const bool face_flip=false, const bool face_rotation=false, const RefinementCase< dim - 1 > &face_refinement_case=RefinementCase< dim - 1 >::isotropic_refinement)
virtual std::size_t memory_consumption() const override
Definition: fe_abf.cc:612
const unsigned int dofs_per_face
Definition: fe_base.h:275
static std::vector< unsigned int > get_dpo_vector(const unsigned int degree)
Definition: fe_abf.cc:458
void reinit_restriction_and_prolongation_matrices(const bool isotropic_restriction_only=false, const bool isotropic_prolongation_only=false)
Definition: fe.cc:276
static ::ExceptionBase & ExcNotImplemented()
unsigned int n_components(const DoFHandler< dim, spacedim > &dh)
void initialize_support_points(const unsigned int rt_degree)
Definition: fe_abf.cc:151
double weight(const unsigned int i) const
virtual bool has_support_on_face(const unsigned int shape_index, const unsigned int face_index) const override
Definition: fe_abf.cc:491
static ::ExceptionBase & ExcInternalError()
static DataSetDescriptor face(const unsigned int face_no, const bool face_orientation, const bool face_flip, const bool face_rotation, const unsigned int n_quadrature_points)
Definition: quadrature.cc:1139