tria.h

// ---------------------------------------------------------------------
//
// Copyright (C) 2008 - 2019 by the deal.II authors
//
// This file is part of the deal.II library.
//
// The deal.II library is free software; you can use it, redistribute
// it, and/or modify it under the terms of the GNU Lesser General
// Public License as published by the Free Software Foundation; either
// version 2.1 of the License, or (at your option) any later version.
// The full text of the license can be found in the file LICENSE.md at
// the top level directory of deal.II.
//
// ---------------------------------------------------------------------

#ifndef dealii_distributed_tria_h
#define dealii_distributed_tria_h


#include <deal.II/base/config.h>

#include <deal.II/base/smartpointer.h>
#include <deal.II/base/subscriptor.h>
#include <deal.II/base/template_constraints.h>

#include <deal.II/distributed/p4est_wrappers.h>
#include <deal.II/distributed/tria_base.h>

#include <deal.II/grid/tria.h>

#include <boost/range/iterator_range.hpp>

#include <functional>
#include <list>
#include <set>
#include <tuple>
#include <type_traits>
#include <utility>
#include <vector>

#ifdef DEAL_II_WITH_MPI
#  include <mpi.h>
#endif

#ifdef DEAL_II_WITH_P4EST
#  include <p4est.h>
#  include <p4est_connectivity.h>
#  include <p4est_ghost.h>
#  include <p8est.h>
#  include <p8est_connectivity.h>
#  include <p8est_ghost.h>
#endif


DEAL_II_NAMESPACE_OPEN

#ifdef DEAL_II_WITH_P4EST

namespace internal
{
  namespace DoFHandlerImplementation
  {
    namespace Policy
    {
      template <typename>
      class ParallelDistributed;
    }
  } // namespace DoFHandlerImplementation
} // namespace internal

namespace FETools
{
  namespace internal
  {
    template <int, int, class>
    class ExtrapolateImplementation;
  }
} // namespace FETools

// forward declaration of the data type for periodic face pairs
namespace GridTools
{
  template <typename CellIterator>
  struct PeriodicFacePair;
}

namespace parallel
{
  namespace distributed
  {
    template <int dim, int spacedim = dim>
    class Triangulation : public ::parallel::Triangulation<dim, spacedim>
    {
    public:
      using cell_iterator =
        typename ::Triangulation<dim, spacedim>::cell_iterator;

      using active_cell_iterator =
        typename ::Triangulation<dim, spacedim>::active_cell_iterator;

      using CellStatus =
        typename ::Triangulation<dim, spacedim>::CellStatus;

      enum Settings
      {
        default_setting = 0x0,
        mesh_reconstruction_after_repartitioning = 0x1,
        construct_multigrid_hierarchy = 0x2,
        no_automatic_repartitioning = 0x4
      };



      Triangulation(
        MPI_Comm mpi_communicator,
        const typename ::Triangulation<dim, spacedim>::MeshSmoothing
                       smooth_grid = (::Triangulation<dim, spacedim>::none),
        const Settings settings    = default_setting);

      virtual ~Triangulation() override;

      virtual void
      clear() override;

      virtual void
      copy_triangulation(
        const ::Triangulation<dim, spacedim> &other_tria) override;

      virtual void
      create_triangulation(const std::vector<Point<spacedim>> &vertices,
                           const std::vector<CellData<dim>> &  cells,
                           const SubCellData &subcelldata) override;

      virtual void
      execute_coarsening_and_refinement() override;

      virtual bool
      prepare_coarsening_and_refinement() override;

      void
      repartition();

      void
      communicate_locally_moved_vertices(
        const std::vector<bool> &vertex_locally_moved);


      virtual bool
      has_hanging_nodes() const override;

      virtual std::size_t
      memory_consumption() const override;

      virtual std::size_t
      memory_consumption_p4est() const;

      void
      write_mesh_vtk(const std::string &file_basename) const;

