particles.h

Go to the documentation of this file.

#ifndef PARTICLES_H
#define PARTICLES_H
 
#include <algorithm>
#include <functional>
#include <iomanip>
#include <iostream>
#include <quadgrid_config.h>
#include <json.hpp>
#include <map>
#include <quadgrid_cpp.h>
#include <string>


 
struct
particles_t {

  using idx_t = quadgrid_t<vector_t<real_t>>::idx_t;
 
  idx_t num_particles;    
  vector_t<real_t> x;  
  vector_t<real_t> y;  

  std::map<std::string, vector_t<idx_t>> iprops;

  std::map<std::string, vector_t<real_t>> dprops;
 
  vector_t<real_t> M; 
  std::map<idx_t, vector_t<idx_t>> grd_to_ptcl;   
  vector_t<idx_t> ptcl_to_grd;                    
  vector_t<idx_t> ptcl_grd_color;                 
  const quadgrid_t<vector_t<real_t>>& grid;       
 
  #ifdef USE_THRUST
  device_vector_t<idx_t> device_ptcl_to_grd;
  device_vector_t<real_t> device_x, device_y;
 
  std::map<std::string, device_vector_t<real_t>> device_dprops;
  std::map<std::string, device_vector_t<idx_t>> device_iprops;
 
  device_vector_t<real_t> device_grid_M; //mass matrix
  std::map<std::string, device_vector_t<real_t>> device_grid_vars;
 
  #endif

  enum class
  output_format : idx_t {
    csv = 0,           
 
    octave_ascii = 1,  
 
    json = 2           
  };

  enum class
  cell_color : idx_t {
    red   = 0,
    green = 1,
    blue  = 2,
    black = 3,
    all_colors = 4
  };
 
  

  real_t
  default_x_generator ();


  real_t
  default_y_generator ();


  template<output_format fmt>
  void
  print (std::ostream & os) const {
    os << "output format not implementd" << std::endl;
  }


  particles_t (idx_t n, const quadgrid_t<vector_t<real_t>>& grid_)
    : num_particles(n), grid(grid_) { }


  particles_t (const nlohmann::json &j,
               const quadgrid_t<vector_t<real_t>>& grid_)
    :  grid(grid_)
  {
    j["dprops"].get_to<std::map<std::string, vector_t<real_t>>> (dprops);
    j["iprops"].get_to<std::map<std::string, vector_t<int>>> (iprops);
    j["x"].get_to<vector_t<real_t>> (x);
    j["y"].get_to<vector_t<real_t>> (y);
    j["num_particles"].get_to<idx_t> (num_particles);
  }


  particles_t (idx_t n, const vector_t<std::string>& ipropnames,
               const vector_t<std::string>& dpropnames,
               const quadgrid_t<vector_t<real_t>>& grid_);


  particles_t (idx_t n, const vector_t<std::string>& ipropnames,
               const vector_t<std::string>& dpropnames,
               const quadgrid_t<vector_t<real_t>>& grid_,
               const vector_t<real_t> & xv,
               const vector_t<real_t> & yv);


  particles_t (idx_t n, const vector_t<std::string>& ipropnames,
               const vector_t<std::string>& dpropnames,
               const quadgrid_t<vector_t<real_t>>& grid_,
               std::function<real_t ()> xgen,
               std::function<real_t ()> ygen);


  void
  init_props (const vector_t<std::string>& ipropnames,
              const vector_t<std::string>& dpropnames);


  void
  remove_in_region (std::function<bool (real_t, real_t)> fun) {
    vector_t<idx_t> vin{};
    for (idx_t i = 0; i < num_particles; ++i) {
      if (fun (x[i], y[i])) {
        vin.push_back (i);
      }
    }
 
    for (auto &dprop : dprops) {
      for (auto id = vin.rbegin (); id != vin.rend (); ++id) {
        dprop.second.erase (dprop.second.begin () + (*id));
      }
    }
 
    for (auto &iprop : iprops) {
      for (auto id = vin.rbegin (); id != vin.rend (); ++id) {
        iprop.second.erase (iprop.second.begin () + (*id));
      }
    }
 
    for (auto id = vin.rbegin (); id != vin.rend (); ++id) {
      x.erase (x.begin () + (*id));
      y.erase (y.begin () + (*id));
      --num_particles;
    }
  };

  void
  memcpy_host_to_device ();

  void
  memcpy_device_to_host ();


  void
  init_particle_mesh ();

  //  `grd_to_ptcl`, if not needed.
  // backend is the policy
  template<class backend>
  void
  update_ptcl_to_grd ();
  
  //policies for update_ptcl_to_grd()
  class 
  update_ptcl_to_grd_host;
  
  class
  update_ptcl_to_grd_device;
  
 
  void
  mark_by_cell_color ();

  void
  reorder (vector_t<idx_t> &);


  void
  init_particle_positions (std::function<real_t ()> xgentr,
                           std::function<real_t ()> ygentr);


  void
  build_mass ();

  real_t &
  dp (const std::string & name, idx_t ii) {
    return dprops.at (name) [ii];
  }

  const real_t &
  dp (const std::string & name, idx_t ii) const {
    return dprops.at (name) [ii];
  }

  idx_t &
  ip (const std::string & name, idx_t ii) {
    return iprops.at (name) [ii];
  }

  const idx_t &
  ip (const std::string & name, idx_t ii) const {
    return iprops.at (name) [ii];
  }
 
  static
  const std::string &
  getkey(std::map<std::string, device_vector_t<real_t>> const &varnames,
         std::size_t ivar)  {
    return std::next (varnames.begin (), ivar)->first;
  };
 
  static
  const std::string &
  getkey(vector_t<std::string> const &varnames,
         std::size_t ivar)  {
    return varnames[ivar];
  };
 
  static
  const char*
  getkey(std::initializer_list<const char *> const &varnames,
         std::size_t ivar)  {
    return *(std::next (varnames.begin (), ivar));
  };


  void
  p2g (std::map<std::string, device_vector_t<real_t>> & vars,
       bool apply_mass = false)  {
    p2g (vars, vars, vars, apply_mass);
  }


  template<typename GT, typename PT>
  void
  p2g (std::map<std::string, device_vector_t<real_t>> & vars,
       PT const & pvarnames,
       GT const & gvarnames,
       bool apply_mass = false) ;
 
  template<typename str>
  void
  p2g (std::map<std::string, device_vector_t<real_t>> & vars,
       std::initializer_list<str> const & pvarnames,
       std::initializer_list<str> const & gvarnames,
       bool apply_mass = false) ;
 
  //p2g with custom helper
  template<typename UnaryFunction>
  void
  p2g (std::map<std::string, device_vector_t<real_t>> & vars,
       UnaryFunction helper);
 
  template<typename GT, typename PT>
  void
  p2g_print (std::map<std::string, vector_t<real_t>> & vars,
       PT const & pvarnames,
       GT const & gvarnames,
       bool apply_mass = false) const;
 
  template<typename str>
  void
  p2g_print (std::map<std::string, vector_t<real_t>> & vars,
       std::initializer_list<str> const & pvarnames,
       std::initializer_list<str> const & gvarnames,
       bool apply_mass = false) const;
 
  template<typename GT, typename PT>
  void
  p2gd (std::map<std::string, device_vector_t<real_t>> & vars,
        PT const & pxvarnames,
        PT const & pyvarnames,
        std::string const &area,
        GT const & gvarnames,
        bool apply_mass = false);
 
  template<typename str>
  void
  p2gd (std::map<std::string, device_vector_t<real_t>> & vars,
        std::initializer_list<str> const & pxvarnames,
        std::initializer_list<str> const & pyvarnames,
        std::string const & area,
        std::initializer_list<str> const & gvarnames,
        bool apply_mass = false);
 
  //p2gd with custom helper
  template<typename UnaryFunction>
  void 
  p2gd (std::map<std::string, device_vector_t<real_t>> & vars,
        UnaryFunction helper);
 
  void
  g2p (const std::map<std::string, device_vector_t<real_t>>& vars,
       bool apply_mass = false) {
    g2p (vars, vars, vars, apply_mass);
  }
 
  template<typename str>
  void
  g2p (const std::map<std::string, device_vector_t<real_t>>& vars,
       std::initializer_list<str> const & gvarnames,
       std::initializer_list<str> const & pvarnames,
       bool apply_mass = false);
 
  template<typename GT, typename PT>
  void
  g2p (const std::map<std::string, device_vector_t<real_t>>& vars,
       GT const & gvarnames,
       PT const & pvarnames,
       bool apply_mass = false);
 
  //g2p with custom helper
  template<typename UnaryFunction>
  void 
  g2p (std::map<std::string, device_vector_t<real_t>> & vars,
       UnaryFunction helper);
 
  template<typename GT, typename PT>
  void
  g2pd (const std::map<std::string, device_vector_t<real_t>>& vars,
        GT const & gvarnames,
        PT const & pxvarnames,
        PT const & pyvarnames,
        bool apply_mass = false);
 
  template<typename str>
  void
  g2pd (const std::map<std::string, device_vector_t<real_t>>& vars,
        std::initializer_list<str> const & gvarnames,
        std::initializer_list<str> const &pxvarnames,
        std::initializer_list<str> const & pyvarnames,
        bool apply_mass = false);
 
  //g2pd with custom helper
  template<typename UnaryFunction>
  void 
  g2pd (std::map<std::string, device_vector_t<real_t>> & vars,
        UnaryFunction helper);
 
};

void
to_json (nlohmann::json &j, const particles_t &p);

template<typename P2G_t, typename COORD_t>
class
ptcl_to_grd_update_t {
 
  using idx_t=particles_t::idx_t;
  P2G_t ptcl_to_grd;
  const COORD_t x;
  const COORD_t y;
  const real_t hx;
  const real_t hy;
  const idx_t nrows;
  
public :
  ptcl_to_grd_update_t (P2G_t ptcl_to_grd_,
             const COORD_t x_, const COORD_t y_,
             const real_t hx_, const real_t hy_, const idx_t nrows_)
    : ptcl_to_grd(ptcl_to_grd_), x(x_), y(y_), hx(hx_), hy(hy_), nrows(nrows_) {};
 
    HOST DEVICE 
    void operator() (particles_t::idx_t ii) {
    ptcl_to_grd[ii] = quadgrid_t<COORD_t>::sub2gind (static_cast<idx_t> (std::floor (y[ii] / hy)),
                             static_cast<idx_t> (std::floor (x[ii] / hx)),
                             nrows);
  }
};
 
 
#include "particles_imp.h"
 
#endif /* PARTICLES_H */