Reference documentation for deal.II version 9.5.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\}}\)
Loading...
Searching...
No Matches
mg_transfer_matrix_free.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2016 - 2023 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
55template <int dim, typename Number>
57 : public MGLevelGlobalTransfer<LinearAlgebra::distributed::Vector<Number>>
58{
59public:
65
71
75 virtual ~MGTransferMatrixFree() override = default;
76
80 void
82
86 void
87 clear();
88
104 void
105 build(const DoFHandler<dim, dim> &dof_handler,
106 const std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>
107 &external_partitioners =
108 std::vector<std::shared_ptr<const Utilities::MPI::Partitioner>>());
109
114 void
115 build(const DoFHandler<dim, dim> &dof_handler,
116 const std::function<void(const unsigned int,
119
134 virtual void
136 const unsigned int to_level,
138 const LinearAlgebra::distributed::Vector<Number> &src) const override;
139
140 virtual void
142 const unsigned int to_level,
144 const LinearAlgebra::distributed::Vector<Number> &src) const override;
145
164 virtual void
166 const unsigned int from_level,
168 const LinearAlgebra::distributed::Vector<Number> &src) const override;
169
184 template <typename BlockVectorType2, int spacedim>
185 void
187 const DoFHandler<dim, spacedim> & dof_handler,
189 const BlockVectorType2 & src) const;
190
195
199 std::size_t
200 memory_consumption() const;
201
202private:
208 unsigned int fe_degree;
209
215
220 unsigned int n_components;
221
227 unsigned int n_child_cell_dofs;
228
239 std::vector<std::vector<unsigned int>> level_dof_indices;
240
245 std::vector<std::vector<std::pair<unsigned int, unsigned int>>>
247
252 std::vector<unsigned int> n_owned_level_cells;
253
259
264
276 std::vector<AlignedVector<VectorizedArray<Number>>> weights_on_refined;
277
283 std::vector<std::vector<std::vector<unsigned short>>> dirichlet_indices;
284
293
297 template <int degree>
298 void
300 const unsigned int to_level,
303
307 template <int degree>
308 void
309 do_restrict_add(const unsigned int from_level,
312};
313
314
315
321template <int dim, typename Number, typename TransferType>
323 : public MGTransferBase<LinearAlgebra::distributed::BlockVector<Number>>
324{
325public:
328 {}
329
344 virtual void
346 const unsigned int to_level,
348 const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
349
350 virtual void
352 const unsigned int to_level,
354 const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
355
374 virtual void
376 const unsigned int from_level,
378 const LinearAlgebra::distributed::BlockVector<Number> &src) const override;
379
390 template <typename BlockVectorType2, int spacedim>
391 void
393 const DoFHandler<dim, spacedim> & dof_handler,
395 const BlockVectorType2 & src) const;
396
400 template <typename BlockVectorType2, int spacedim>
401 void
403 const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
405 const BlockVectorType2 & src) const;
406
410 template <typename BlockVectorType2, int spacedim>
411 void
413 const DoFHandler<dim, spacedim> &dof_handler,
414 BlockVectorType2 & dst,
416 const;
417
421 template <typename BlockVectorType2, int spacedim>
422 void
424 const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
425 BlockVectorType2 & dst,
427 const;
428
433 static const bool supports_dof_handler_vector = true;
434
435protected:
440 virtual const TransferType &
441 get_matrix_free_transfer(const unsigned int b) const = 0;
442
447 const bool same_for_all;
448};
449
450
451
465template <int dim, typename Number>
468 Number,
469 MGTransferMatrixFree<dim, Number>>
470{
471public:
477
482 MGTransferBlockMatrixFree(const MGConstrainedDoFs &mg_constrained_dofs);
483
488 const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
489
493 virtual ~MGTransferBlockMatrixFree() override = default;
494
498 void
499 initialize_constraints(const MGConstrainedDoFs &mg_constrained_dofs);
500
504 void
506 const std::vector<MGConstrainedDoFs> &mg_constrained_dofs);
507
511 void
512 clear();
513
517 void
518 build(const DoFHandler<dim, dim> &dof_handler);
519
523 void
524 build(const std::vector<const DoFHandler<dim, dim> *> &dof_handler);
525
529 std::size_t
530 memory_consumption() const;
531
532protected:
534 get_matrix_free_transfer(const unsigned int b) const override;
535
536private:
540 std::vector<MGTransferMatrixFree<dim, Number>> matrix_free_transfer_vector;
541};
542
543
547//------------------------ templated functions -------------------------
548#ifndef DOXYGEN
549
550
551template <int dim, typename Number>
552template <typename BlockVectorType2, int spacedim>
553void
555 const DoFHandler<dim, spacedim> & dof_handler,
557 const BlockVectorType2 & src) const
558{
559 const unsigned int min_level = dst.min_level();
560 const unsigned int max_level = dst.max_level();
561
562 Assert(max_level == dof_handler.get_triangulation().n_global_levels() - 1,
564 max_level, dof_handler.get_triangulation().n_global_levels() - 1));
565
566 const FiniteElement<dim, spacedim> &fe = dof_handler.get_fe();
567
568 for (unsigned int level = min_level; level <= max_level; ++level)
569 if (dst[level].size() != dof_handler.n_dofs(level) ||
570 dst[level].locally_owned_size() !=
572 dst[level].reinit(this->vector_partitioners[level]);
573
574 // copy fine level vector to active cells in MG hierarchy
575 this->copy_to_mg(dof_handler, dst, src, true);
576
577 // FIXME: maybe need to store hanging nodes constraints per level?
578 // MGConstrainedDoFs does NOT keep this info right now, only periodicity
579 // constraints...
580
581 // do the transfer from level to level-1:
582 dst[max_level].update_ghost_values();
583
584 for (unsigned int level = max_level; level > min_level; --level)
585 {
586 // auxiliary vector which always has ghost elements
587 const LinearAlgebra::distributed::Vector<Number> *input = nullptr;
589 if (dst[level].get_partitioner().get() ==
590 this->vector_partitioners[level].get())
591 input = &dst[level];
592 else
593 {
594 ghosted_fine.reinit(this->vector_partitioners[level]);
595 ghosted_fine.copy_locally_owned_data_from(dst[level]);
596 ghosted_fine.update_ghost_values();
597 input = &ghosted_fine;
598 }
599
600 std::vector<Number> dof_values_coarse(fe.n_dofs_per_cell());
601 Vector<Number> dof_values_fine(fe.n_dofs_per_cell());
603 std::vector<types::global_dof_index> dof_indices(fe.n_dofs_per_cell());
604 for (const auto &cell : dof_handler.cell_iterators_on_level(level - 1))
605 if (cell->is_locally_owned_on_level())
606 {
607 // if we get to a cell without children (== active), we can
608 // skip it as there values should be already set by the
609 // equivalent of copy_to_mg()
610 if (cell->is_active())
611 continue;
612
613 std::fill(dof_values_coarse.begin(), dof_values_coarse.end(), 0.);
614 for (unsigned int child = 0; child < cell->n_children(); ++child)
615 {
616 cell->child(child)->get_mg_dof_indices(dof_indices);
617
619 this->mg_constrained_dofs, level, dof_indices);
620
621 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
622 dof_values_fine(i) = (*input)(dof_indices[i]);
623 fe.get_restriction_matrix(child, cell->refinement_case())
624 .vmult(tmp, dof_values_fine);
625 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
626 if (fe.restriction_is_additive(i))
627 dof_values_coarse[i] += tmp[i];
628 else if (tmp(i) != 0.)
629 dof_values_coarse[i] = tmp[i];
630 }
631 cell->get_mg_dof_indices(dof_indices);
632 for (unsigned int i = 0; i < fe.n_dofs_per_cell(); ++i)
633 if (dof_handler.locally_owned_mg_dofs(level - 1).is_element(
634 dof_indices[i]))
635 dst[level - 1](dof_indices[i]) = dof_values_coarse[i];
636 }
637
638 dst[level - 1].update_ghost_values();
639 }
640}
641
642
643
644template <int dim, typename Number, typename TransferType>
645template <typename BlockVectorType2, int spacedim>
646void
648 const DoFHandler<dim, spacedim> & dof_handler,
650 const BlockVectorType2 & src) const
651{
652 Assert(same_for_all,
654 "This object was initialized with support for usage with one "
655 "DoFHandler for each block, but this method assumes that "
656 "the same DoFHandler is used for all the blocks!"));
657 const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(src.n_blocks(),
658 &dof_handler);
659
660 copy_to_mg(mg_dofs, dst, src);
661}
662
663
664
665template <int dim, typename Number, typename TransferType>
666template <typename BlockVectorType2, int spacedim>
667void
669 const std::vector<const DoFHandler<dim, spacedim> *> & dof_handler,
671 const BlockVectorType2 & src) const
672{
673 const unsigned int n_blocks = src.n_blocks();
674 AssertDimension(dof_handler.size(), n_blocks);
675
676 if (n_blocks == 0)
677 return;
678
679 const unsigned int min_level = dst.min_level();
680 const unsigned int max_level = dst.max_level();
681
682 for (unsigned int level = min_level; level <= max_level; ++level)
683 if (dst[level].n_blocks() != n_blocks)
684 dst[level].reinit(n_blocks);
685
686 // FIXME: this a quite ugly as we need a temporary object:
688 min_level, max_level);
689
690 for (unsigned int b = 0; b < n_blocks; ++b)
691 {
692 const unsigned int data_block = same_for_all ? 0 : b;
693 get_matrix_free_transfer(data_block)
694 .copy_to_mg(*dof_handler[b], dst_non_block, src.block(b));
695
696 for (unsigned int l = min_level; l <= max_level; ++l)
697 dst[l].block(b) = dst_non_block[l];
698 }
699
700 for (unsigned int level = min_level; level <= max_level; ++level)
701 dst[level].collect_sizes();
702}
703
704template <int dim, typename Number, typename TransferType>
705template <typename BlockVectorType2, int spacedim>
706void
708 const DoFHandler<dim, spacedim> &dof_handler,
709 BlockVectorType2 & dst,
711 const
712{
713 Assert(same_for_all,
715 "This object was initialized with support for usage with one "
716 "DoFHandler for each block, but this method assumes that "
717 "the same DoFHandler is used for all the blocks!"));
718 const std::vector<const DoFHandler<dim, spacedim> *> mg_dofs(dst.n_blocks(),
719 &dof_handler);
720
721 copy_from_mg(mg_dofs, dst, src);
722}
723
724template <int dim, typename Number, typename TransferType>
725template <typename BlockVectorType2, int spacedim>
726void
728 const std::vector<const DoFHandler<dim, spacedim> *> &dof_handler,
729 BlockVectorType2 & dst,
731 const
732{
733 const unsigned int n_blocks = dst.n_blocks();
734 AssertDimension(dof_handler.size(), n_blocks);
735
736 if (n_blocks == 0)
737 return;
738
739 const unsigned int min_level = src.min_level();
740 const unsigned int max_level = src.max_level();
741
742 for (unsigned int l = min_level; l <= max_level; ++l)
743 AssertDimension(src[l].n_blocks(), dst.n_blocks());
744
745 // FIXME: this a quite ugly as we need a temporary object:
747 min_level, max_level);
748
749 for (unsigned int b = 0; b < n_blocks; ++b)
750 {
751 for (unsigned int l = min_level; l <= max_level; ++l)
752 {
753 src_non_block[l].reinit(src[l].block(b));
754 src_non_block[l] = src[l].block(b);
755 }
756 const unsigned int data_block = same_for_all ? 0 : b;
757 get_matrix_free_transfer(data_block)
758 .copy_from_mg(*dof_handler[b], dst.block(b), src_non_block);
759 }
760}
761
762
763
764#endif // DOXYGEN
765
766
768
769#endif
const FiniteElement< dim, spacedim > & get_fe(const types::fe_index index=0) const
const IndexSet & locally_owned_mg_dofs(const unsigned int level) const
const Triangulation< dim, spacedim > & get_triangulation() const
types::global_dof_index n_dofs() const
unsigned int n_dofs_per_cell() const
virtual const FullMatrix< double > & get_restriction_matrix(const unsigned int child, const RefinementCase< dim > &refinement_case=RefinementCase< dim >::isotropic_refinement) const
bool restriction_is_additive(const unsigned int index) const
size_type n_elements() const
Definition index_set.h:1816
bool is_element(const size_type index) const
Definition index_set.h:1776
void copy_locally_owned_data_from(const Vector< Number2, MemorySpace > &src)
void reinit(const size_type size, const bool omit_zeroing_entries=false)
SmartPointer< const MGConstrainedDoFs, MGLevelGlobalTransfer< LinearAlgebra::distributed::Vector< Number > > > mg_constrained_dofs
std::function< void(const unsigned int, LinearAlgebra::distributed::Vector< Number > &)> initialize_dof_vector
virtual void prolongate(const unsigned int to_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 BlockVectorType2 &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
void copy_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::BlockVector< Number > > &dst, const BlockVectorType2 &src) const
virtual const TransferType & get_matrix_free_transfer(const unsigned int b) const =0
void copy_from_mg(const DoFHandler< dim, spacedim > &dof_handler, BlockVectorType2 &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number > > &src) const
virtual void restrict_and_add(const unsigned int from_level, LinearAlgebra::distributed::BlockVector< Number > &dst, const LinearAlgebra::distributed::BlockVector< Number > &src) const override
MGTransferBlockMatrixFreeBase(const bool same_for_all)
void copy_from_mg(const std::vector< const DoFHandler< dim, spacedim > * > &dof_handler, BlockVectorType2 &dst, const MGLevelObject< LinearAlgebra::distributed::BlockVector< Number > > &src) const
std::vector< MGTransferMatrixFree< dim, Number > > matrix_free_transfer_vector
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
void interpolate_to_mg(const DoFHandler< dim, spacedim > &dof_handler, MGLevelObject< LinearAlgebra::distributed::Vector< Number > > &dst, const BlockVectorType2 &src) const
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 > >())
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
void update_ghost_values() const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
unsigned int level
Definition grid_out.cc:4618
IteratorRange< cell_iterator > cell_iterators_on_level(const unsigned int level) const
#define DeclException0(Exception0)
Definition exceptions.h:465
static ::ExceptionBase & ExcNoProlongation()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
std::enable_if_t< IsBlockVector< VectorType >::value, void > collect_sizes(VectorType &vector)
Definition operators.h:96
std::enable_if_t< IsBlockVector< VectorType >::value, unsigned int > n_blocks(const VectorType &vector)
Definition operators.h:50
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)