Reference documentation for deal.II version Git 193422c69f 2020-07-08 17:07:46 +0200
\(\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\}}\)
mg_transfer_component.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2003 - 2020 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/function.h>
18 #include <deal.II/base/logstream.h>
19 
21 #include <deal.II/dofs/dof_tools.h>
22 
23 #include <deal.II/fe/fe.h>
24 
25 #include <deal.II/grid/tria.h>
27 
32 #include <deal.II/lac/vector.h>
33 
36 #include <deal.II/multigrid/mg_transfer_component.templates.h>
37 
38 #include <algorithm>
39 #include <iostream>
40 #include <numeric>
41 
43 
44 
45 namespace
46 {
79  template <int dim, typename number, int spacedim>
80  void
81  reinit_vector_by_components(
82  const ::DoFHandler<dim, spacedim> & mg_dof,
84  const std::vector<bool> & sel,
85  const std::vector<unsigned int> & target_comp,
86  std::vector<std::vector<types::global_dof_index>> &ndofs)
87  {
88  std::vector<bool> selected = sel;
89  std::vector<unsigned int> target_component = target_comp;
90  const unsigned int ncomp = mg_dof.get_fe(0).n_components();
91 
92  // If the selected and
93  // target_component have size 0,
94  // they must be replaced by default
95  // values.
96  //
97  // Since we already made copies
98  // directly after this function was
99  // called, we use the arguments
100  // directly.
101  if (target_component.size() == 0)
102  {
103  target_component.resize(ncomp);
104  for (unsigned int i = 0; i < ncomp; ++i)
105  target_component[i] = i;
106  }
107 
108  // If selected is an empty vector,
109  // all components are selected.
110  if (selected.size() == 0)
111  {
112  selected.resize(target_component.size());
113  std::fill_n(selected.begin(), ncomp, false);
114  for (const unsigned int component : target_component)
115  selected[component] = true;
116  }
117 
118  Assert(selected.size() == target_component.size(),
119  ExcDimensionMismatch(selected.size(), target_component.size()));
120 
121  // Compute the number of blocks needed
122  const unsigned int n_selected =
123  std::accumulate(selected.begin(), selected.end(), 0u);
124 
125  if (ndofs.size() == 0)
126  {
127  std::vector<std::vector<types::global_dof_index>> new_dofs(
128  mg_dof.get_triangulation().n_levels(),
129  std::vector<types::global_dof_index>(target_component.size()));
130  std::swap(ndofs, new_dofs);
131  MGTools::count_dofs_per_block(mg_dof, ndofs, target_component);
132  }
133 
134  for (unsigned int level = v.min_level(); level <= v.max_level(); ++level)
135  {
136  v[level].reinit(n_selected, 0);
137  unsigned int k = 0;
138  for (unsigned int i = 0;
139  i < selected.size() && (k < v[level].n_blocks());
140  ++i)
141  {
142  if (selected[i])
143  {
144  v[level].block(k++).reinit(ndofs[level][i]);
145  }
146  v[level].collect_sizes();
147  }
148  }
149  }
150 
151 
178  template <int dim, typename number, int spacedim>
179  void
180  reinit_vector_by_components(
181  const ::DoFHandler<dim, spacedim> & mg_dof,
183  const ComponentMask & component_mask,
184  const std::vector<unsigned int> & target_component,
185  std::vector<std::vector<types::global_dof_index>> &ndofs)
186  {
187  Assert(component_mask.represents_n_components(target_component.size()),
188  ExcMessage("The component mask does not have the correct size."));
189 
190  unsigned int selected_block = 0;
191  for (unsigned int i = 0; i < target_component.size(); ++i)
192  if (component_mask[i])
193  selected_block = target_component[i];
194 
195  if (ndofs.size() == 0)
196  {
197  std::vector<std::vector<types::global_dof_index>> new_dofs(
198  mg_dof.get_triangulation().n_levels(),
199  std::vector<types::global_dof_index>(target_component.size()));
200  std::swap(ndofs, new_dofs);
201  MGTools::count_dofs_per_block(mg_dof, ndofs, target_component);
202  }
203 
204  for (unsigned int level = v.min_level(); level <= v.max_level(); ++level)
205  {
206  v[level].reinit(ndofs[level][selected_block]);
207  }
208  }
209 } // namespace
210 
211 
212 template <typename number>
213 template <int dim, class InVector, int spacedim>
214 void
216  const DoFHandler<dim, spacedim> &mg_dof_handler,
218  const InVector & src) const
219 {
220  dst = 0;
221 
222  Assert(sizes.size() == mg_dof_handler.get_triangulation().n_levels(),
223  ExcMatricesNotBuilt());
224 
225  reinit_vector_by_components(
226  mg_dof_handler, dst, mg_component_mask, mg_target_component, sizes);
227 
228  // traverse the grid top-down
229  // (i.e. starting with the most
230  // refined grid). this way, we can
231  // always get that part of one
232  // level of the output vector which
233  // corresponds to a region which is
234  // more refined, by restriction of
235  // the respective vector on the
236  // next finer level, which we then
237  // already have built.
238 
239  for (unsigned int level = mg_dof_handler.get_triangulation().n_levels();
240  level != 0;)
241  {
242  --level;
243 
244  using IT = std::vector<
245  std::pair<types::global_dof_index, unsigned int>>::const_iterator;
246  for (IT i = copy_to_and_from_indices[level].begin();
247  i != copy_to_and_from_indices[level].end();
248  ++i)
249  dst[level](i->second) = src(i->first);
250  }
251 }
252 
253 
254 
255 template <int dim, int spacedim>
256 void
258  const DoFHandler<dim, spacedim> &mg_dof)
259 {
260  build(mg_dof);
261 }
262 
263 
264 
265 template <int dim, int spacedim>
266 void
268 {
269  // Fill target component with
270  // standard values (identity) if it
271  // is empty
272  if (target_component.size() == 0)
273  {
274  target_component.resize(mg_dof.get_fe(0).n_components());
275  for (unsigned int i = 0; i < target_component.size(); ++i)
276  target_component[i] = i;
277  }
278  else
279  {
280  // otherwise, check it for consistency
281  Assert(target_component.size() == mg_dof.get_fe(0).n_components(),
282  ExcDimensionMismatch(target_component.size(),
283  mg_dof.get_fe(0).n_components()));
284 
285  for (unsigned int i = 0; i < target_component.size(); ++i)
286  {
287  AssertIndexRange(i, target_component.size());
288  }
289  }
290  // Do the same for the multilevel
291  // components. These may be
292  // different.
293  if (mg_target_component.size() == 0)
294  {
295  mg_target_component.resize(mg_dof.get_fe(0).n_components());
296  for (unsigned int i = 0; i < mg_target_component.size(); ++i)
297  mg_target_component[i] = target_component[i];
298  }
299  else
300  {
301  Assert(mg_target_component.size() == mg_dof.get_fe(0).n_components(),
303  mg_dof.get_fe(0).n_components()));
304 
305  for (unsigned int i = 0; i < mg_target_component.size(); ++i)
306  {
308  }
309  }
310 
311  const FiniteElement<dim> &fe = mg_dof.get_fe();
312 
313  // Effective number of components
314  // is the maximum entry in
315  // mg_target_component. This
316  // assumes that the values in that
317  // vector don't have holes.
318  const unsigned int n_components =
319  *std::max_element(mg_target_component.begin(), mg_target_component.end()) +
320  1;
321  const unsigned int dofs_per_cell = fe.dofs_per_cell;
322  const unsigned int n_levels = mg_dof.get_triangulation().n_levels();
323 
325  ExcMessage("Component mask has wrong size."));
326 
327  // Compute the lengths of all blocks
328  sizes.resize(n_levels);
329  for (unsigned int l = 0; l < n_levels; ++l)
330  sizes[l].resize(n_components);
331 
333 
334  // Fill some index vectors
335  // for later use.
337  // Compute start indices from sizes
338  for (auto &level_components : mg_component_start)
339  {
341  for (types::global_dof_index &level_component_start : level_components)
342  {
343  const types::global_dof_index t = level_component_start;
344  level_component_start = k;
345  k += t;
346  }
347  }
348 
349  component_start = DoFTools::count_dofs_per_fe_block(mg_dof, target_component);
350 
352  for (types::global_dof_index &first_index : component_start)
353  {
354  const types::global_dof_index t = first_index;
355  first_index = k;
356  k += t;
357  }
358 
359  // Build index vectors for
360  // copy_to_mg and
361  // copy_from_mg. These vectors must
362  // be prebuilt, since the
363  // get_dof_indices functions are
364  // too slow
365 
366  copy_to_and_from_indices.resize(n_levels);
367 
368  // Building the prolongation matrices starts here!
369 
370  // reset the size of the array of
371  // matrices. call resize(0) first,
372  // in order to delete all elements
373  // and clear their memory. then
374  // repopulate these arrays
375  //
376  // note that on resize(0), the
377  // shared_ptr class takes care of
378  // deleting the object it points to
379  // by itself
380  prolongation_matrices.resize(0);
381  prolongation_sparsities.resize(0);
382  prolongation_matrices.reserve(n_levels - 1);
383  prolongation_sparsities.reserve(n_levels - 1);
384 
385  for (unsigned int i = 0; i < n_levels - 1; ++i)
386  {
389  }
390 
391  // two fields which will store the
392  // indices of the multigrid dofs
393  // for a cell and one of its children
394  std::vector<types::global_dof_index> dof_indices_parent(dofs_per_cell);
395  std::vector<types::global_dof_index> dof_indices_child(dofs_per_cell);
396 
397  // for each level: first build the
398  // sparsity pattern of the matrices
399  // and then build the matrices
400  // themselves. note that we only
401  // need to take care of cells on
402  // the coarser level which have
403  // children
404  for (unsigned int level = 0; level < n_levels - 1; ++level)
405  {
406  // reset the dimension of the
407  // structure. note that for
408  // the number of entries per
409  // row, the number of parent
410  // dofs coupling to a child dof
411  // is necessary. this, is the
412  // number of degrees of freedom
413  // per cell
414  prolongation_sparsities[level]->reinit(n_components, n_components);
415  for (unsigned int i = 0; i < n_components; ++i)
416  for (unsigned int j = 0; j < n_components; ++j)
417  if (i == j)
418  prolongation_sparsities[level]->block(i, j).reinit(
419  sizes[level + 1][i], sizes[level][j], dofs_per_cell + 1);
420  else
421  prolongation_sparsities[level]->block(i, j).reinit(
422  sizes[level + 1][i], sizes[level][j], 0);
423 
424  prolongation_sparsities[level]->collect_sizes();
425 
426  for (typename DoFHandler<dim, spacedim>::cell_iterator cell =
427  mg_dof.begin(level);
428  cell != mg_dof.end(level);
429  ++cell)
430  if (cell->has_children())
431  {
432  cell->get_mg_dof_indices(dof_indices_parent);
433 
434  for (unsigned int child = 0; child < cell->n_children(); ++child)
435  {
436  // set an alias to the
437  // prolongation matrix for
438  // this child
439  const FullMatrix<double> &prolongation =
440  mg_dof.get_fe().get_prolongation_matrix(
441  child, cell->refinement_case());
442 
443  cell->child(child)->get_mg_dof_indices(dof_indices_child);
444 
445  // now tag the entries in the
446  // matrix which will be used
447  // for this pair of parent/child
448  for (unsigned int i = 0; i < dofs_per_cell; ++i)
449  for (unsigned int j = 0; j < dofs_per_cell; ++j)
450  if (prolongation(i, j) != 0)
451  {
452  const unsigned int icomp =
453  fe.system_to_component_index(i).first;
454  const unsigned int jcomp =
455  fe.system_to_component_index(j).first;
456  if ((icomp == jcomp) && mg_component_mask[icomp])
458  dof_indices_child[i], dof_indices_parent[j]);
459  }
460  }
461  }
462  prolongation_sparsities[level]->compress();
463 
465  // now actually build the matrices
466  for (typename DoFHandler<dim, spacedim>::cell_iterator cell =
467  mg_dof.begin(level);
468  cell != mg_dof.end(level);
469  ++cell)
470  if (cell->has_children())
471  {
472  cell->get_mg_dof_indices(dof_indices_parent);
473 
474  for (unsigned int child = 0; child < cell->n_children(); ++child)
475  {
476  // set an alias to the
477  // prolongation matrix for
478  // this child
479  const FullMatrix<double> &prolongation =
480  mg_dof.get_fe().get_prolongation_matrix(
481  child, cell->refinement_case());
482 
483  cell->child(child)->get_mg_dof_indices(dof_indices_child);
484 
485  // now set the entries in the
486  // matrix
487  for (unsigned int i = 0; i < dofs_per_cell; ++i)
488  for (unsigned int j = 0; j < dofs_per_cell; ++j)
489  if (prolongation(i, j) != 0)
490  {
491  const unsigned int icomp =
492  fe.system_to_component_index(i).first;
493  const unsigned int jcomp =
494  fe.system_to_component_index(j).first;
495  if ((icomp == jcomp) && mg_component_mask[icomp])
497  dof_indices_child[i],
498  dof_indices_parent[j],
499  prolongation(i, j));
500  }
501  }
502  }
503  }
504  // impose boundary conditions
505  // but only in the column of
506  // the prolongation matrix
507  // TODO: this way is not very efficient
508 
509  if (boundary_indices.size() != 0)
510  {
511  std::vector<std::vector<types::global_dof_index>> dofs_per_component(
512  mg_dof.get_triangulation().n_levels(),
513  std::vector<types::global_dof_index>(n_components));
514 
516  dofs_per_component,
518  for (unsigned int level = 0; level < n_levels - 1; ++level)
519  {
520  if (boundary_indices[level].size() == 0)
521  continue;
522 
523  for (unsigned int iblock = 0; iblock < n_components; ++iblock)
524  for (unsigned int jblock = 0; jblock < n_components; ++jblock)
525  if (iblock == jblock)
526  {
527  const types::global_dof_index n_dofs =
528  prolongation_matrices[level]->block(iblock, jblock).m();
529  for (types::global_dof_index i = 0; i < n_dofs; ++i)
530  {
533  ->block(iblock, jblock)
534  .begin(i),
536  ->block(iblock, jblock)
537  .end(i);
538  for (; begin != end; ++begin)
539  {
540  const types::global_dof_index column_number =
541  begin->column();
542 
543  // convert global indices into local ones
544  const BlockIndices block_indices_coarse(
545  dofs_per_component[level]);
546  const types::global_dof_index global_j =
547  block_indices_coarse.local_to_global(iblock,
548  column_number);
549 
550  std::set<types::global_dof_index>::const_iterator
551  found_dof = boundary_indices[level].find(global_j);
552 
553  const bool is_boundary_index =
554  (found_dof != boundary_indices[level].end());
555 
556  if (is_boundary_index)
557  {
559  ->block(iblock, jblock)
560  .set(i, column_number, 0);
561  }
562  }
563  }
564  }
565  }
566  }
567 }
568 
569 
570 
571 template <typename number>
572 template <int dim, int spacedim>
573 void
575  const DoFHandler<dim, spacedim> & /*dof*/,
576  const DoFHandler<dim, spacedim> & mg_dof,
577  unsigned int select,
578  unsigned int mg_select,
579  const std::vector<unsigned int> & t_component,
580  const std::vector<unsigned int> & mg_t_component,
581  const std::vector<std::set<types::global_dof_index>> &bdry_indices)
582 {
583  build(mg_dof, select, mg_select, t_component, mg_t_component, bdry_indices);
584 }
585 
586 
587 
588 template <typename number>
589 template <int dim, int spacedim>
590 void
592  const DoFHandler<dim, spacedim> & mg_dof,
593  unsigned int select,
594  unsigned int mg_select,
595  const std::vector<unsigned int> & t_component,
596  const std::vector<unsigned int> & mg_t_component,
597  const std::vector<std::set<types::global_dof_index>> &bdry_indices)
598 {
599  const FiniteElement<dim> &fe = mg_dof.get_fe();
600  unsigned int ncomp = mg_dof.get_fe(0).n_components();
601 
602  target_component = t_component;
603  mg_target_component = mg_t_component;
604  boundary_indices = bdry_indices;
605 
606  selected_component = select;
607  mg_selected_component = mg_select;
608 
609  {
610  std::vector<bool> tmp(ncomp, false);
611  for (unsigned int c = 0; c < ncomp; ++c)
612  if (t_component[c] == selected_component)
613  tmp[c] = true;
614  component_mask = ComponentMask(tmp);
615  }
616 
617  {
618  std::vector<bool> tmp(ncomp, false);
619  for (unsigned int c = 0; c < ncomp; ++c)
620  if (mg_t_component[c] == mg_selected_component)
621  tmp[c] = true;
623  }
624 
625  // If components are renumbered,
626  // find the first original
627  // component corresponding to the
628  // target component.
629  for (unsigned int i = 0; i < target_component.size(); ++i)
630  {
631  if (target_component[i] == select)
632  {
633  selected_component = i;
634  break;
635  }
636  }
637 
638  for (unsigned int i = 0; i < mg_target_component.size(); ++i)
639  {
640  if (mg_target_component[i] == mg_select)
641  {
642  mg_selected_component = i;
643  break;
644  }
645  }
646 
648 
649  interface_dofs.resize(mg_dof.get_triangulation().n_levels());
650  for (unsigned int l = 0; l < mg_dof.get_triangulation().n_levels(); ++l)
651  {
652  interface_dofs[l].clear();
653  interface_dofs[l].set_size(mg_dof.n_dofs(l));
654  }
655  MGTools::extract_inner_interface_dofs(mg_dof, interface_dofs);
656 
657  // use a temporary vector to create the
658  // relation between global and level dofs
659  std::vector<types::global_dof_index> temp_copy_indices;
660  std::vector<types::global_dof_index> global_dof_indices(fe.dofs_per_cell);
661  std::vector<types::global_dof_index> level_dof_indices(fe.dofs_per_cell);
662  for (int level = mg_dof.get_triangulation().n_levels() - 1; level >= 0;
663  --level)
664  {
667  mg_dof.begin_active(level);
668  const typename DoFHandler<dim, spacedim>::active_cell_iterator level_end =
669  mg_dof.end_active(level);
670 
671  temp_copy_indices.resize(0);
672  temp_copy_indices.resize(mg_dof.n_dofs(level),
674 
675  // Compute coarse level right hand side
676  // by restricting from fine level.
677  for (; level_cell != level_end; ++level_cell)
678  {
679  // get the dof numbers of
680  // this cell for the global
681  // and the level-wise
682  // numbering
683  level_cell->get_dof_indices(global_dof_indices);
684  level_cell->get_mg_dof_indices(level_dof_indices);
685 
686  for (unsigned int i = 0; i < fe.dofs_per_cell; ++i)
687  {
688  const unsigned int component =
689  fe.system_to_component_index(i).first;
690  if (component_mask[component] &&
691  !interface_dofs[level].is_element(level_dof_indices[i]))
692  {
693  const types::global_dof_index level_start =
695  const types::global_dof_index global_start =
696  component_start[target_component[component]];
697  temp_copy_indices[level_dof_indices[i] - level_start] =
698  global_dof_indices[i] - global_start;
699  }
700  }
701  }
702 
703  // write indices from vector into the map from
704  // global to level dofs
705  const types::global_dof_index n_active_dofs =
706  std::count_if(temp_copy_indices.begin(),
707  temp_copy_indices.end(),
708  [](const types::global_dof_index index) {
709  return index != numbers::invalid_dof_index;
710  });
711  copy_to_and_from_indices[level].resize(n_active_dofs);
712  types::global_dof_index counter = 0;
713  for (types::global_dof_index i = 0; i < temp_copy_indices.size(); ++i)
714  if (temp_copy_indices[i] != numbers::invalid_dof_index)
715  copy_to_and_from_indices[level][counter++] =
716  std::pair<types::global_dof_index, unsigned int>(
717  temp_copy_indices[i], i);
718  Assert(counter == n_active_dofs, ExcInternalError());
719  }
720 }
721 
722 
723 
724 // explicit instantiations
725 #include "mg_transfer_component.inst"
726 
727 
void build(const DoFHandler< dim, spacedim > &dof_handler)
const Triangulation< dim, spacedim > & get_triangulation() const
std::vector< std::shared_ptr< BlockSparseMatrix< double > > > prolongation_matrices
void count_dofs_per_block(const DoFHandler< dim, spacedim > &dof_handler, std::vector< std::vector< types::global_dof_index >> &dofs_per_block, std::vector< unsigned int > target_block={})
Definition: mg_tools.cc:1158
cell_iterator begin(const unsigned int level=0) const
std::vector< std::set< types::global_dof_index > > boundary_indices
cell_iterator end() const
#define AssertIndexRange(index, range)
Definition: exceptions.h:1628
bool represents_n_components(const unsigned int n) const
std::vector< types::global_dof_index > count_dofs_per_fe_block(const DoFHandler< dim, spacedim > &dof, const std::vector< unsigned int > &target_block=std::vector< unsigned int >())
Definition: dof_tools.cc:1998
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
active_cell_iterator begin_active(const unsigned int level=0) const
static ::ExceptionBase & ExcMessage(std::string arg1)
void build_matrices(const DoFHandler< dim, spacedim > &dof, const DoFHandler< dim, spacedim > &mg_dof, unsigned int selected, unsigned int mg_selected, const std::vector< unsigned int > &target_component=std::vector< unsigned int >(), const std::vector< unsigned int > &mg_target_component=std::vector< unsigned int >(), const std::vector< std::set< types::global_dof_index >> &boundary_indices=std::vector< std::set< types::global_dof_index >>())
#define Assert(cond, exc)
Definition: exceptions.h:1403
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
std::vector< types::global_dof_index > component_start
active_cell_iterator end_active(const unsigned int level) const
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:363
std::vector< std::vector< std::pair< types::global_dof_index, unsigned int > > > copy_to_and_from_indices
unsigned int level
Definition: grid_out.cc:4355
types::global_dof_index n_dofs() const
VectorType::value_type * end(VectorType &V)
size_type local_to_global(const unsigned int block, const size_type index) const
void build(const DoFHandler< dim, spacedim > &dof, unsigned int selected, unsigned int mg_selected, const std::vector< unsigned int > &target_component=std::vector< unsigned int >(), const std::vector< unsigned int > &mg_target_component=std::vector< unsigned int >(), const std::vector< std::set< types::global_dof_index >> &boundary_indices=std::vector< std::set< types::global_dof_index >>())
const unsigned int dofs_per_cell
Definition: fe_base.h:280
std::vector< std::shared_ptr< BlockSparsityPattern > > prolongation_sparsities
void do_copy_to_mg(const DoFHandler< dim, spacedim > &mg_dof, MGLevelObject< Vector< number >> &dst, const InVector &src) const
std::vector< unsigned int > mg_target_component
void swap(MemorySpaceData< Number, MemorySpace > &, MemorySpaceData< Number, MemorySpace > &)
Definition: memory_space.h:102
std::pair< unsigned int, unsigned int > system_to_component_index(const unsigned int index) const
Definition: fe.h:3062
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:362
VectorType::value_type * begin(VectorType &V)
void extract_inner_interface_dofs(const DoFHandler< dim, spacedim > &mg_dof_handler, std::vector< IndexSet > &interface_dofs)
Definition: mg_tools.cc:1526
typename ActiveSelector::cell_iterator cell_iterator
Definition: dof_handler.h:371
std::vector< std::vector< types::global_dof_index > > mg_component_start
typename ActiveSelector::active_cell_iterator active_cell_iterator
Definition: dof_handler.h:343
const types::global_dof_index invalid_dof_index
Definition: types.h:206
void build_matrices(const DoFHandler< dim, spacedim > &dof, const DoFHandler< dim, spacedim > &mg_dof)
std::vector< std::vector< types::global_dof_index > > sizes
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
static ::ExceptionBase & ExcInternalError()