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