mesh_adaptor.h

Go to the documentation of this file.

//
// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
// Copyright (C) 2024 by the ryujin authors
//
 
#pragma once
 
#include "mpi_ensemble.h"
#include "offline_data.h"
 
#include <deal.II/base/parameter_acceptor.h>
#include <deal.II/base/smartpointer.h>
 
#include <random>
 
namespace ryujin
{
  enum class AdaptationStrategy {
    global_refinement,
 
    random_adaptation,
 
    kelly_estimator,
  };
 
  enum class MarkingStrategy {
    fixed_number,
    fixed_fraction,
  };
 
  enum class TimePointSelectionStrategy {
    fixed_time_points,
 
    simulation_cycle,
  };
} // namespace ryujin
 
#ifndef DOXYGEN
DECLARE_ENUM(
    ryujin::AdaptationStrategy,
    LIST({ryujin::AdaptationStrategy::global_refinement, "global refinement"},
         {ryujin::AdaptationStrategy::random_adaptation, "random adaptation"},
         {ryujin::AdaptationStrategy::kelly_estimator, "kelly estimator"}, ));
 
DECLARE_ENUM(ryujin::MarkingStrategy,
             LIST({ryujin::MarkingStrategy::fixed_number, "fixed number"},
                  {ryujin::MarkingStrategy::fixed_fraction,
                   "fixed fraction"}, ));
 
DECLARE_ENUM(ryujin::TimePointSelectionStrategy,
             LIST({ryujin::TimePointSelectionStrategy::fixed_time_points,
                   "fixed time points"},
                  {ryujin::TimePointSelectionStrategy::simulation_cycle,
                   "simulation cycle"}, ));
#endif
 
namespace ryujin
{
  template <typename Description, int dim, typename Number = double>
  class MeshAdaptor final : public dealii::ParameterAcceptor
  {
  public:
 
    using HyperbolicSystem = typename Description::HyperbolicSystem;
    using ParabolicSystem = typename Description::ParabolicSystem;
 
    using Triangulation = typename Discretization<dim>::Triangulation;
 
    using View =
        typename Description::template HyperbolicSystemView<dim, Number>;
 
    static constexpr auto problem_dimension = View::problem_dimension;
 
    using StateVector = typename View::StateVector;
    using InitialPrecomputedVector = typename View::InitialPrecomputedVector;
    using ScalarVector = Vectors::ScalarVector<Number>;
 
 
 
    MeshAdaptor(const MPIEnsemble &mpi_ensemble,
                const OfflineData<dim, Number> &offline_data,
                const HyperbolicSystem &hyperbolic_system,
                const ParabolicSystem &parabolic_system,
                const InitialPrecomputedVector &initial_precomputed,
                const ScalarVector &alpha,
                const std::string &subsection = "/MeshAdaptor");
 
    void prepare(const Number t);
 
    void analyze(const StateVector &state_vector,
                 const Number t,
                 unsigned int cycle);
 
    ACCESSOR_READ_ONLY(need_mesh_adaptation)
 
    
    void mark_cells_for_coarsening_and_refinement(
        Triangulation &triangulation) const;
 
  private:
 
    AdaptationStrategy adaptation_strategy_;
    std::uint_fast64_t random_adaptation_mersenne_twister_seed_;
 
    MarkingStrategy marking_strategy_;
    double refinement_fraction_;
    double coarsening_fraction_;
    unsigned int min_refinement_level_;
    unsigned int max_refinement_level_;
    unsigned int max_num_cells_;
 
    TimePointSelectionStrategy time_point_selection_strategy_;
    std::vector<Number> adaptation_time_points_;
    unsigned int adaptation_cycle_interval_;
 
    std::vector<std::string> kelly_quantities_;
 
 
 
    const MPIEnsemble &mpi_ensemble_;
 
    dealii::SmartPointer<const OfflineData<dim, Number>> offline_data_;
    dealii::SmartPointer<const HyperbolicSystem> hyperbolic_system_;
    dealii::SmartPointer<const ParabolicSystem> parabolic_system_;
 
    bool need_mesh_adaptation_;
 
    mutable dealii::Vector<float> indicators_;
 
    /* random adaptation: */
 
    void compute_random_indicators() const;
 
    mutable std::mt19937_64 mersenne_twister_;
 
    /* Kelly estimator: */
 
    void populate_kelly_quantities(const StateVector &state_vector) const;
    void compute_kelly_indicators() const;
 
    const InitialPrecomputedVector &initial_precomputed_;
    const ScalarVector &alpha_;
 
    mutable std::vector<ScalarVector> kelly_components_;
  };
 
} // namespace ryujin