      unsigned int
      get_checksum() const;

      void
      save(const std::string &filename) const;

      void
      load(const std::string &filename, const bool autopartition = true);

      unsigned int
      register_data_attach(
        const std::function<std::vector<char>(const cell_iterator &,
                                              const CellStatus)> &pack_callback,
        const bool returns_variable_size_data);

      void
      notify_ready_to_unpack(
        const unsigned int handle,
        const std::function<void(
          const cell_iterator &,
          const CellStatus,
          const boost::iterator_range<std::vector<char>::const_iterator> &)>
          &unpack_callback);

      const std::vector<types::global_dof_index> &
      get_p4est_tree_to_coarse_cell_permutation() const;

      const std::vector<types::global_dof_index> &
      get_coarse_cell_to_p4est_tree_permutation() const;

      const typename ::internal::p4est::types<dim>::forest *
      get_p4est() const;

      virtual void
      add_periodicity(
        const std::vector<::GridTools::PeriodicFacePair<cell_iterator>> &)
        override;


    private:
      virtual void
      update_number_cache() override;

      Settings settings;

      bool triangulation_has_content;

      typename ::internal::p4est::types<dim>::connectivity *connectivity;

      typename ::internal::p4est::types<dim>::forest *parallel_forest;

      typename ::internal::p4est::types<dim>::ghost *parallel_ghost;

      struct CellAttachedData
      {
        unsigned int n_attached_data_sets;

        unsigned int n_attached_deserialize;

        using pack_callback_t = std::function<std::vector<char>(
          typename Triangulation<dim, spacedim>::cell_iterator,
          CellStatus)>;

        std::vector<pack_callback_t> pack_callbacks_fixed;
        std::vector<pack_callback_t> pack_callbacks_variable;
      };

      CellAttachedData cell_attached_data;

      using quadrant_cell_relation_t = typename std::tuple<
        typename ::internal::p4est::types<dim>::quadrant *,
        CellStatus,
        cell_iterator>;

      std::vector<quadrant_cell_relation_t> local_quadrant_cell_relations;

      void
      update_quadrant_cell_relations();

      class DataTransfer
      {
      public:
        DataTransfer(MPI_Comm mpi_communicator);

        void
        pack_data(
          const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
          const std::vector<typename CellAttachedData::pack_callback_t>
            &pack_callbacks_fixed,
          const std::vector<typename CellAttachedData::pack_callback_t>
            &pack_callbacks_variable);

        void
        execute_transfer(
          const typename ::internal::p4est::types<dim>::forest
            *parallel_forest,
          const typename ::internal::p4est::types<dim>::gloidx
            *previous_global_first_quadrant);

        void
        unpack_cell_status(
          std::vector<quadrant_cell_relation_t> &quad_cell_relations) const;

        void
        unpack_data(
          const std::vector<quadrant_cell_relation_t> &quad_cell_relations,
          const unsigned int                           handle,
          const std::function<void(
            const typename ::Triangulation<dim, spacedim>::cell_iterator
              &,
            const typename ::Triangulation<dim, spacedim>::CellStatus &,
            const boost::iterator_range<std::vector<char>::const_iterator> &)>
            &unpack_callback) const;

        void
        save(const typename ::internal::p4est::types<dim>::forest
               *                parallel_forest,
             const std::string &filename) const;

        void
        load(const typename ::internal::p4est::types<dim>::forest
               *                parallel_forest,
             const std::string &filename,
             const unsigned int n_attached_deserialize_fixed,
             const unsigned int n_attached_deserialize_variable);

        void
        clear();

      private:
        MPI_Comm mpi_communicator;

        bool variable_size_data_stored;

        std::vector<unsigned int> sizes_fixed_cumulative;

        std::vector<char> src_data_fixed;
        std::vector<char> dest_data_fixed;

        std::vector<int>  src_sizes_variable;
        std::vector<int>  dest_sizes_variable;
        std::vector<char> src_data_variable;
        std::vector<char> dest_data_variable;
      };

