Reference documentation for deal.II version GIT fac6d58a2d 2022-05-23 17:55:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
mg_transfer_matrix_free.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2016 - 2022 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_mg_transfer_matrix_free_h
17 #define dealii_mg_transfer_matrix_free_h
18 
19 #include <deal.II/base/config.h>
20 
23 
25 
28 
33 
34 
36 
37 
40 
53 template <int dim, typename Number>
55  : public MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>
56 {
57 public:
63 
69 
73  virtual ~MGTransferMatrixFree() override = default;
74 
78  void
80 
84  void
85  clear();
86 
102  void
103  build(const DoFHandler<dim, dim> &dof_handler,
104  const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
105  &external_partitioners =
106  std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>());
107 
112  void
113  build(const DoFHandler<dim, dim> &dof_handler,
114  const std::function<void(const unsigned int,
117 
132  virtual void
133  prolongate(
134  const unsigned int to_level,
136  const LinearAlgebra::distributed::Vector<Number> &src) const override;
137 
138  virtual void
140  const unsigned int to_level,
142  const LinearAlgebra::distributed::Vector<Number> &src) const override;
143 
162  virtual void
164  const unsigned int from_level,
166  const LinearAlgebra::distributed::Vector<Number> &src) const override;
167 
182  template <typename Number2, int spacedim>
183  void
185  const DoFHandler<dim, spacedim> & dof_handler,
188 
193 
197  std::size_t
198  memory_consumption() const;
199 
200 private:
206  unsigned int fe_degree;
207 
213 
218  unsigned int n_components;
219 
225  unsigned int n_child_cell_dofs;
226 
237  std::vector<std::vector<unsigned int>> level_dof_indices;
238 
243  std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
245 
250  std::vector<unsigned int> n_owned_level_cells;
251 
257 
262 
274  std::vector<AlignedVector<VectorizedArray<Number>>> weights_on_refined;
275 
281  std::vector<std::vector<std::vector<unsigned short>>> dirichlet_indices;
282 
291 
295  template <int degree>
296  void
298  const unsigned int to_level,
301 
305  template <int degree>
306  void
307  do_restrict_add(const unsigned int from_level,
310 };
311 
312 
313 
319 template <int dim, typename Number, typename TransferType>
321  : public MGTransferBase<LinearAlgebra::distributed::BlockVector<Number>>
322 {
323 public:
326  {}
327 
342  virtual void
344  const unsigned int to_level,
346  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
347 
348  virtual void
350  const unsigned int to_level,
352  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
353 
372  virtual void
374  const unsigned int from_level,
376  const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
377 
388  template <typename Number2, int spacedim>
389  void
391  const DoFHandler<dim, spacedim> & dof_handler,
394 
398  template <typename Number2, int spacedim>
399  void
401  const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
404 
408  template <typename Number2, int spacedim>
409  void
411  const DoFHandler<dim, spacedim> & dof_handler,
414  const;
415 
419  template <typename Number2, int spacedim>
420  void
422  const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
425  const;
426 
431  static const bool supports_dof_handler_vector = true;
432 
433 protected:
438  virtual const TransferType &
439  get_matrix_free_transfer(const unsigned int b) const = 0;
440 
445  const bool same_for_all;
446 };
447 
448 
449 
463 template <int dim, typename Number>
465  : public MGTransferBlockMatrixFreeBase<dim,
466  Number,
467  MGTransferMatrixFree<dim, Number>>
468 {
469 public:
475 
480  MGTransferBlockMatrixFree(const MGConstrainedDoFs &mg_constrained_dofs);
481 
486  const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
487 
491  virtual ~MGTransferBlockMatrixFree() override = default;
492 
496  void
497  initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs);
498 
502  void
504  const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
505 
509  void
510  clear();
511 
515  void
516  build(const DoFHandler<dim, dim> &dof_handler);
517 
521  void
522  build(const std::vector<const DoFHandler<dim, dim> *> &dof_handler);
523 
527  std::size_t
528  memory_consumption() const;
529 
530 protected:
532  get_matrix_free_transfer(const unsigned int b) const override;
533 
534 private:
538  std::vector<MGTransferMatrixFree<dim, Number>> matrix_free_transfer_vector;
539 };
540 
541 
545 //------------------------ templated functions -------------------------
546 #ifndef DOXYGEN
547 
548 
549 template <int dim, typename Number>
550 template <typename Number2, int spacedim>
551 void
553  const DoFHandler<dim, spacedim> & dof_handler,
556 {
557  const unsigned int min_level = dst.min_level();
558  const unsigned int max_level = dst.max_level();
559 
560  Assert(max_level == dof_handler.get_triangulation().n_global_levels() - 1,
562  max_level, dof_handler.get_triangulation().n_global_levels() - 1));
563 
564  const FiniteElement<dim, spacedim> &fe = dof_handler.get_fe();
565 
566  for (unsigned int level = min_level; level <= max_level; ++level)
567  if (dst[level].size() != dof_handler.n_dofs(level) ||
568  dst[level].locally_owned_size() !=
569  dof_handler.locally_owned_mg_dofs(level).n_elements())
570  dst[level].reinit(this->vector_partitioners[level]);
571 
572  // copy fine level vector to active cells in MG hierarchy
573  this->copy_to_mg(dof_handler, dst, src, true);
574 
575  // FIXME: maybe need to store hanging nodes constraints per level?
576  // MGConstrainedDoFs does NOT keep this info right now, only periodicity
577  // constraints...
578 
579  // do the transfer from level to level-1:
580  dst[max_level].update_ghost_values();
581 
582  for (unsigned int level = max_level; level > min_level; --level)
583  {
584  // auxiliary vector which always has ghost elements
585  const LinearAlgebra::distributed::Vector<Number> *input = nullptr;
587  if (dst[level].get_partitioner().get() ==
588  this->vector_partitioners[level].get())
589  input = &dst[level];
590  else
591  {
592  ghosted_fine.reinit(this->vector_partitioners[level]);
593  ghosted_fine.copy_locally_owned_data_from(dst[level]);
594  ghosted_fine.update_ghost_values();
595  input = &ghosted_fine;
596  }
597 
598  std::vector<Number> dof_values_coarse(fe.n_dofs_per_cell());
599  Vector<Number> dof_values_fine(fe.n_dofs_per_cell());
601  std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
602  for (const auto &cell : dof_handler.cell_iterators_on_level(level - 1))
603  if (cell->is_locally_owned_on_level())
604  {
605  // if we get to a cell without children (== active), we can
606  // skip it as there values should be already set by the
607  // equivalent of copy_to_mg()
608  if (cell->is_active())
609  continue;
610 
611  std::fill(dof_values_coarse.begin(), dof_values_coarse.end(), 0.);
612  for (unsigned int child = 0; child < cell->n_children(); ++child)
613  {
614  cell->child(child)->get_mg_dof_indices(dof_indices);
615 
617  this->mg_constrained_dofs, level, dof_indices);
618 
619  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
620  dof_values_fine(i) = (*input)(dof_indices[i]);
621  fe.get_restriction_matrix(child, cell->refinement_case())
622  .vmult(tmp, dof_values_fine);
623  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
624  if (fe.restriction_is_additive(i))
625  dof_values_coarse[i] += tmp[i];
626  else if (tmp(i) != 0.)
627  dof_values_coarse[i] = tmp[i];
628  }
629  cell->get_mg_dof_indices(dof_indices);
630  for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
631  if (dof_handler.locally_owned_mg_dofs(level - 1).is_element(
632  dof_indices[i]))
633  dst[level - 1](dof_indices[i]) = dof_values_coarse[i];
634  }
635 
636  dst[level - 1].update_ghost_values();
637  }
638 }
639 
640 
641 
642 template <int dim, typename Number, typename TransferType>
643 template <typename Number2, int spacedim>
644 void
646  const DoFHandler<dim, spacedim> & dof_handler,
649 {
650  Assert(same_for_all,
651  ExcMessage(
652  "This object was initialized with support for usage with one "
653  "DoFHandler for each block, but this method assumes that "
654  "the same DoFHandler is used for all the blocks!"));
655  const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(src.n_blocks(),
656  &dof_handler);
657 
658  copy_to_mg(mg_dofs, dst, src);
659 }
660 
661 
662 
663 template <int dim, typename Number, typename TransferType>
664 template <typename Number2, int spacedim>
665 void
667  const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
670 {
671  const unsigned int n_blocks = src.n_blocks();
672  AssertDimension(dof_handler.size(), n_blocks);
673 
674  if (n_blocks == 0)
675  return;
676 
677  const unsigned int min_level = dst.min_level();
678  const unsigned int max_level = dst.max_level();
679 
680  for (unsigned int level = min_level; level <= max_level; ++level)
681  if (dst[level].n_blocks() != n_blocks)
682  dst[level].reinit(n_blocks);
683 
684  // FIXME: this a quite ugly as we need a temporary object:
686  min_level, max_level);
687 
688  for (unsigned int b = 0; b < n_blocks; ++b)
689  {
690  const unsigned int data_block = same_for_all ? 0 : b;
691  get_matrix_free_transfer(data_block)
692  .copy_to_mg(*dof_handler[b], dst_non_block, src.block(b));
693 
694  for (unsigned int l = min_level; l <= max_level; ++l)
695  dst[l].block(b) = dst_non_block[l];
696  }
697 
698  for (unsigned int level = min_level; level <= max_level; ++level)
699  dst[level].collect_sizes();
700 }
701 
702 template <int dim, typename Number, typename TransferType>
703 template <typename Number2, int spacedim>
704 void
706  const DoFHandler<dim, spacedim> & dof_handler,
709  const
710 {
711  Assert(same_for_all,
712  ExcMessage(
713  "This object was initialized with support for usage with one "
714  "DoFHandler for each block, but this method assumes that "
715  "the same DoFHandler is used for all the blocks!"));
716  const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(dst.n_blocks(),
717  &dof_handler);
718 
719  copy_from_mg(mg_dofs, dst, src);
720 }
721 
722 template <int dim, typename Number, typename TransferType>
723 template <typename Number2, int spacedim>
724 void
726  const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
729  const
730 {
731  const unsigned int n_blocks = dst.n_blocks();
732  AssertDimension(dof_handler.size(), n_blocks);
733 
734  if (n_blocks == 0)
735  return;
736 
737  const unsigned int min_level = src.min_level();
738  const unsigned int max_level = src.max_level();
739 
740  for (unsigned int l = min_level; l <= max_level; ++l)
741  AssertDimension(src[l].n_blocks(), dst.n_blocks());
742 
743  // FIXME: this a quite ugly as we need a temporary object:
745  min_level, max_level);
746 
747  for (unsigned int b = 0; b < n_blocks; ++b)
748  {
749  for (unsigned int l = min_level; l <= max_level; ++l)
750  {
751  src_non_block[l].reinit(src[l].block(b));
752  src_non_block[l] = src[l].block(b);
753  }
754  const unsigned int data_block = same_for_all ? 0 : b;
755  get_matrix_free_transfer(data_block)
756  .copy_from_mg(*dof_handler[b], dst.block(b), src_non_block);
757  }
758 }
759 
760 
761 
762 #endif // DOXYGEN
763 
764 
766 
767 #endif
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
const FiniteElement< dim, spacedim > & get_fe(const unsigned int index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
unsigned int n_dofs_per_cell() const
bool restriction_is_additive(const unsigned int index) const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
size_type n_elements() const
Definition: index_set.h:1834
bool is_element(const size_type index) const
Definition: index_set.h:1767
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
Definition: mg_transfer.h:586
std::function< void(const unsigned int, LinearAlgebra::distributed::Vector< Number > &)> initialize_dof_vector
Definition: mg_transfer.h:619
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, LinearAlgebra::distributed::BlockVector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &src) const
virtual const TransferType & get_matrix_free_transfer(const unsigned int b) const =0
void copy_from_mg(const std::vector< const DoFHandler< dim, spacedim > * > &dof_handler, LinearAlgebra::distributed::BlockVector< Number2 > &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &src) const
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
void copy_to_mg(const std::vector< const DoFHandler< dim, spacedim > * > &dof_handler, MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &dst, const LinearAlgebra::distributed::BlockVector< Number2 > &src) const
MGTransferBlockMatrixFreeBase(const bool same_for_all)
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::BlockVector< Number >> &dst, const LinearAlgebra::distributed::BlockVector< Number2 > &src) const
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
std::size_t memory_consumption() const
MGTransferBlockMatrixFree()=default
virtual ~MGTransferBlockMatrixFree() override=default
const MGTransferMatrixFree< dim, Number > & get_matrix_free_transfer(const unsigned int b) const override
void build(const DoFHandler< dim, dim > &dof_handler)
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
void do_restrict_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
void do_prolongate_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const
virtual ~MGTransferMatrixFree() override=default
std::vector< std::vector< unsigned int > > level_dof_indices
virtual void prolongate(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
void initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs)
std::vector< std::vector< std::vector< unsigned short > > > dirichlet_indices
MGLevelObject< std::shared_ptr< const Utilities::MPI::Partitioner > > vector_partitioners
void build(const DoFHandler< dim, dim > &dof_handler, const std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >> &external_partitioners=std::vector< std::shared_ptr< const Utilities::MPI::Partitioner >>())
void interpolate_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::Vector< Number >> &dst, const LinearAlgebra::distributed::Vector< Number2 > &src) const
AlignedVector< VectorizedArray< Number > > evaluation_data
std::vector< std::vector< std::pair< unsigned int, unsigned int > > > parent_child_connect
virtual void prolongate_and_add(const unsigned int to_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
std::size_t memory_consumption() const
std::vector< AlignedVector< VectorizedArray< Number > > > weights_on_refined
std::vector< unsigned int > n_owned_level_cells
AlignedVector< VectorizedArray< Number > > prolongation_matrix_1d
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::Vector< Number > &dst, const LinearAlgebra::distributed::Vector< Number > &src) const override
virtual unsigned int n_global_levels() const
Definition: vector.h:109
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:416
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:417
unsigned int level
Definition: grid_out.cc:4606
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
Definition: exceptions.h:464
static ::ExceptionBase & ExcNoProlongation()
unsigned int n_blocks() const
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
#define Assert(cond, exc)
Definition: exceptions.h:1473
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1667
BlockType & block(const unsigned int i)
static ::ExceptionBase & ExcMessage(std::string arg1)
void copy_locally_owned_data_from(const Vector< Number2, MemorySpace > &src)
void reinit(const size_type size, const bool omit_zeroing_entries=false)
void update_ghost_values() const
std::enable_if< IsBlockVector< VectorType >::value, unsigned int >::type n_blocks(const VectorType &vector)
Definition: operators.h:50
std::enable_if< IsBlockVector< VectorType >::value, void >::type collect_sizes(VectorType &vector)
Definition: operators.h:98
Tensor< 2, dim, Number > l(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
void resolve_identity_constraints(const MGConstrainedDoFs *mg_constrained_dofs, const unsigned int level, std::vector< types::global_dof_index > &dof_indices)