ERF
Energy Research and Forecasting: An Atmospheric Modeling Code
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
ERFPC.H
Go to the documentation of this file.
1 #ifndef ERF_PC_H_
2 #define ERF_PC_H_
3 
4 #ifdef ERF_USE_PARTICLES
5 
6 #include <string>
7 #include <AMReX_Particles.H>
8 
9 struct ERFParticlesIntIdxAoS
10 {
11  enum {
12  k = 0,
13  ncomps
14  };
15 };
16 
17 struct ERFParticlesRealIdxAoS
18 {
19  enum {
20  ncomps = 0
21  };
22 };
23 
24 struct ERFParticlesIntIdxSoA
25 {
26  enum {
27  ncomps = 0
28  };
29 };
30 
31 struct ERFParticlesRealIdxSoA
32 {
33  enum {
34  vx = 0,
35  vy,
36  vz,
37  mass,
38  temperature,
39  ncomps
40  };
41 };
42 
43 namespace ERFParticleInitializations
44 {
45  /* list of particle initializations */
46  const std::string init_box_uniform = "box";
47 }
48 
49 namespace ERFParticleNames
50 {
51  const std::string tracers = "tracer_particles";
52 }
53 
54 struct ERFParticlesAssignor
55 {
56  template <typename P>
57  AMREX_GPU_HOST_DEVICE
58  amrex::IntVect operator() ( P const& p,
59  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& plo,
60  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& dxi,
61  const amrex::Box& domain ) const noexcept
62  {
63  amrex::IntVect iv(
64  AMREX_D_DECL( int(amrex::Math::floor((p.pos(0)-plo[0])*dxi[0])),
65  int(amrex::Math::floor((p.pos(1)-plo[1])*dxi[1])),
66  p.idata(ERFParticlesIntIdxAoS::k) ) );
67  iv[0] += domain.smallEnd()[0];
68  iv[1] += domain.smallEnd()[1];
69  return iv;
70  }
71 };
72 
73 template <typename P>
74 AMREX_GPU_HOST_DEVICE AMREX_FORCE_INLINE
75 void update_location_idata ( P& a_p,
76  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& a_plo,
77  amrex::GpuArray<amrex::Real,AMREX_SPACEDIM> const& a_dxi,
78  const amrex::Array4<amrex::Real const>& a_height_arr )
79 {
80  amrex::IntVect iv( int(amrex::Math::floor((a_p.pos(0)-a_plo[0])*a_dxi[0])),
81  int(amrex::Math::floor((a_p.pos(1)-a_plo[1])*a_dxi[1])),
82  a_p.idata(ERFParticlesIntIdxAoS::k) );
83 
84  if (a_height_arr) {
85  amrex::Real lx = (a_p.pos(0)-a_plo[0])*a_dxi[0] - static_cast<amrex::Real>(iv[0]);
86  amrex::Real ly = (a_p.pos(1)-a_plo[1])*a_dxi[1] - static_cast<amrex::Real>(iv[1]);
87  auto zlo = a_height_arr(iv[0] ,iv[1] ,iv[2] ) * (1.0-lx) * (1.0-ly) +
88  a_height_arr(iv[0]+1,iv[1] ,iv[2] ) * lx * (1.0-ly) +
89  a_height_arr(iv[0] ,iv[1]+1,iv[2] ) * (1.0-lx) * ly +
90  a_height_arr(iv[0]+1,iv[1]+1,iv[2] ) * lx * ly;
91  auto zhi = a_height_arr(iv[0] ,iv[1] ,iv[2]+1) * (1.0-lx) * (1.0-ly) +
92  a_height_arr(iv[0]+1,iv[1] ,iv[2]+1) * lx * (1.0-ly) +
93  a_height_arr(iv[0] ,iv[1]+1,iv[2]+1) * (1.0-lx) * ly +
94  a_height_arr(iv[0]+1,iv[1]+1,iv[2]+1) * lx * ly;
95 
96  if (a_p.pos(2) > zhi) {
97  a_p.idata(ERFParticlesIntIdxAoS::k) += 1;
98  } else if (a_p.pos(2) <= zlo) {
99  a_p.idata(ERFParticlesIntIdxAoS::k) -= 1;
100  }
101  }
102 }
103 
104 class ERFPC : public amrex::ParticleContainer< ERFParticlesRealIdxAoS::ncomps, // AoS real attributes
105  ERFParticlesIntIdxAoS::ncomps, // AoS integer attributes
106  ERFParticlesRealIdxSoA::ncomps, // SoA real attributes
107  ERFParticlesIntIdxSoA::ncomps, // SoA integer attributes
108  amrex::DefaultAllocator,
109  ERFParticlesAssignor >
110 {
111  public:
112 
113  /*! Constructor */
114  ERFPC ( amrex::ParGDBBase* a_gdb,
115  const std::string& a_name = "particles" )
116  : amrex::ParticleContainer< ERFParticlesRealIdxAoS::ncomps, // AoS real attributes
117  ERFParticlesIntIdxAoS::ncomps, // AoS integer attributes
118  ERFParticlesRealIdxSoA::ncomps, // SoA real attributes
119  ERFParticlesIntIdxSoA::ncomps, // SoA integer attributes
120  amrex::DefaultAllocator,
121  ERFParticlesAssignor> (a_gdb)
122  {
123  BL_PROFILE("ERFPCPC::ERFPC()");
124  m_name = a_name;
125  readInputs();
126  }
127 
128  /*! Constructor */
129  ERFPC ( const amrex::Geometry& a_geom,
130  const amrex::DistributionMapping& a_dmap,
131  const amrex::BoxArray& a_ba,
132  const std::string& a_name = "particles" )
133  : amrex::ParticleContainer< ERFParticlesRealIdxAoS::ncomps, // AoS real attributes
134  ERFParticlesIntIdxAoS::ncomps, // AoS real attributes
135  ERFParticlesRealIdxSoA::ncomps, // SoA real attributes
136  ERFParticlesIntIdxSoA::ncomps, // SoA integer attributes
137  amrex::DefaultAllocator,
138  ERFParticlesAssignor> ( a_geom, a_dmap, a_ba )
139  {
140  BL_PROFILE("ERFPCPC::ERFPC()");
141  m_name = a_name;
142  readInputs();
143  }
144 
145  /*! Initialize particles in domain */
146  virtual void InitializeParticles (const amrex::Real time, const std::unique_ptr<amrex::MultiFab>& a_ptr = nullptr);
147 
148  /*! Evolve particles for one time step */
149  virtual void EvolveParticles ( int,
150  amrex::Real,
151  amrex::Vector<amrex::Vector<amrex::MultiFab>>&,
152  const amrex::Vector<std::unique_ptr<amrex::MultiFab>>& );
153 
154  /*! Get real-type particle attribute names */
155  virtual amrex::Vector<std::string> varNames () const
156  {
157  BL_PROFILE("ERFPCPC::varNames()");
158  return {AMREX_D_DECL("xvel","yvel","zvel"),"mass","temperature"};
159  }
160 
161  /*! Get real-type particle attribute names */
162  virtual amrex::Vector<std::string> meshPlotVarNames () const
163  {
164  BL_PROFILE("ERFPCPC::varNames()");
165  return {"mass_density"};
166  }
167 
168  /*! Uses midpoint method to advance particles using flow velocity. */
169  virtual void AdvectWithFlow ( amrex::MultiFab*,
170  int,
171  amrex::Real,
172  const std::unique_ptr<amrex::MultiFab>& );
173 
174  /*! Uses midpoint method to advance particles falling under gravity. */
175  virtual void AdvectWithGravity ( int,
176  amrex::Real,
177  const std::unique_ptr<amrex::MultiFab>& );
178 
179  /*! Interpolates flow temperature to particles */
180  virtual void ComputeTemperature ( const amrex::MultiFab&,
181  int,
182  amrex::Real,
183  const std::unique_ptr<amrex::MultiFab>& );
184 
185  /*! Compute mass density */
186  virtual void massDensity ( amrex::MultiFab&, const int&, const int& a_comp = 0) const;
187 
188  /*! Compute mesh variable from particles */
189  virtual void computeMeshVar( const std::string& a_var_name,
190  amrex::MultiFab& a_mf,
191  const int a_lev) const
192  {
193  if (a_var_name == "mass_density") {
194  massDensity( a_mf, a_lev );
195  } else {
196  a_mf.setVal(0.0);
197  }
198  }
199 
200  /*! Specify if particles should advect with flow */
201  inline void setAdvectWithFlow (bool a_flag)
202  {
203  BL_PROFILE("ERFPCPC::setAdvectWithFlow()");
204  m_advect_w_flow = a_flag;
205  }
206  /*! Specify if particles fall under gravity */
207  inline void setAdvectWithGravity (bool a_flag)
208  {
209  BL_PROFILE("ERFPCPC::setAdvectWithGravity()");
210  m_advect_w_gravity = a_flag;
211  }
212 
213  // the following functions should ideally be private or protected, but need to be
214  // public due to CUDA extended lambda capture rules
215 
216  /*! Default particle initialization */
217  void initializeParticlesUniformDistributionInBox (const std::unique_ptr<amrex::MultiFab>& a_ptr,
218  const amrex::RealBox& particle_box);
219 
220  protected:
221 
222  bool m_advect_w_flow; /*!< advect with flow velocity */
223  bool m_advect_w_gravity; /*!< advect under gravitational force */
224 
225  amrex::RealBox m_particle_box; /*!< box within which to place particles */
226 
227  std::string m_name; /*!< name of this particle species */
228 
229  std::string m_initialization_type; /*!< initial particle distribution type */
230  int m_ppc_init; /*!< initial number of particles per cell */
231 
232  amrex::Real m_inject_start_time; /*!< particles will only be initialized/injected if time >= m_inject_start_time */
233 
234  bool m_initialized = false; /*!< have the particles been initialized */
235 
236  bool m_stable_redistribute; /*!< use stable redistribute for deterministic simulations */
237 
238  /*! read inputs from file */
239  virtual void readInputs ();
240 
241  private:
242 
243  bool place_randomly_in_cells; /*!< place particles at random positions? */
244 };
245 
246 #endif
247 #endif
Definition: ERF_ConsoleIO.cpp:12