CD_ItoSolver.H

Go to the documentation of this file.

/*
 * SPDX-FileCopyrightText: 2021-2026 SINTEF Energy Research
 *
 * SPDX-License-Identifier: GPL-3.0-or-later
 */
 
#ifndef CD_ITOSOLVER_H
#define CD_ITOSOLVER_H
 
// Std includes
#include <optional>
 
// Our includes
#include <CD_AmrMesh.H>
#include <CD_ComputationalGeometry.H>
#include <CD_ItoSpecies.H>
#include <CD_ItoParticle.H>
#include <CD_ItoMergeParticle.H>
#include <CD_ParticleContainer.H>
#include <CD_EBRepresentation.H>
#include <CD_EBIntersection.H>
#include <CD_CellInfo.H>
#include <CD_ParticleManagement.H>
#include <CD_NamespaceHeader.H>
 
class ItoSolver
{
public:
  enum class WhichContainer
  {
    Bulk,
    EB,
    Domain,
    Source,
    Covered,
    Scratch
  };
 
  ItoSolver();
 
  virtual ~ItoSolver();
 
  ItoSolver(const ItoSolver&) = delete;
 
  ItoSolver&
  operator=(const ItoSolver&) = delete;
 
  ItoSolver(ItoSolver&&) = default;
 
  ItoSolver&
  operator=(ItoSolver&&) = default;
 
  virtual void
  setParticleCellMerger(const ParticleManagement::ParticleMerger<ItoMergeParticle>& a_particleCellMerger) noexcept;
 
  virtual std::string
  getName() const;
 
  virtual std::string
  getRealm() const;
 
  virtual void
  setRealm(const std::string& a_realm);
 
  virtual void
  parseOptions();
 
  virtual void
  parseRuntimeOptions();
 
  virtual EBIntersection
  getIntersectionAlgorithm() const noexcept;
 
  virtual void
  initialData();
 
  virtual void
  generateParticlesFromDensity(ParticleContainer<ItoParticle>&              a_particles,
                               const std::function<Real(const RealVect x)>& a_densityFunc,
                               const Vector<int>&                           a_maxParticlesPerCell) const noexcept;
 
  virtual void
  generateParticlesFromDensity(ParticleContainer<ItoParticle>&              a_particles,
                               const std::function<Real(const RealVect x)>& a_densityFunc,
                               int                                          a_maxParticlesPerCell) const noexcept;
 
  virtual void
  regrid(int a_lmin, int a_oldFinestLevel, int a_newFinestLevel);
 
  virtual void
  allocate();
 
  virtual void
  clear(WhichContainer a_container);
 
