1 /* ----------------------------------------------------------------------
2    SPARTA - Stochastic PArallel Rarefied-gas Time-accurate Analyzer
3    http://sparta.sandia.gov
4    Steve Plimpton, sjplimp@sandia.gov, Michael Gallis, magalli@sandia.gov
5    Sandia National Laboratories
6 
7    Copyright (2014) Sandia Corporation.  Under the terms of Contract
8    DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
9    certain rights in this software.  This software is distributed under
10    the GNU General Public License.
11 
12    See the README file in the top-level SPARTA directory.
13 ------------------------------------------------------------------------- */
14 
15 #ifndef SPARTA_UPDATE_KOKKOS_H
16 #define SPARTA_UPDATE_KOKKOS_H
17 
18 #include "update.h"
19 #include "kokkos_type.h"
20 #include "particle.h"
21 #include "grid_kokkos.h"
22 #include "domain_kokkos.h"
23 #include "kokkos_copy.h"
24 #include "surf_collide_diffuse_kokkos.h"
25 #include "surf_collide_specular_kokkos.h"
26 #include "surf_collide_vanish_kokkos.h"
27 #include "surf_collide_piston_kokkos.h"
28 #include "surf_collide_transparent_kokkos.h"
29 #include "compute_boundary_kokkos.h"
30 #include "compute_surf_kokkos.h"
31 
32 namespace SPARTA_NS {
33 
34 #define KOKKOS_SURF_COLL_TYPES 4
35 #define KOKKOS_MAX_SURF_COLL_PER_TYPE 2
36 #define KOKKOS_TOT_SURF_COLL 6
37 #define KOKKOS_MAX_BLIST 2
38 #define KOKKOS_MAX_SLIST 2
39 
40 struct s_UPDATE_REDUCE {
41   int ntouch_one,nexit_one,nboundary_one,
42       entryexit,ncomm_one,
43       nscheck_one,nscollide_one,nreact_one,nstuck,
44       naxibad,error_flag;
45   KOKKOS_INLINE_FUNCTION
s_UPDATE_REDUCEs_UPDATE_REDUCE46   s_UPDATE_REDUCE() {
47     ntouch_one    = 0;
48     nexit_one     = 0;
49     nboundary_one = 0;
50     ncomm_one     = 0;
51     nscheck_one   = 0;
52     nscollide_one = 0;
53     nreact_one    = 0;
54     nstuck        = 0;
55     naxibad       = 0;
56   }
57 
58   KOKKOS_INLINE_FUNCTION
59   void operator+=(const s_UPDATE_REDUCE &rhs) {
60     ntouch_one    += rhs.ntouch_one   ;
61     nexit_one     += rhs.nexit_one    ;
62     nboundary_one += rhs.nboundary_one;
63     ncomm_one     += rhs.ncomm_one    ;
64     nscheck_one   += rhs.nscheck_one  ;
65     nscollide_one += rhs.nscollide_one;
66     nreact_one    += rhs.nreact_one   ;
67     nstuck        += rhs.nstuck       ;
68     naxibad       += rhs.naxibad      ;
69   }
70 
71   KOKKOS_INLINE_FUNCTION
72   void operator+=(const volatile s_UPDATE_REDUCE &rhs) volatile {
73     ntouch_one    += rhs.ntouch_one   ;
74     nexit_one     += rhs.nexit_one    ;
75     nboundary_one += rhs.nboundary_one;
76     ncomm_one     += rhs.ncomm_one    ;
77     nscheck_one   += rhs.nscheck_one  ;
78     nscollide_one += rhs.nscollide_one;
79     nreact_one    += rhs.nreact_one   ;
80     nstuck        += rhs.nstuck       ;
81     naxibad       += rhs.naxibad      ;
82   }
83 };
84 typedef struct s_UPDATE_REDUCE UPDATE_REDUCE;
85 
86 template<int DIM, int SURF, int ATOMIC_REDUCTION>
87 struct TagUpdateMove{};
88 
89 class UpdateKokkos : public Update {
90  public:
91   typedef UPDATE_REDUCE value_type;
92 
93   DAT::tdual_int_1d k_mlist;
94   DAT::tdual_int_1d k_mlist_small;
95   //DAT::t_int_1d d_mlist_small;
96   //HAT::t_int_scalar h_mlist_small;
97   //int* mlist_small;
98 
99   UpdateKokkos(class SPARTA *);
100   ~UpdateKokkos();
101   void init();
102   void setup();
103   void run(int);
104 
105   template<int DIM, int SURF, int ATOMIC_REDUCTION>
106   KOKKOS_INLINE_FUNCTION
107   void operator()(TagUpdateMove<DIM,SURF,ATOMIC_REDUCTION>, const int&) const;
108 
109   template<int DIM, int SURF, int ATOMIC_REDUCTION>
110   KOKKOS_INLINE_FUNCTION
111   void operator()(TagUpdateMove<DIM,SURF,ATOMIC_REDUCTION>, const int&, UPDATE_REDUCE&) const;
112 
113  private:
114 
115   double dt;
116 
117   t_cell_1d d_cells;
118   t_sinfo_1d d_sinfo;
119   t_pcell_1d d_pcells;
120 
121   Kokkos::Crs<int, DeviceType, void, int> d_csurfs;
122   Kokkos::Crs<int, DeviceType, void, int> d_csplits;
123   Kokkos::Crs<int, DeviceType, void, int> d_csubs;
124 
125   t_line_1d d_lines;
126   t_tri_1d d_tris;
127 
128   t_particle_1d d_particles;
129 
130   KKCopy<GridKokkos> grid_kk_copy;
131   KKCopy<DomainKokkos> domain_kk_copy;
132 
133   int sc_type_list[KOKKOS_TOT_SURF_COLL];
134   int sc_map[KOKKOS_TOT_SURF_COLL];
135   KKCopy<SurfCollideSpecularKokkos> sc_kk_specular_copy[KOKKOS_MAX_SURF_COLL_PER_TYPE];
136   KKCopy<SurfCollideDiffuseKokkos> sc_kk_diffuse_copy[KOKKOS_MAX_SURF_COLL_PER_TYPE];
137   KKCopy<SurfCollideVanishKokkos> sc_kk_vanish_copy[KOKKOS_MAX_SURF_COLL_PER_TYPE];
138   KKCopy<SurfCollidePistonKokkos> sc_kk_piston_copy[KOKKOS_MAX_SURF_COLL_PER_TYPE];
139   KKCopy<SurfCollideTransparentKokkos> sc_kk_transparent_copy[KOKKOS_MAX_SURF_COLL_PER_TYPE];
140   KKCopy<ComputeBoundaryKokkos> blist_active_copy[KOKKOS_MAX_BLIST];
141   KKCopy<ComputeSurfKokkos> slist_active_copy[KOKKOS_MAX_SLIST];
142 
143 
144   typedef Kokkos::DualView<int[12], DeviceType::array_layout, DeviceType> tdual_int_12;
145   typedef tdual_int_12::t_dev t_int_12;
146   typedef tdual_int_12::t_host t_host_int_12;
147   t_int_12 d_scalars;
148   t_host_int_12 h_scalars;
149 
150   DAT::t_int_scalar d_ntouch_one;
151   HAT::t_int_scalar h_ntouch_one;
152 
153   DAT::t_int_scalar d_nexit_one;
154   HAT::t_int_scalar h_nexit_one;
155 
156   DAT::t_int_scalar d_nboundary_one;
157   HAT::t_int_scalar h_nboundary_one;
158 
159   DAT::t_int_scalar d_nmigrate;
160   HAT::t_int_scalar h_nmigrate;
161 
162   DAT::t_int_scalar d_entryexit;
163   HAT::t_int_scalar h_entryexit;
164 
165   DAT::t_int_scalar d_ncomm_one;
166   HAT::t_int_scalar h_ncomm_one;
167 
168   DAT::t_int_scalar d_nscheck_one;
169   HAT::t_int_scalar h_nscheck_one;
170 
171   DAT::t_int_scalar d_nscollide_one;
172   HAT::t_int_scalar h_nscollide_one;
173 
174   DAT::t_int_scalar d_nreact_one;
175   HAT::t_int_scalar h_nreact_one;
176 
177   DAT::t_int_scalar d_nstuck;
178   HAT::t_int_scalar h_nstuck;
179 
180   DAT::t_int_scalar d_naxibad;
181   HAT::t_int_scalar h_naxibad;
182 
183   DAT::t_int_scalar d_error_flag;
184   HAT::t_int_scalar h_error_flag;
185 
186   void bounce_set(bigint);
187 
188   // remap x and v components into axisymmetric plane
189   // input x at end of linear move (x = xold + dt*v)
190   // change x[1] = sqrt(x[1]^2 + x[2]^2), x[2] = 0.0
191   // change vy,vz by rotation into axisymmetric plane
192 
193   KOKKOS_INLINE_FUNCTION
axi_remap(double * x,double * v)194   void axi_remap(double *x, double *v) const {
195     double ynew = x[1];
196     double znew = x[2];
197     x[1] = sqrt(ynew*ynew + znew*znew);
198     x[2] = 0.0;
199     double rn = ynew / x[1];
200     double wn = znew / x[1];
201     double vy = v[1];
202     double vz = v[2];
203     v[1] = vy*rn + vz*wn;
204     v[2] = -vy*wn + vz*rn;
205   };
206 
207   typedef void (UpdateKokkos::*FnPtr)();
208   FnPtr moveptr;             // ptr to move method
209   template < int, int > void move();
210 
211   //
212   //int perturbflag;
213   typedef void (UpdateKokkos::*FnPtr2)(double, double *, double *) const;
214   FnPtr2 moveperturb;        // ptr to moveperturb method
215   //
216   //// variants of moveperturb method
217   //// adjust end-of-move x,v due to perturbation on straight-line advection
218 
219   KOKKOS_INLINE_FUNCTION
220   int split3d(int, double*) const;
221 
222   KOKKOS_INLINE_FUNCTION
223   int split2d(int, double*) const;
224 
225   int gravity_3d_flag,gravity_2d_flag;
226 
227   // variants of moveperturb method
228   // adjust end-of-move x,v due to perturbation on straight-line advection
229 
230   KOKKOS_INLINE_FUNCTION
gravity2d(double dt,double * x,double * v)231   void gravity2d(double dt, double *x, double *v) const {
232     double dtsq = 0.5*dt*dt;
233     x[0] += dtsq*gravity[0];
234     x[1] += dtsq*gravity[1];
235     v[0] += dt*gravity[0];
236     v[1] += dt*gravity[1];
237   };
238 
239   KOKKOS_INLINE_FUNCTION
gravity3d(double dt,double * x,double * v)240   void gravity3d(double dt, double *x, double *v) const {
241     double dtsq = 0.5*dt*dt;
242     x[0] += dtsq*gravity[0];
243     x[1] += dtsq*gravity[1];
244     x[2] += dtsq*gravity[2];
245     v[0] += dt*gravity[0];
246     v[1] += dt*gravity[1];
247     v[2] += dt*gravity[2];
248   };
249 };
250 
251 }
252 
253 #endif
254 
255 /* ERROR/WARNING messages:
256 
257 */
258