XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
grid_field_pack.hpp
Go to the documentation of this file.
1 #ifndef GRID_FIELD_PACK_HPP
2 #define GRID_FIELD_PACK_HPP
3 #include "space_settings.hpp"
4 #include "NamelistReader.hpp"
5 
6 #include "globals.hpp"
7 #include "grid.hpp"
8 #include "linear_weights.hpp"
9 #include "local_fields.hpp"
10 #include "push_controls.hpp"
11 #include "gyro_radius.hpp"
12 #include "grid_field.hpp"
13 #include "label.hpp"
14 #include "pack.hpp"
15 
16 // Base struct for GridFieldPack so that non-GPU code doesn't need to be templated just because there is a GridFieldPack argument
18  virtual ~GridFieldPackGeneric() = default; // Virtual destructor
19 };
20 
21 using GridFieldPackPtr = std::unique_ptr<GridFieldPackGeneric>;
22 
23 template<class Device, typename GFPT>
25  using pack_type = GFPT;
26 
27  bool turb_efield; // Ideally this would not be here but it is tedious to pass into the gather otherwise
28 
29  int phi_offset;
31 
32  // Views stored in Pack
34 
35  // Const access to pack contents for inside kernels (can modify data but not structure)
36  template<Label FN>
37  KOKKOS_INLINE_FUNCTION const auto& get() const{
38  return pack.template get<FN>();
39  }
40 
41  // Non-const access to pack contents (can modify structure)
42  template<Label FN>
43  inline auto& get(){
44  return pack.template get<FN>();
45  }
46 
47  // This will be a ScalarGridField in case of XGCa and a grid field with field-following
48  // interpolation in case of XGC1
50 
52 
53  GridFieldPack(bool turb_efield_in)
54  : turb_efield(turb_efield_in),
55  phi_offset(0),
56  node_offset(0) {}
57 
58  // Device functions:
59 
60  KOKKOS_INLINE_FUNCTION void phi_from_para(SimdVector& vec, int i_simd, const SimdVector& B, double Bmag) const{
61  vec.phi[i_simd]=(vec.phi[i_simd]*Bmag - vec.r[i_simd]*B.r[i_simd] - vec.z[i_simd]*B.z[i_simd]) / B.phi[i_simd];
62  }
63 
64 #ifndef XGCA
65  template<typename FieldT>
66  KOKKOS_INLINE_FUNCTION void gather_deltaf_conv(int i_simd, int node, int i_phi, double wp, const FieldCorrection& corr, const SimdGridWeights<Order::One, PIT_GLOBAL>& grid_wts, LocalFields<GFPT>& fld, const FieldT& labeled_field) const {
67  if constexpr(FieldT::name == Label::E00){
68  get<Label::E00>()(node).gather(fld.template get<Label::E00>(),i_simd, wp, grid_wts.phi_wts, corr);
69  get<Label::ddpotdt>()(i_phi,node).gather(fld.template get<Label::ddpotdt>(),i_simd, wp, grid_wts.phi_wts);
70  }
71  }
72 #endif
73 
74  KOKKOS_INLINE_FUNCTION void gather_all_fields(const PushControls &push_controls, int i_simd, int i_node, double wp, const FieldCorrection& corr, const SimdGridWeights<Order::One, PIT_GLOBAL>& grid_wts, const SimdGyroWeights<DriftKin>& rho_wts, LocalFields<GFPT>& fld) const {
75  const int node = i_node - node_offset;
76 #ifdef XGCA
77  get<Label::E>()(node,0).gather(fld.template get<Label::E>(),i_simd, wp, corr); // irho = 0
78  // Loop voltage
79  //loop_voltage(node).gather(fld.loop_voltage,i_simd,wp);
80 #else
81  const int i_phi = (grid_wts.phi_wts.phi.i[i_simd] - phi_offset)%get<Label::E>().f.extent(0); //should be grid.nplanes (?);
82  get<Label::E>()(i_phi,node).gather(fld.template get<Label::E>(),i_simd, wp, grid_wts.phi_wts, corr);
83 #ifdef EXPLICIT_EM
84  get<Label::Ah>()(i_phi,node).gather(fld.template get<Label::Ah>(),i_simd, wp, grid_wts.phi_wts);
85  get<Label::dAh>()(i_phi,node).gather(fld.template get<Label::dAh>(),i_simd, wp, grid_wts.phi_wts, corr);
86 
87  if(push_controls.em_mixed_variable) {
88  get<Label::dAs>()(i_phi,node).gather(fld.template get<Label::dAs>(),i_simd, wp, grid_wts.phi_wts, corr);
89  get<Label::As>()(i_phi,node).gather(fld.template get<Label::As>(),i_simd, wp, grid_wts.phi_wts);
90  }
91  if(push_controls.em_mixed_variable && (push_controls.em_pullback_mode==4)) {
92  get<Label::Epar_em>()(i_phi,node).gather(fld.template get<Label::Epar_em>(),i_simd,wp,grid_wts.phi_wts);
93  }
94 #endif
95  // Loop voltage
96  //loop_voltage(node).gather(fld.loop_voltage,i_simd,wp,grid_wts.phi_wts);
97 
98 
99  // Could replace with new function like: if constexpr(pack.is_present<Label::E00>())
100  fld.pack.for_each_labeled([&](auto& labeled_field) {
101  gather_deltaf_conv(i_simd, node, i_phi, wp, corr, grid_wts, fld, labeled_field);
102  });
103 #endif
104  }
105 
106  // Gather all fields (ion)
107  KOKKOS_INLINE_FUNCTION void gather_all_fields(const PushControls &push_controls, int i_simd, int node, double wp, const FieldCorrection& corr, const SimdGridWeights<Order::One, PIT_GLOBAL>& grid_wts, const SimdGyroWeights<GyroKin>& rho_wts, LocalFields<GFPT>& fld) const {
108  int irho = rho_wts.irho(i_simd);
109 #ifdef XGCA
110  get<Label::E>()(node,irho).gather(fld.template get<Label::E>(),i_simd, wp, rho_wts, corr, get<Label::E>()(node,irho+1));
111 # ifdef SONIC_GK
112  get<Label::dEr_B2>()(node,0).gather(fld.template get<Label::dEr_B2>(),i_simd, wp, corr); // irho = 0
113  get<Label::dEz_B2>()(node,0).gather(fld.template get<Label::dEz_B2>(),i_simd, wp, corr); // irho = 0
114  get<Label::du2_E>()(node,0).gather(fld.template get<Label::du2_E>(),i_simd, wp, corr); // irho = 0
115 # endif
116  // Loop voltage
117  //loop_voltage(node).gather(fld.loop_voltage,i_simd,wp);
118 #else
119  get<Label::E>()(node,irho).gather(fld.template get<Label::E>(),i_simd, wp, rho_wts, grid_wts.phi_wts, corr, get<Label::E>()(node,irho+1));
120 #ifdef EXPLICIT_EM
121  get<Label::Ah>()(node,irho).gather(fld.template get<Label::Ah>(),i_simd, wp, rho_wts, grid_wts.phi_wts, get<Label::Ah>()(node,irho+1));
122  get<Label::dAh>()(node,irho).gather(fld.template get<Label::dAh>(),i_simd, wp, rho_wts, grid_wts.phi_wts, corr, get<Label::dAh>()(node,irho+1));
123 
124  if(push_controls.em_mixed_variable) {
125  get<Label::dAs>()(node,irho).gather(fld.template get<Label::dAs>(),i_simd, wp, rho_wts, grid_wts.phi_wts, corr, get<Label::dAs>()(node,irho+1));
126  get<Label::As>()(node,irho).gather(fld.template get<Label::As>(),i_simd, wp, rho_wts, grid_wts.phi_wts, get<Label::As>()(node,irho+1));
127  }
128  if(push_controls.em_mixed_variable && (push_controls.em_pullback_mode==4)) {
129  get<Label::Epar_em>()(node,irho).gather(fld.template get<Label::Epar_em>(),i_simd,wp,rho_wts, grid_wts.phi_wts, get<Label::Epar_em>()(node,irho+1));
130  }
131 #endif
132  // Loop voltage
133  //loop_voltage(node).gather(fld.loop_voltage,i_simd,wp,grid_wts.phi_wts);
134 #endif
135  }
136 
137 
148  template<KinType PT>
149  KOKKOS_INLINE_FUNCTION void fields_at_point(const PushControls &push_controls, const Grid<Device> &grid, const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights<Order::One, PIT_GLOBAL>& grid_wts, const SimdGyroWeights<PT>& rho_wts, LocalFields<GFPT>& fld) const
150  {
151  // Initialize local fields to 0
152  fld.pack.for_each([&](auto& field) { field.zero(); });
153 
154  for (int i_simd = 0; i_simd<SIMD_SIZE; i_simd++){
155  // If particle is inactive, set itr_work=1, phi=0D0, iphi=0
156  int itr_work = grid_wts.is_valid(i_simd) ? grid_wts.itr[i_simd] : 1;
157 
158  // Field basis correction
159  FieldCorrection field_correction(gradpsi, i_simd);
160 
161  // Loop over the 3 vertices of the grid triangle and add up the contributions from each
162  for (int ip = 0; ip<3; ip++){
163  int node = grid_wts.is_valid(i_simd) ? grid.get_node_index(itr_work, ip) : 0;
164  double wp = grid_wts.is_valid(i_simd) ? grid_wts.p[ip][i_simd] : 0.0; // all these trinary ops are probably unnecessary
165 
166  field_correction.set(i_simd, grid.uses_rz_basis(node), gradpsi);
167 
168  gather_all_fields(push_controls, i_simd, node, wp, field_correction, grid_wts, rho_wts, fld);
169  }
170 
171 #ifdef NEOCLASSICAL_TEST
172  fld.template get<Label::E>().phi[i_simd]=0.0;
173 #else
174 #ifdef XGC1
175  if(turb_efield){
176  // Correct phi component of E-field
177  double Bmag = B.magnitude(i_simd);
178 
179  // Get phi from para
180  phi_from_para(fld.template get<Label::E>(), i_simd, B, Bmag);
181 #ifdef EXPLICIT_EM
182  phi_from_para(fld.template get<Label::dAh>(), i_simd, B, Bmag);
183  phi_from_para(fld.template get<Label::dAs>(), i_simd, B, Bmag);
184 #endif
185  }else{
186  fld.template get<Label::E>().phi[i_simd]=0.0;
187 #ifdef EXPLICIT_EM
188  fld.template get<Label::dAh>().phi[i_simd]=0.0;
189  fld.template get<Label::dAs>().phi[i_simd]=0.0;
190 #endif
191  }
192 #endif
193 #endif
194  }
195  }
196 };
197 #endif
bool turb_efield
Definition: grid_field_pack.hpp:27
KOKKOS_INLINE_FUNCTION void set(int i_simd, double basis, const SimdVector2D &gradpsi)
Definition: field.hpp:35
KOKKOS_INLINE_FUNCTION void gather_all_fields(const PushControls &push_controls, int i_simd, int node, double wp, const FieldCorrection &corr, const SimdGridWeights< Order::One, PIT_GLOBAL > &grid_wts, const SimdGyroWeights< GyroKin > &rho_wts, LocalFields< GFPT > &fld) const
Definition: grid_field_pack.hpp:107
Simd< double > r
Definition: simd.hpp:150
Definition: simd.hpp:149
std::unique_ptr< GridFieldPackGeneric > GridFieldPackPtr
Definition: grid_field_pack.hpp:21
GridFieldPack(bool turb_efield_in)
Definition: grid_field_pack.hpp:53
KOKKOS_INLINE_FUNCTION int irho(int i_simd) const
Definition: gyro_radius.hpp:109
Definition: grid_weights.hpp:47
Definition: push_controls.hpp:9
int node_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:30
Definition: grid.hpp:21
pack_type pack
Definition: local_fields.hpp:36
Definition: grid_field_pack.hpp:24
int phi_offset
Offset for phi_ff field decomposition.
Definition: grid_field_pack.hpp:29
KOKKOS_INLINE_FUNCTION int uses_rz_basis(const int inode) const
Definition: grid.tpp:252
KOKKOS_INLINE_FUNCTION void fields_at_point(const PushControls &push_controls, const Grid< Device > &grid, const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights< Order::One, PIT_GLOBAL > &grid_wts, const SimdGyroWeights< PT > &rho_wts, LocalFields< GFPT > &fld) const
Definition: grid_field_pack.hpp:149
Definition: gyro_radius.hpp:116
Definition: gyro_radius.hpp:24
Definition: gyro_radius.hpp:84
KOKKOS_INLINE_FUNCTION int get_node_index(int triangle_index, int tri_vertex_index) const
Definition: grid.tpp:158
Definition: field.hpp:23
Simd< double > phi
Definition: simd.hpp:152
Simd< double > z
Definition: simd.hpp:151
GFPT pack_type
Definition: grid_field_pack.hpp:25
GridFieldPack()
Definition: grid_field_pack.hpp:51
Definition: simd.hpp:139
KOKKOS_INLINE_FUNCTION void gather_all_fields(const PushControls &push_controls, int i_simd, int i_node, double wp, const FieldCorrection &corr, const SimdGridWeights< Order::One, PIT_GLOBAL > &grid_wts, const SimdGyroWeights< DriftKin > &rho_wts, LocalFields< GFPT > &fld) const
Definition: grid_field_pack.hpp:74
GridField< Device, VarType::Scalar, PIT_GLOBAL, TorType::OnePlane, KinType::DriftKin > loop_voltage
Definition: grid_field_pack.hpp:49
KOKKOS_INLINE_FUNCTION void phi_from_para(SimdVector &vec, int i_simd, const SimdVector &B, double Bmag) const
Definition: grid_field_pack.hpp:60
KOKKOS_INLINE_FUNCTION double magnitude(const int i_simd) const
Definition: simd.hpp:176
virtual ~GridFieldPackGeneric()=default
Definition: grid_field_pack.hpp:17
KOKKOS_INLINE_FUNCTION void gather_deltaf_conv(int i_simd, int node, int i_phi, double wp, const FieldCorrection &corr, const SimdGridWeights< Order::One, PIT_GLOBAL > &grid_wts, LocalFields< GFPT > &fld, const FieldT &labeled_field) const
Definition: grid_field_pack.hpp:66
pack_type pack
Definition: grid_field_pack.hpp:33