  virtual void
  clear(ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  depositConductivity(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  depositConductivity(EBAMRCellData&                  a_phi,
                      ParticleContainer<ItoParticle>& a_particles,
                      DepositionType                  a_deposition,
                      CoarseFineDeposition            a_coarseFineDeposition) const;
 
  virtual void
  depositDiffusivity(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  depositDiffusivity(EBAMRCellData&                  a_phi,
                     ParticleContainer<ItoParticle>& a_particles,
                     DepositionType                  a_deposition,
                     CoarseFineDeposition            a_coarseFineDeposition) const;
 
  virtual void
  depositEnergyDensity(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  depositEnergyDensity(EBAMRCellData&                  a_phi,
                       ParticleContainer<ItoParticle>& a_particles,
                       DepositionType                  a_deposition,
                       CoarseFineDeposition            a_coarseFineDeposition) const;
 
  virtual void
  computeAverageMobility(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  computeAverageDiffusion(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  computeAverageEnergy(EBAMRCellData& a_phi, ParticleContainer<ItoParticle>& a_particles) const;
 
  virtual void
  depositParticles();
 
  virtual void
  depositParticles(WhichContainer a_container);
 
  void
  depositWeightNGP(LevelData<EBCellFAB>&                 a_output,
                   const ParticleContainer<ItoParticle>& a_particles,
                   int                                   a_level) const noexcept;
 
  template <typename Gather>
  void
  depositGatheredNGP(LevelData<EBCellFAB>&                 a_output,
                     const ParticleContainer<ItoParticle>& a_particles,
                     int                                   a_level,
                     Gather                                a_gather) const noexcept;
 
  void
  coarsenAndFillGhosts(EBAMRCellData& a_phi) const;
 
  template <typename P, typename Traits>
  void
  depositWeight(EBAMRCellData&                      a_phi,
                const ParticleContainer<P, Traits>& a_particles,
                DepositionType                      a_deposition,
                CoarseFineDeposition                a_coarseFineDeposition) const;
 
  template <typename Gather>
  void
  depositGathered(EBAMRCellData&                        a_phi,
                  const ParticleContainer<ItoParticle>& a_particles,
                  DepositionType                        a_deposition,
                  CoarseFineDeposition                  a_coarseFineDeposition,
                  Gather                                a_gather) const;
 
  virtual void
  removeCoveredParticles(EBRepresentation a_representation, Real a_tol);
 
  virtual void
  removeCoveredParticles(WhichContainer a_container, EBRepresentation a_representation, Real a_tol);
 
  virtual void
  removeCoveredParticles(ParticleContainer<ItoParticle>& a_particles,
                         EBRepresentation                a_representation,
                         Real                            a_tol) const;
 
  virtual void
  transferCoveredParticles(EBRepresentation a_representation, Real a_tol);
 
  virtual void
  transferCoveredParticles(WhichContainer   a_containerFrom,
                           WhichContainer   a_containerTo,
                           EBRepresentation a_representation,
                           Real             a_tol);
 
  virtual void
  transferCoveredParticles(ParticleContainer<ItoParticle>& a_particlesFrom,
                           ParticleContainer<ItoParticle>& a_particlesTo,
                           EBRepresentation                a_representation,
                           Real                            a_tol) const;
 
  virtual void
  intersectParticles(
    const EBIntersection                                               a_ebIntersection,
    const bool                                                         a_deleteParticles,
    const std::function<void(ParticleSoA<ItoParticle>&, std::size_t)>& a_nonDeletionModifier =
      [](ParticleSoA<ItoParticle>&, std::size_t) -> void {
      return;
    });
 
  virtual void
  intersectParticles(
    const WhichContainer                                               a_particles,
    const WhichContainer                                               a_ebParticles,
    const WhichContainer                                               a_domainParticles,
    const EBIntersection                                               a_ebIntersection,
    const bool                                                         a_deleteParticles,
    const std::function<void(ParticleSoA<ItoParticle>&, std::size_t)>& a_nonDeletionModifier =
      [](ParticleSoA<ItoParticle>&, std::size_t) -> void {
      return;
    });
 
  virtual void
  intersectParticles(
    ParticleContainer<ItoParticle>&                                    a_particles,
    ParticleContainer<ItoParticle>&                                    a_ebParticles,
    ParticleContainer<ItoParticle>&                                    a_domainParticles,
    const EBIntersection                                               a_ebIntersection,
    const bool                                                         a_deleteParticles,
    const std::function<void(ParticleSoA<ItoParticle>&, std::size_t)>& a_nonDeletionModifier =
      [](ParticleSoA<ItoParticle>&, std::size_t) -> void {
      return;
    });
 
  virtual void
  computeLoads(Vector<long int>& a_loads, const DisjointBoxLayout& a_dbl, int a_level);
 
  virtual int
  getNumberOfPlotVariables() const;
 
  virtual Vector<std::string>
  getPlotVariableNames() const;
 
  virtual void
  writePlotData(LevelData<EBCellFAB>& a_output,
                int&                  a_comp,
                const std::string&    a_outputRealm,
                int                   a_level) const noexcept;
 
  virtual unsigned long long
  getNumParticles(WhichContainer a_whichContainer, bool a_localOnly) const;
 
  virtual void
  registerOperators() const;
 
  virtual void
  setComputationalGeometry(const RefCountedPtr<ComputationalGeometry>& a_computationalGeometry);
 
  virtual void
  setAmr(const RefCountedPtr<AmrMesh>& a_amr);
 
  virtual void
  setPhase(phase::which_phase a_phase);
 
  virtual void
  setVerbosity(int a_verbosity);
 
  virtual void
  setTime(int a_step, Real a_time, Real a_dt);
 
#ifdef CH_USE_HDF5
  virtual void
  writeCheckpointLevel(HDF5Handle& a_handle, int a_level) const;
#endif
 
#ifdef CH_USE_HDF5
  virtual void
  readCheckpointLevel(HDF5Handle& a_handle, int a_level);
#endif
 
  virtual void
  setSpecies(const RefCountedPtr<ItoSpecies>& a_species);
 
  virtual void
  preRegrid(int a_lbase, int a_oldFinestLevel);
 
  virtual bool
  isMobile() const;
 
  virtual bool
  isDiffusive() const;
 
  virtual ParticleContainer<ItoParticle>&
  getParticles(WhichContainer a_container);
 
  virtual const ParticleContainer<ItoParticle>&
  getParticles(WhichContainer a_container) const;
 
  virtual const RefCountedPtr<ItoSpecies>&
  getSpecies() const;
 
  virtual EBAMRCellData&
  getPhi();
 
  virtual EBAMRCellData&
  getVelocityFunction();
 
  virtual const EBAMRCellData&
  getVelocityFunction() const;
 
  virtual EBAMRCellData&
  getDiffusionFunction();
 
  virtual const EBAMRCellData&
  getDiffusionFunction() const;
 
  virtual EBAMRCellData&
  getMobilityFunction();
 
  virtual const EBAMRCellData&
  getMobilityFunction() const;
 
  virtual void
  setDiffusionFunction(Real a_diffusionCoefficient);
 
  virtual void
  setVelocityFunction(const RealVect& a_velocity);
 
  virtual void
  setParticleMobility(Real a_mobility);
 
  virtual void
  setParticleDiffusion(Real a_diffCo);
 
  virtual void
  interpolateVelocities();
 
  virtual void
  interpolateVelocities(int a_lvl, const DataIndex& a_dit);
 
  virtual void
  interpolateMobilities();
 
  virtual void
  updateMobilities();
 
  virtual void
  interpolateDiffusion();
 
  virtual void
  updateDiffusion();
 
  virtual void
  organizeParticlesByCell(WhichContainer a_container);
 
  virtual void
  organizeParticlesByPatch(WhichContainer a_container);
 
  virtual void
  makeSuperparticles(WhichContainer a_container);
 
  virtual void
  makeSuperparticles(WhichContainer a_container, int a_particlesPerCell);
 
  virtual void
  makeSuperparticles(WhichContainer a_container, const Vector<int>& a_particlesPerCell);
 
  virtual void
  makeSuperparticles(WhichContainer                                a_container,
                     const Vector<int>&                            a_particlesPerCell,
                     const ParticleManagement::ParticleMergeMethod a_method);
 
  virtual const Vector<int>&
  getParticlesPerCell() const noexcept;
 
  virtual ParticleContainer<ItoMergeParticle>&
  getRegridParticles() noexcept;
 
  virtual void
  remap();
 
  virtual void
  remap(WhichContainer a_container);
 
  virtual void
  remapAll();
 
  inline RealVect
  randomGaussian() const;
 
  inline RealVect
  randomDirection() const;
 
  virtual Real
  computeDt() const;
 
  virtual Real
  computeDt(int a_lvl) const;
 
  virtual Real
  computeDt(int a_lvl, const DataIndex& a_dit) const;
 
  virtual Real
  computeHopDt(Real a_maxCellsToMove) const;
 
  virtual Real
  computeHopDt(Real a_maxCellsToMove, int a_lvl) const;
 
  virtual Real
  computeHopDt(Real a_maxCellsToMove, int a_lvl, const DataIndex& a_dit) const;
 
  virtual Real
  computeAdvectiveDt() const;
 
  virtual Real
  computeAdvectiveDt(int a_lvl) const;
 
  virtual Real
  computeAdvectiveDt(int a_lvl, const DataIndex& a_dit) const;
 
  virtual Real
  computeDiffusiveDt() const;
 
  virtual Real
  computeDiffusiveDt(int a_lvl) const;
 
  virtual Real
  computeDiffusiveDt(int a_lvl, const DataIndex& a_dit) const;
 
  DepositionType
  getDeposition() const;
 
  CoarseFineDeposition
  getCoarseFineDeposition() const;
 
  phase::which_phase
  getPhase() const;
 
protected:
  static constexpr int m_comp = 0;
 
  static constexpr int m_nComp = 1;
 
  enum class WhichCheckpoint
  {
    Particles,
    Numbers
  };
 
  enum class WhichMobilityInterpolation
  {
    Direct,
    Velocity
  };
 
  WhichCheckpoint m_checkpointing;
 
  WhichMobilityInterpolation m_mobilityInterp;
 
  ParticleManagement::ParticleMerger<ItoMergeParticle> m_particleCellMerger;
 
  ParticleManagement::ParticleMergeMethod m_mergeMethod = ParticleManagement::ParticleMergeMethod::EqualWeightKD;
 
  Vector<int> m_particlesPerCell;
 
  std::optional<ParticleManagement::ParticleMergeMethod> m_regridMergeMethod;
 
  bool m_nnPairIterate = true;
 
  int m_nnPairFallback = 1;
 
  int m_nnPairMaxCellDistance = -1;
 
  int m_nnPairMaxRounds = 3;
 
  Real m_kdSplitWeightLeafDx = 1.0;
 
  int m_restartPPC;
 
  std::string m_realm;
 
  RefCountedPtr<ComputationalGeometry> m_computationalGeometry;
 
  RefCountedPtr<AmrMesh> m_amr;
 
  RefCountedPtr<ItoSpecies> m_species;
 
  phase::which_phase m_phase;
 
  std::string m_name;
 
  std::string m_className;
 
  Real m_normalDistributionTruncation;
 
  Real m_bisectionStep;
 
  int m_verbosity;
 
  int m_timeStep;
 
  Real m_time;
 
  Real m_dt;
 
  bool m_forceIrregDepositionNGP;
 
  bool m_forceIrregInterpolationNGP;
 
  bool m_forceHaloNGP;
 
  bool m_useRedistribution;
 
  bool m_blendConservation;
 
  bool m_isDiffusive;
 
  bool m_isMobile;
 
  bool m_plotPhi;
 
  bool m_plotVelocity;
 
  bool m_plotDiffCo;
 
  bool m_plotParticles;
 
  bool m_plotParticlesEB;
 
  bool m_plotParticlesDomain;
 
  bool m_plotParticlesSource;
 
  bool m_plotParticlesCovered;
 
  bool m_plotEnergyDensity;
 
  bool m_plotAverageEnergy;
 
  DepositionType m_deposition;
 
  CoarseFineDeposition m_coarseFineDeposition;
 
  DepositionType m_plotDeposition;
 
  EBAMRCellData m_phi;
 
  EBAMRCellData m_mobilityFunction;
 
  EBAMRCellData m_velocityFunction;
 
  EBAMRCellData m_diffusionFunction;
 
  EBAMRFAB m_kdMergeCellHistogram;
 
  EBAMRFAB m_kdMergeLeafQuota;
 
  mutable EBAMRIVData m_depositionNC;
 
  mutable EBAMRIVData m_massDiff;
 
  std::map<WhichContainer, ParticleContainer<ItoParticle>> m_particleContainers;
 
  ParticleContainer<ItoMergeParticle> m_regridParticles;
 
  EBIntersection m_intersectionAlg;
 
  void
  parseRNG();
 
  void
  parseTruncation();
 
  void
  parseDeposition();
 
  void
  parseParticleMerger();
 
  void
  mergeLite(ParticleContainer<ItoMergeParticle>&    a_particles,
            const Vector<int>&                      a_particlesPerCell,
            ParticleManagement::ParticleMergeMethod a_method);
 
  void
  applyCellMerger(ParticleContainer<ItoMergeParticle>&                        a_particles,
                  const Vector<int>&                                          a_particlesPerCell,
                  const ParticleManagement::ParticleMerger<ItoMergeParticle>& a_merger);
 
  void
  invalidateParticleIDs(ParticleContainer<ItoMergeParticle>& a_particles) const;
 
  void
  extractIntoMergeContainer(WhichContainer a_container, ParticleContainer<ItoMergeParticle>& a_merge);
 
  void
  splitFromMergeContainer(ParticleContainer<ItoMergeParticle>& a_merge, int a_numParticlesPerCellThresh);
 
  void
  rebuildFromMergeContainer(WhichContainer a_container, ParticleContainer<ItoMergeParticle>& a_merge);
 
  enum class NnPairSearchBackend
  {
    Tree, 
    Hash  
  };
 
  void
  mergeEqualWeightKD(ParticleContainer<ItoMergeParticle>& a_particles, const Vector<int>& a_particlesPerCell);
 
  void
  mergeReinitialize(ParticleContainer<ItoMergeParticle>& a_particles, const Vector<int>& a_particlesPerCell);
 
  void
  mergeReinitializeBVH(ParticleContainer<ItoMergeParticle>& a_particles, const Vector<int>& a_particlesPerCell);
 
  void
  mergeNnSfc(ParticleContainer<ItoMergeParticle>& a_particles, const Vector<int>& a_particlesPerCell);
 
  template <NnPairSearchBackend Backend>
  void
  mergeNnPairSearch(ParticleContainer<ItoMergeParticle>& merge, const Vector<int>& a_particlesPerCell);
 
  void
  mergeNnPairOneCell(ParticleContainer<ItoMergeParticle>& merge, const Vector<int>& a_particlesPerCell);
 
  void
  mergeKDSkinNn(ParticleContainer<ItoMergeParticle>& merge, const Vector<int>& a_particlesPerCell);
 
  void
  mergeKDImpl(ParticleContainer<ItoMergeParticle>& merge,
              const Vector<int>&                   a_particlesPerCell,
              bool                                 a_enableBoundaryCarve);
 
  void
  parseIntersectionEB();
 
  void
  parseVerbosity();
 
  void
  parsePlotVariables();
 
  void
  parseDiffusionHop();
 
  void
  parseRedistribution();
 
  void
  parseDivergenceComputation();
 
  void
  parseCheckpointing();
 
#ifdef CH_USE_HDF5
  virtual void
  writeCheckPointLevelParticles(HDF5Handle& a_handle, int a_level) const;
#endif
 
#ifdef CH_USE_HDF5
  virtual void
  writeCheckPointLevelFluid(HDF5Handle& a_handle, int a_level) const;
#endif
 
#ifdef CH_USE_HDF5
  virtual void
  readCheckpointLevelParticles(HDF5Handle& a_handle, int a_level);
#endif
 
#ifdef CH_USE_HDF5
  virtual void
  readCheckpointLevelFluid(HDF5Handle& a_handle, int a_level);
#endif
 
  void
  drawNewParticles(const LevelData<EBCellFAB>& a_particlesPerCell, int a_level, int a_newPPC);
 
  virtual void
  redistributeAMR(EBAMRCellData& a_phi) const;
 
  virtual void
  depositNonConservative(EBAMRIVData& a_depositionNC, const EBAMRCellData& a_depositionKappaC) const;
 
  virtual void
  depositHybrid(EBAMRCellData& a_depositionH, EBAMRIVData& a_massDifference, const EBAMRIVData& a_depositionNC) const;
 
  virtual void
  interpolateMobilities(int a_lvl, const DataIndex& a_dit, const EBCellFAB& a_velocityMagnitude) noexcept;
 
  virtual void
  interpolateMobilitiesDirect(int a_lvl, const DataIndex& a_dit) noexcept;
 
  virtual void
  interpolateMobilitiesVelocity(int a_lvl, const DataIndex& a_dit, const EBCellFAB& a_velocityMagnitude) noexcept;
 
  virtual void
  updateMobilities(int a_level, const DataIndex& a_dit);
 
  virtual void
  interpolateDiffusion(int a_lvl, const DataIndex& a_dit);
 
  virtual void
  updateDiffusion(int a_level, const DataIndex& a_dit);
 
  virtual void
  writeData(LevelData<EBCellFAB>& a_output,
            int&                  a_comp,
            const EBAMRCellData&  a_data,
            const std::string&    a_outputRealm,
            int                   a_level,
            bool                  a_interpToCentroids,
            bool                  a_interpGhost) const noexcept;
};
 
#include <CD_NamespaceFooter.H>
 
#include <CD_ItoSolverImplem.H>
 
#endif