LocalDtAcousticSolver.cpp

#include <scheme/LocalDtAcousticSolver.hpp>

#include <language/utils/InterpolateItemValue.hpp>
#include <mesh/ItemValueUtils.hpp>
#include <mesh/ItemValueVariant.hpp>
#include <mesh/MeshFaceBoundary.hpp>
#include <mesh/MeshFlatNodeBoundary.hpp>
#include <mesh/MeshNodeBoundary.hpp>
#include <mesh/SubItemValuePerItemVariant.hpp>
#include <scheme/DirichletBoundaryConditionDescriptor.hpp>
#include <scheme/DiscreteFunctionP0.hpp>
#include <scheme/DiscreteFunctionUtils.hpp>
#include <scheme/ExternalBoundaryConditionDescriptor.hpp>
#include <scheme/FixedBoundaryConditionDescriptor.hpp>
#include <scheme/IBoundaryConditionDescriptor.hpp>
#include <scheme/IDiscreteFunctionDescriptor.hpp>
#include <scheme/SymmetryBoundaryConditionDescriptor.hpp>
#include <utils/Socket.hpp>

#include <variant>
#include <vector>

template <size_t Dimension>
class LocalDtAcousticSolverHandler::LocalDtAcousticSolver final
  : public LocalDtAcousticSolverHandler::ILocalDtAcousticSolver
{
 private:
  using Rdxd = TinyMatrix<Dimension>;
  using Rd   = TinyVector<Dimension>;

  using MeshType     = Mesh<Connectivity<Dimension>>;
  using MeshDataType = MeshData<Dimension>;

  using DiscreteScalarFunction = DiscreteFunctionP0<Dimension, const double>;
  using DiscreteVectorFunction = DiscreteFunctionP0<Dimension, const Rd>;

  class FixedBoundaryCondition;
  class PressureBoundaryCondition;
  class SymmetryBoundaryCondition;
  class VelocityBoundaryCondition;
  class ExternalFSIVelocityBoundaryCondition;

  using BoundaryCondition = std::variant<FixedBoundaryCondition,
                                         PressureBoundaryCondition,
                                         SymmetryBoundaryCondition,
                                         VelocityBoundaryCondition,
                                         ExternalFSIVelocityBoundaryCondition>;

  using BoundaryConditionList = std::vector<BoundaryCondition>;

  NodeValuePerCell<const Rdxd>
  _computeGlaceAjr(const MeshType& mesh, const DiscreteScalarFunction& rhoc) const
  {
    MeshDataType& mesh_data = MeshDataManager::instance().getMeshData(mesh);

    const NodeValuePerCell<const Rd> Cjr = mesh_data.Cjr();
    const NodeValuePerCell<const Rd> njr = mesh_data.njr();

    NodeValuePerCell<Rdxd> Ajr{mesh.connectivity()};

    parallel_for(
      mesh.numberOfCells(), PUGS_LAMBDA(CellId j) {
        const size_t& nb_nodes = Ajr.numberOfSubValues(j);
        const double& rhoc_j   = rhoc[j];
        for (size_t r = 0; r < nb_nodes; ++r) {
          Ajr(j, r) = tensorProduct(rhoc_j * Cjr(j, r), njr(j, r));
        }
      });

    return Ajr;