Program Listing for File domain.h

↰ Return to documentation for file (include/immersx/core/domain.h)

// ---------------------------------------------------------------------
//
// Copyright (C) 2026 by Luca Heltai
//
// This file is part of the ImmersX application, based on the deal.II
// library.
//
// The ImmersX application is free software; you can use
// it, redistribute it, and/or modify it under the terms of the Apache-2.0
// License WITH LLVM-exception (the "License") as published by the Free
// Software Foundation; either version 3.0 of the License, or (at your option)
// any later version.
//
// ---------------------------------------------------------------------

#ifndef immersx_domain_h
#define immersx_domain_h

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

#include <deal.II/base/observer_pointer.h>
#include <deal.II/base/parameter_acceptor.h>
#include <deal.II/base/point.h>

#include <deal.II/distributed/fully_distributed_tria.h>
#include <deal.II/distributed/repartitioning_policy_tools.h>
#include <deal.II/distributed/tria.h>
#include <deal.II/distributed/tria_base.h>

#include <immersx/core/integration_support.h>
#include <immersx/coupling/tensor_product_lift.h>
#include <immersx/io/imported_finite_element_fields.h>
#include <mpi.h>

#include <functional>
#include <memory>
#include <string>
#include <vector>

namespace ImmersX
{
  template <int reduced_dim, int surface_dim, int spacedim, int n_components>
  class LiftedIntegrationSupport;

  inline std::string
  domain_subsection(const std::string &subsection)
  {
    if (subsection.empty())
      return "/Domain/";

    std::string normalized = subsection;
    if (normalized.front() != '/')
      normalized.insert(normalized.begin(), '/');
    if (normalized.back() != '/')
      normalized.push_back('/');
    return normalized + "Domain/";
  }

  template <int dim, int spacedim = dim>
  class DomainParameters : public dealii::ParameterAcceptor
  {
  public:
    explicit DomainParameters(const std::string &subsection = "/Domain/");

    std::string  name_of_grid       = "hyper_cube";
    std::string  arguments_for_grid = "-1: 1: false";
    double       grid_scale         = 1.;
    unsigned int initial_refinement = 0;
    std::string  triangulation_type = "distributed";

    std::string mesh_input = "auto";

    double surface_mesh_size = 1.;

    double surface_feature_angle = 40.;

    std::string surface_volume_cells = "tetrahedra";

    bool read_vtk_fields = false;

    std::string vtk_boundary_id_field;

    std::string vtk_material_id_field;
  };


  template <int dim, int spacedim = dim>
  class Domain
  {
  public:
    using ImportedFields = ImportedFiniteElementFields<dim, spacedim>;

    using DistributedTriangulation =
      dealii::parallel::distributed::Triangulation<dim, spacedim>;
    using FullyDistributedTriangulation =
      dealii::parallel::fullydistributed::Triangulation<dim, spacedim>;
    using TriangulationType =
      dealii::parallel::TriangulationBase<dim, spacedim>;

    explicit Domain(const DomainParameters<dim, spacedim> &parameters,
                    MPI_Comm communicator = MPI_COMM_WORLD);

    Domain(const DomainParameters<dim, spacedim> &parameters,
           TriangulationType                     &triangulation);

    ~Domain();

    Domain(const Domain &) = delete;
    Domain &
    operator=(const Domain &) = delete;

    void
    make_grid();

    void
    prepare_grid(bool serial = false);

    void
    prepare_grid(const dealii::Triangulation<dim, spacedim> &coarse_mesh);

    void
    finalize_grid();

    void
    finalize_grid(
      const dealii::RepartitioningPolicyTools::Base<dim, spacedim> &policy);

    bool
    grid_is_prepared() const;

    dealii::Triangulation<dim, spacedim> &
    refinement_triangulation();

    TriangulationType &
    triangulation();

    const TriangulationType &
    triangulation() const;

    bool
    uses_fully_distributed_triangulation() const;

    MPI_Comm
    communicator() const;

    IntegrationSupport<dim, spacedim>
    support() const
    {
      return IntegrationSupport<dim, spacedim>(triangulation());
    }

    IntegrationSupport<dim, spacedim>
    support(const dealii::types::material_id material_id) const
    {
      return IntegrationSupport<dim, spacedim>(triangulation(), material_id);
    }

    template <int surface_dim, int n_components>
    LiftedIntegrationSupport<dim, surface_dim, spacedim, n_components>
    lift(const TensorProductLift<dim, surface_dim, spacedim, n_components>
           &descriptor) const;

    template <int surface_dim, int n_components>
    LiftedIntegrationSupport<dim, surface_dim, spacedim, n_components>
    lift(const TensorProductLift<dim, surface_dim, spacedim, n_components>
                                          &descriptor,
         FixedThicknessEvaluator<spacedim> thickness_evaluator) const;

    BoundaryIntegrationSupport<dim, spacedim>
    boundary_support() const
    {
      return BoundaryIntegrationSupport<dim, spacedim>(triangulation());
    }

    BoundaryIntegrationSupport<dim, spacedim>
    boundary_support(const dealii::types::boundary_id boundary_id) const
    {
      return BoundaryIntegrationSupport<dim, spacedim>(triangulation(),
                                                       boundary_id);
    }

    bool
    has_vtk_fields() const;

    const ImportedFields &
    vtk_fields() const;

    typename ImportedFields::FieldView
    field(const std::string &name, unsigned int component = 0) const;

  private:
    void
    finalize_grid_impl(
      const dealii::RepartitioningPolicyTools::Base<dim, spacedim> *policy);

    static std::unique_ptr<TriangulationType>
    make_triangulation_storage(
      const DomainParameters<dim, spacedim> &parameters,
      MPI_Comm                               communicator);

    void
    prepare_grid_impl(bool                                        serial,
                      const dealii::Triangulation<dim, spacedim> *coarse_mesh);

    void
    import_vtk_fields_if_requested();

    struct SerialFields;

    const DomainParameters<dim, spacedim>     &parameters_;
    MPI_Comm                                   communicator_;
    std::unique_ptr<TriangulationType>         triangulation_storage_;
    dealii::ObserverPointer<TriangulationType> triangulation_;
    std::unique_ptr<ImportedFields>            vtk_fields_;
    std::unique_ptr<dealii::Triangulation<dim, spacedim>> serial_triangulation_;
    std::unique_ptr<SerialFields>                         serial_fields_;
    bool grid_prepared_ = false;
  };
} // namespace ImmersX

#include <immersx/core/lifted_integration_support.h>

#endif // immersx_domain_h