metric_to_mesh_generator.h

#ifndef __METRIC_TO_MESH_GENERATOR_H__
#define __METRIC_TO_MESH_GENERATOR_H__

#include <deal.II/fe/fe_q.h>
#include <deal.II/fe/fe_system.h>
#include <deal.II/fe/mapping_fe_field.h>

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

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

#include <deal.II/lac/vector.h>
#include <deal.II/lac/la_parallel_vector.templates.h>

#include <deal.II/dofs/dof_handler.h>

namespace PHiLiP {
#if PHILIP_DIM==1 // dealii::parallel::distributed::Triangulation<dim> does not work for 1D
template <int dim, int nstate, typename real, typename MeshType = dealii::Triangulation<dim>>
#else
template <int dim, int nstate, typename real, typename MeshType = dealii::parallel::distributed::Triangulation<dim>>
#endif
class MetricToMeshGenerator {

    using VectorType = dealii::LinearAlgebra::distributed::Vector<double>; 
    using DoFHandlerType = dealii::DoFHandler<dim>; 

public:
    MetricToMeshGenerator(
        std::shared_ptr<dealii::MappingFEField<dim,dim,VectorType,DoFHandlerType>> _volume_nodes_mapping,
        std::shared_ptr<MeshType> _triangulation);

    ~MetricToMeshGenerator() = default;
    
    void generate_mesh_from_cellwise_metric(const std::vector<dealii::Tensor<2, dim, real>> &cellwise_optimal_metric);

    std::string get_generated_mesh_filename() const;

private:
    void reinit();
    
    void interpolate_metric_to_vertices(const std::vector<dealii::Tensor<2, dim, real>> &cellwise_optimal_metric);

    void write_geo_file() const;

    void write_pos_file() const;
    
    
    std::shared_ptr<dealii::MappingFEField<dim,dim,VectorType,DoFHandlerType>> volume_nodes_mapping;

    std::shared_ptr<MeshType> triangulation;
    
    dealii::DoFHandler<dim> dof_handler_vertices;

    const dealii::FE_Q<dim> fe_q;

    const dealii::FESystem<dim> fe_system;
    
    std::array<VectorType, dim*dim> optimal_metric_at_vertices;

    std::vector<dealii::Point<dim>> all_vertices;

    MPI_Comm mpi_communicator;
    
    dealii::ConditionalOStream pcout;

    int mpi_rank;

    int n_mpi;

    const std::string filename;
    const std::string filename_pos;
    const std::string filename_geo;
    const std::string filename_msh;
    
    dealii::IndexSet locally_owned_vertex_dofs; 
    dealii::IndexSet ghost_vertex_dofs; 
    dealii::IndexSet locally_relevant_vertex_dofs; 
    
}; // class ends

} // PHiLiP namespace
#endif