grid.h

/*  ______________________________________________________________________
 *
 *  ExaDG - High-Order Discontinuous Galerkin for the Exa-Scale
 *
 *  Copyright (C) 2021 by the ExaDG authors
 *
 *  This program is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  This program is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program. If not, see <https://www.gnu.org/licenses/>.
 *  ______________________________________________________________________
 */
 
#ifndef EXADG_GRID_GRID_H_
#define EXADG_GRID_GRID_H_
 
// deal.II
#include <deal.II/grid/grid_tools.h>
#include <deal.II/grid/tria.h>
 
// ExaDG
#include <exadg/grid/mapping_dof_vector.h>
 
namespace ExaDG
{
template<int dim>
class Grid
{
public:
  typedef typename std::vector<
    dealii::GridTools::PeriodicFacePair<typename dealii::Triangulation<dim>::cell_iterator>>
    PeriodicFacePairs;

  std::shared_ptr<dealii::Triangulation<dim>> triangulation;

  PeriodicFacePairs periodic_face_pairs;

  std::vector<std::shared_ptr<dealii::Triangulation<dim> const>> coarse_triangulations;

  std::vector<PeriodicFacePairs> coarse_periodic_face_pairs;
};

template<int dim, typename Number>
class MultigridMappings
{
public:
  MultigridMappings(std::shared_ptr<dealii::Mapping<dim>> mapping,
                    std::shared_ptr<dealii::Mapping<dim>> mapping_coarse_levels)
    : mapping_fine_level(mapping),
      mapping_coarse_levels(mapping_coarse_levels),
      degree_coarse_mappings(1)
  {
  }

  MultigridMappings(std::shared_ptr<MappingDoFVector<dim, Number>> mapping_dof_vector,
                    unsigned int const                             degree_coarse_mappings)
    : mapping_dof_vector_fine_level(mapping_dof_vector),
      degree_coarse_mappings(degree_coarse_mappings)
  {
  }

  void
  initialize_coarse_mappings(Grid<dim> const & grid, unsigned int const n_h_levels)
  {
    // we need to initialize mappings for coarse levels only if we have a mapping of type
    // MappingDoFVector
    if(mapping_dof_vector_fine_level.get())
    {
      if(n_h_levels > 1)
      {
        mapping_dof_vector_coarse_levels.resize(n_h_levels - 1);
 
        lambda_initialize_coarse_mappings(grid.triangulation, grid.coarse_triangulations);
      }
    }
    else // standard dealii::Mapping
    {
      AssertThrow(mapping_fine_level.get(),
                  dealii::ExcMessage("Fine-level mapping is uninitialized."));
 
      // when using standard dealii::Mapping's, there is nothing to do, i.e. we assume that all the
      // mappings have been set up prior to calling the constructor of this class.
    }
  }

  dealii::Mapping<dim> const &
  get_mapping(unsigned int const h_level, unsigned int const n_h_levels) const
  {
    if(mapping_dof_vector_fine_level.get()) // ExaDG::MappingDoFVector
    {
      // fine level
      if(h_level == n_h_levels - 1)
      {
        return *(mapping_dof_vector_fine_level->get_mapping());
      }
      else // coarse levels
      {
        AssertThrow(h_level < mapping_dof_vector_coarse_levels.size(),
                    dealii::ExcMessage("Vector of coarse mappings seems to have incorrect size."));
 
        return *(mapping_dof_vector_coarse_levels[h_level]->get_mapping());
      }
    }
    else // standard dealii::Mapping
    {
      // use mapping_fine_level on the fine level or on all levels if mapping_coarse_levels is
      // uninitialized
      if(h_level == n_h_levels - 1 or not(mapping_coarse_levels.get()))
      {
        AssertThrow(mapping_fine_level.get(),
                    dealii::ExcMessage("Fine-level mapping is uninitialized."));
 
        return *mapping_fine_level;
      }
      else // coarse levels
      {
        AssertThrow(mapping_coarse_levels.get(),
                    dealii::ExcMessage("mapping_coarse_levels is uninitialized."));
 
        return *mapping_coarse_levels;
      }
    }
  }

  std::function<void(
    std::shared_ptr<dealii::Triangulation<dim> const> const &              fine_triangulation,
    std::vector<std::shared_ptr<dealii::Triangulation<dim> const>> const & coarse_triangulations)>
    lambda_initialize_coarse_mappings =
      [&](std::shared_ptr<dealii::Triangulation<dim> const> const & fine_triangulation,
          std::vector<std::shared_ptr<dealii::Triangulation<dim> const>> const &
            coarse_triangulations) {
        AssertThrow(
          mapping_dof_vector_fine_level.get(),
          dealii::ExcMessage(
            "Coarse mappings can not be initialized because fine level mapping is invalid."));
 
        MappingTools::initialize_coarse_mappings<dim, Number>(mapping_dof_vector_coarse_levels,
                                                              degree_coarse_mappings,
                                                              mapping_dof_vector_fine_level,
                                                              fine_triangulation,
                                                              coarse_triangulations);
      };
 
private:
  std::shared_ptr<dealii::Mapping<dim>> mapping_fine_level;

  std::shared_ptr<dealii::Mapping<dim>> mapping_coarse_levels;

  std::shared_ptr<MappingDoFVector<dim, Number>> mapping_dof_vector_fine_level;

  std::vector<std::shared_ptr<MappingDoFVector<dim, Number>>> mapping_dof_vector_coarse_levels;

  unsigned int const degree_coarse_mappings;
};
 
} // namespace ExaDG
 
#endif /* EXADG_GRID_GRID_H_ */