FloodFiller.h

#ifndef INCLUDE_FLOODFILLER_H_
#define INCLUDE_FLOODFILLER_H_
 
#include <deal.II/base/quadrature.h>
#include <deal.II/dofs/dof_handler.h>
#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/fe_values.h>
#include <deal.II/grid/tria_accessor.h>
#include <deal.II/grid/tria_iterator.h>
#include <deal.II/lac/la_parallel_vector.h>
#include <deal.II/matrix_free/fe_evaluation.h>
 
template <int dim>
class GrainSet
{
public:
  void
  setGrainIndex(unsigned int _grain_index)
  {
    grain_index = _grain_index;
  };
 
  [[nodiscard]] unsigned int
  getGrainIndex() const
  {
    return grain_index;
  };
 
  void
  setOrderParameterIndex(unsigned int _order_parameter_index)
  {
    order_parameter_index = _order_parameter_index;
  };
 
  [[nodiscard]] unsigned int
  getOrderParameterIndex() const
  {
    return order_parameter_index;
  };
 
  void
  addVertexList(std::vector<dealii::Point<dim>> _vertices)
  {
    list_of_vertices.push_back(_vertices);
  };
 
  std::vector<std::vector<dealii::Point<dim>>>
  getVertexList() const
  {
    return list_of_vertices;
  };
 
private:
  unsigned int grain_index;
 
  unsigned int order_parameter_index;
 
  std::vector<std::vector<dealii::Point<dim>>> list_of_vertices;
};
 
template <int dim, int degree>
class FloodFiller
{
public:
  FloodFiller(dealii::FESystem<dim> &_fe, dealii::QGaussLobatto<dim> _quadrature)
    : quadrature(_quadrature)
    , num_quad_points(_quadrature.size())
    , dofs_per_cell(_fe.dofs_per_cell)
    , fe(&_fe) {};
 
  void
  calcGrainSets(dealii::FESystem<dim>                              &finite_element,
                dealii::DoFHandler<dim>                            &dof_handler,
                dealii::LinearAlgebra::distributed::Vector<double> *solution_field,
                double                                              threshold_lower,
                double                                              threshold_upper,
                int                                                 min_id,
                unsigned int                                        order_parameter_index,
                std::vector<GrainSet<dim>>                         &grain_sets);
 
protected:
  template <typename T>
  void
  recursiveFloodFill(T                                                   cell,
                     T                                                   cell_end,
                     dealii::LinearAlgebra::distributed::Vector<double> *solution_field,
                     double                                              threshold_lower,
                     double                                              threshold_upper,
                     int                                                 min_id,
                     unsigned int                                       &grain_index,
                     std::vector<GrainSet<dim>>                         &grain_sets,
                     bool                                               &grain_assigned);
 
  void
  createGlobalGrainSetList(std::vector<GrainSet<dim>> &grain_sets) const;
 
  void
  mergeSplitGrains(std::vector<GrainSet<dim>> &grain_sets) const;
 
  dealii::QGaussLobatto<dim> quadrature;
 
  const unsigned int num_quad_points;
 
  const unsigned int dofs_per_cell;
 
  dealii::FESystem<dim> *fe;
};
 
#endif