      DataTransfer data_transfer;

      std::vector<types::global_dof_index>
        coarse_cell_to_p4est_tree_permutation;
      std::vector<types::global_dof_index>
        p4est_tree_to_coarse_cell_permutation;

      typename ::internal::p4est::types<dim>::tree *
      init_tree(const int dealii_coarse_cell_index) const;

      void
      setup_coarse_cell_to_p4est_tree_permutation();

      void
      copy_new_triangulation_to_p4est(std::integral_constant<int, 2>);
      void
      copy_new_triangulation_to_p4est(std::integral_constant<int, 3>);

      void
      copy_local_forest_to_triangulation();

      std::vector<unsigned int>
      get_cell_weights() const;

      virtual std::map<unsigned int, std::set<::types::subdomain_id>>
      compute_vertices_with_ghost_neighbors() const override;

      std::vector<bool>
      mark_locally_active_vertices_on_level(const int level) const;

      template <typename>
      friend class ::internal::DoFHandlerImplementation::Policy::
        ParallelDistributed;

      template <int, int, class>
      friend class ::FETools::internal::ExtrapolateImplementation;
    };


    template <int spacedim>
    class Triangulation<1, spacedim>
      : public ::parallel::Triangulation<1, spacedim>
    {
    public:
      enum Settings
      {
        default_setting                          = 0x0,
        mesh_reconstruction_after_repartitioning = 0x1,
        construct_multigrid_hierarchy            = 0x2
      };

      Triangulation(
        MPI_Comm mpi_communicator,
        const typename ::Triangulation<1, spacedim>::MeshSmoothing
                       smooth_grid = (::Triangulation<1, spacedim>::none),
        const Settings settings    = default_setting);

      virtual ~Triangulation() override;

      const std::vector<types::global_dof_index> &
      get_p4est_tree_to_coarse_cell_permutation() const;

      void
      communicate_locally_moved_vertices(
        const std::vector<bool> &vertex_locally_moved);

      void
      load(const std::string &filename, const bool autopartition = true);

      void
      save(const std::string &filename) const;

      unsigned int
      register_data_attach(
        const std::function<std::vector<char>(
          const typename ::Triangulation<1, spacedim>::cell_iterator &,
          const typename ::Triangulation<1, spacedim>::CellStatus)>
          &        pack_callback,
        const bool returns_variable_size_data);

      void
      notify_ready_to_unpack(
        const unsigned int handle,
        const std::function<void(
          const typename ::Triangulation<1, spacedim>::cell_iterator &,
          const typename ::Triangulation<1, spacedim>::CellStatus,
          const boost::iterator_range<std::vector<char>::const_iterator> &)>
          &unpack_callback);

      std::vector<types::global_dof_index>
        coarse_cell_to_p4est_tree_permutation;
      std::vector<types::global_dof_index>
        p4est_tree_to_coarse_cell_permutation;


      // TODO: The following variable should really be private, but it is used
      // in dof_handler_policy.cc ...
      Settings settings;

      virtual std::map<unsigned int, std::set<::types::subdomain_id>>
      compute_vertices_with_ghost_neighbors() const override;

      virtual std::map<unsigned int, std::set<::types::subdomain_id>>
      compute_level_vertices_with_ghost_neighbors(
        const unsigned int level) const;

      virtual std::vector<bool>
      mark_locally_active_vertices_on_level(const unsigned int level) const;
    };
  } // namespace distributed
} // namespace parallel


#else // DEAL_II_WITH_P4EST

namespace parallel
{
  namespace distributed
  {
    template <int dim, int spacedim = dim>
    class Triangulation : public ::parallel::Triangulation<dim, spacedim>
    {
    public:
      Triangulation() = delete;
    };
  } // namespace distributed
} // namespace parallel


#endif


DEAL_II_NAMESPACE_CLOSE

#endif