BaseParticleContainer.hxx

Go to the documentation of this file.

 
#ifndef BASEPARTICLESCONTAINER_HXX
#define BASEPARTICLESCONTAINER_HXX
 
// Include libraries
#include "AMReX_CTOParallelForImpl.H"
#include <AMReX_AmrParticles.H>
#include <AMReX_Particles.H>
#include <cctk.h>
#include <cctk_Arguments.h>
#include <cctk_Parameters.h>
#include <cctk_core.h>
#include <string>
 
namespace GInX {
 
template <typename StructType>
class ParticleIterator
    : public amrex::ParIter<0, 0, StructType::n_attributes, 0> {
public:
  using amrex::ParIter<0, 0, StructType::n_attributes, 0>::ParIter;
  using RealVector = typename amrex::ParIter<
      0, 0, StructType::n_attributes>::ContainerType::RealVector;
 
  const std::array<RealVector, StructType::n_attributes> &GetAttribs() const {
    return this->GetStructOfArrays().GetRealData();
  }
 
  std::array<RealVector, StructType::n_attributes> &GetAttributes() {
    return this->GetStructOfArrays().GetRealData();
  }
 
  const RealVector &GetAttribs(int comp) const {
    return this->GetStructOfArrays().GetRealData(comp);
  }
 
  RealVector &GetAttributes(int comp) {
    return this->GetStructOfArrays().GetRealData(comp);
  }
}; // class ParicleIterator
 
template <typename OtherContainer, typename StructType>
class BaseParticleContainer
    : public amrex::AmrParticleContainer<0, 0, StructType::n_attributes, 0> {
public:
  // Name of the particles
  const std::string name = StructType::name;
  // Number of attributes per each particle
  static constexpr int n_attributes = StructType::n_attributes;
 
  BaseParticleContainer(amrex::AmrCore *amr_core)
      : amrex::AmrParticleContainer<0, 0, StructType::n_attributes, 0>(
            amr_core) {}
 
  virtual ~BaseParticleContainer() = default;
 
  template <typename Function>
  void initialize(Function initializer_function, const CCTK_REAL *real_params,
                  const CCTK_INT *int_params) {
    initializer_function(static_cast<OtherContainer &>(*this), real_params,
                         int_params);
  };
 
  template <typename Function>
  void initialize(Function initializer_function, const amrex::MultiFab &metric,
                  const int &level, const CCTK_REAL *real_params,
                  const CCTK_INT *int_params) {
    initializer_function(static_cast<OtherContainer &>(*this), metric, level,
                         real_params, int_params);
  };
 
  virtual void evolve(const amrex::MultiFab &lapse,
                      const amrex::MultiFab &shift,
                      const amrex::MultiFab &metric,
                      const amrex::MultiFab &curv, const CCTK_REAL &dt,
                      const int &lev) = 0;
 
  void check_banned_zones(const int &level, const CCTK_INT4 &zones,
                          const CCTK_REAL (&x)[10], const CCTK_REAL (&y)[10],
                          const CCTK_REAL (&z)[10],
                          const CCTK_REAL (&radius)[10]) {
 
    if (!zones) {
      return;
    }
 
    for (ParticleIterator<StructType> pti(*this, level); pti.isValid(); ++pti) {
      const int np = pti.numParticles();
      auto *AMREX_RESTRICT particles = &(pti.GetArrayOfStructs()[0]);
 
      auto self = this;
      amrex::ParallelFor(np, [=] AMREX_GPU_DEVICE(int i) noexcept {
        bool out = false;
        for (int check = 0; check < zones; check++) {
          const CCTK_REAL dx = particles[i].pos(0) - x[check];
          const CCTK_REAL dy = particles[i].pos(1) - y[check];
          const CCTK_REAL dz = particles[i].pos(2) - z[check];
          out |= (dx * dx + dy * dy + dz * dz <= radius[check] * radius[check]);
        }
        if (out) {
          particles[i].id() = -1;
          return;
        }
      });
    }
  }
 
  void outputParticlesAscii(const int &it, const int &plot_every,
                            const std::string &out_dir) {
    if (plot_every > 0 && it % plot_every == 0) {
      const std::string &file_name =
          out_dir + "/" + amrex::Concatenate(this->name, it);
      CCTK_VINFO(" Writing ascii file %s", file_name.c_str());
 
      this->WriteAsciiFile(file_name);
    }
  };
 
  void outputParticlesPlot(const int &it, const int &plot_every,
                           const std::string &out_dir) {
    if (plot_every > 0 && it % plot_every == 0) {
      const std::string file_name =
          out_dir + "/" + amrex::Concatenate("plt", it);
      CCTK_VINFO("Writing plot file %s", file_name.c_str());
 
      this->WritePlotFile(file_name, this->name);
    }
  };
 
}; // class BaseParticlesContainer
 
} // namespace GInX
 
#endif // !BASEPARTICLESCONTAINER_HXX