operation_type_erasure.hpp

Go to the documentation of this file.

#pragma once
 
#include <deal.II/base/function.h>
 
#include <deal.II/dofs/dof_handler.h>
 
#include <deal.II/lac/la_parallel_vector.h>
 
#include <meltpooldg/compressible_flow/boundary_conditions.hpp>
#include <meltpooldg/compressible_flow/utils.hpp>
#include <meltpooldg/cut/cut_data.hpp>
 
#include <functional>
#include <memory>
#include <utility>
 
namespace MeltPoolDG::CompressibleFlow
{
  template <int dim, typename number>
  class OperationTypeErasure
  {
  public:
    using VectorType = dealii::LinearAlgebra::distributed::Vector<number>;
 
    OperationTypeErasure() = default;
 
    template <typename OperationType>
    explicit OperationTypeErasure(std::unique_ptr<OperationType> &&operation)
      : operation_pimpl(std::make_unique<OperationModel<OperationType>>(std::move(operation)))
    {}
 
    void
    solve(const number current_time, const number time_step)
    {
      operation_pimpl->solve(current_time, time_step);
    }
 
    void
    distribute_dofs(dealii::DoFHandler<dim> &dof_handler) const
    {
      operation_pimpl->distribute_dofs(dof_handler);
    }
 
    void
    reinit()
    {
      operation_pimpl->reinit();
    }
 
    std::vector<std::complex<number>>
    estimate_jacobian_eigenvalues(const number time_step) const
    {
      return operation_pimpl->estimate_jacobian_eigenvalues(time_step);
    }
 
    void
    attach_output_vectors(GenericDataOut<dim, number> &data_out) const
    {
      operation_pimpl->attach_output_vectors(data_out);
    }
 
    number
    compute_time_step_size(const bool do_print = false) const
    {
      return operation_pimpl->compute_time_step_size(do_print);
    }
 
    void
    set_initial_condition(const dealii::Function<dim> &function)
    {
      operation_pimpl->set_initial_condition(function);
    }
 
    void
    set_boundary_conditions(const std::shared_ptr<SimulationCaseBase<dim, number>> &simulation_case,
                            const std::string                                      &operation_name)
    {
      operation_pimpl->set_boundary_conditions(simulation_case, operation_name);
    }
 
    void
    set_body_force(std::unique_ptr<dealii::Function<dim>> body_force_in)
    {
      operation_pimpl->set_body_force(std::move(body_force_in));
    }
 
    bool
    repartition(VectorType &level_set, const std::function<void()> &distribute_level_set_dofs)
    {
      return operation_pimpl->repartition(level_set, distribute_level_set_dofs);
    }
 
    void
    add_external_force(
      std::shared_ptr<ExternalFlowForce<dim, number>>         external_force_residuum,
      std::shared_ptr<ExternalFlowForceJacobian<dim, number>> external_force_jacobian = nullptr);
 
    const dealii::DoFHandler<dim> &
    get_dof_handler() const
    {
      return operation_pimpl->get_dof_handler();
    }
 
    const VectorType &
    get_solution() const
    {
      return operation_pimpl->get_solution();
    }
 
    VectorType &
    get_solution()
    {
      return operation_pimpl->get_solution();
    }
 
    bool
    is_initialized() const
    {
      return (operation_pimpl != nullptr);
    }
 
    void
    clear()
    {
      operation_pimpl.reset(nullptr);
    }
 
  private:
    struct OperationConcept
    {
      virtual ~OperationConcept() = default;
 
      virtual void
      solve(number current_time, number time_step) = 0;
 
      virtual void
      distribute_dofs(dealii::DoFHandler<dim> &dof_handler) const = 0;
 
      virtual void
      attach_output_vectors(GenericDataOut<dim, number> &data_out) const = 0;
 
      virtual void
      reinit() = 0;
 
      virtual number
      compute_time_step_size(bool do_print = false) const = 0;
 
      virtual void
      set_initial_condition(const dealii::Function<dim> &function) = 0;
 
      virtual void
      set_boundary_conditions(
        const std::shared_ptr<SimulationCaseBase<dim, number>> &simulation_case,
        const std::string                                      &operation_name) = 0;
 
      virtual void
      set_body_force(std::unique_ptr<dealii::Function<dim>> body_force_in) = 0;
 
      virtual bool
      repartition(VectorType                  &level_set,
                  const std::function<void()> &distribute_level_set_dofs) = 0;
 
      virtual const VectorType &
      get_solution() const = 0;
 
      virtual VectorType &
      get_solution() = 0;
 
      virtual const dealii::DoFHandler<dim> &
      get_dof_handler() const = 0;
 
      virtual std::vector<std::complex<number>>
      estimate_jacobian_eigenvalues(const number time_step) const = 0;
    };
 
    template <typename OperationType>
    struct OperationModel final : public OperationConcept
    {
    public:
      explicit OperationModel(std::unique_ptr<OperationType> &&operation_in)
        : operation(std::move(operation_in))
      {}
 
      void
      solve(const number current_time, const number time_step) override
      {
        operation->solve(current_time, time_step);
      }
 
      void
      distribute_dofs(dealii::DoFHandler<dim> &dof_handler) const override
      {
        operation->distribute_dofs(dof_handler);
      }
 
      void
      reinit() override
      {
        operation->reinit();
      }
 
      void
      attach_output_vectors(GenericDataOut<dim, number> &data_out) const override
      {
        operation->attach_output_vectors(data_out);
      }
 
      number
      compute_time_step_size(const bool do_print = false) const override
      {
        return operation->compute_time_step_size(do_print);
      }
 
      void
      set_initial_condition(const dealii::Function<dim> &function) override
      {
        operation->set_initial_condition(function);
      }
 
      void
      set_boundary_conditions(
        const std::shared_ptr<SimulationCaseBase<dim, number>> &simulation_case,
        const std::string                                      &operation_name) override
      {
        operation->set_boundary_conditions(simulation_case, operation_name);
      }
 
      void
      set_body_force(std::unique_ptr<dealii::Function<dim>> body_force_in) override
      {
        operation->set_body_force(std::move(body_force_in));
      }
 
      bool
      repartition(VectorType                  &level_set,
                  const std::function<void()> &distribute_level_set_dofs) override
      {
        // repartitioning is only provided by operations implementing it (e.g. the compressible
        // multiphase operation); further, it does only make sense for cut applications
        if constexpr (requires { operation->repartition(level_set, distribute_level_set_dofs); })
          return operation->repartition(level_set, distribute_level_set_dofs);
        else
          {
            (void)level_set;
            (void)distribute_level_set_dofs;
            AssertThrow(false,
                        dealii::ExcMessage(
                          "Repartitioning is not supported by this compressible flow operation."));
            return false;
          }
      }
 
      const VectorType &
      get_solution() const override
      {
        return operation->get_solution();
      }
 
      VectorType &
      get_solution() override
      {
        return operation->get_solution();
      }
 
      const dealii::DoFHandler<dim> &
      get_dof_handler() const override
      {
        return operation->get_dof_handler();
      }
 
      std::vector<std::complex<number>>
      estimate_jacobian_eigenvalues(const number time_step) const override
      {
        return operation->estimate_jacobian_eigenvalues(time_step);
      }
 
    private:
      std::unique_ptr<OperationType> operation;
    };
 
    std::unique_ptr<OperationConcept> operation_pimpl;
  };
} // namespace MeltPoolDG::CompressibleFlow