XGC1
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
field_aligned_local_fields.hpp
Go to the documentation of this file.
1 #ifndef FIELD_ALIGNED_LOCAL_FIELDS_HPP
2 #define FIELD_ALIGNED_LOCAL_FIELDS_HPP
3 
4 #include "globals.hpp"
5 #include "push_controls.hpp"
6 #include "grid.hpp"
7 #include "gyro_radius.hpp"
8 #include "electric_field.hpp"
9 #include "grid_field_pack.hpp"
10 
11 
12 /* This class manages how to obtain the local grid field data during the push. It either:
13  * - a) performs a gather operation when the method fields_at_point is used, or
14  * - b) performs a gather operation when initialized, stores those values and
15  * returns them when the fields_at_point method is used.
16  * Option b) is used in XGC1 for ions, because the ion location in an intermediate RK step
17  * might not have a valid associated field since it could be outside of the toroidal sector
18  * for which field information is present.
19  * The generic declaration uses Option a) */
20 template<KinType KT, PhiInterpType PIT>
22  // No members
23 
24  public:
25 
26  // Nothing to do in constructor
27  template<class Device>
28  KOKKOS_INLINE_FUNCTION FieldAlignedLocalFields(const Grid<Device> &grid, const PushControls &push_controls, const Species<Device> &species,
30 
31  /* Gather fields if calculating for electrons */
32  template<class Device>
33  KOKKOS_INLINE_FUNCTION void fields_at_point(const GridFieldPack<Device, PIT> fields, const PushControls &push_controls, const Grid<Device> &grid,
34  const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights<Order::One, PIT>& grid_wts,
35  SimdGyroWeights<KT>& rho_wts, LocalFields& fld) const{
36  fields.fields_at_point(push_controls, grid, B,gradpsi,grid_wts,rho_wts, fld);
37  }
38 };
39 
40 /* If using XGC1, then the class specialization on KinType GyroKin stores the field aligned
41  * local field calculated in the constructor, and converts it to cylindrical coordinates in
42  * fields_at_point. */
43 template<>
45  // Simd vectors are hard-coded to (r, z, phi), but here they represent field-aligned coordinates instead
46  // "phi" --> parallel to B-field
48 # ifdef EXPLICIT_EM
49  SimdVector dAh;
50  SimdVector dAs;
51 # endif
52 
53  // Members that don't get adjusted but need to be stored
54 #if defined(DELTAF_CONV) && !defined(EXPLICIT_EM)
55  SimdVector E00;
56  Simd<double> ddpotdt;
57 #endif
58 #ifdef EXPLICIT_EM
59  Simd<double> Ah;
60  Simd<double> As;
61  Simd<double> Epar_em;
62 #endif
63 
64  /* Convert stored electric/magnetic fields to field-aligned coordinates from cylindrical coordinates */
65  KOKKOS_INLINE_FUNCTION void convert_cyl_to_mag(const SimdVector &B, const SimdVector2D &gradpsi, const LocalFields& fld){
66  for (int i_simd = 0; i_simd<SIMD_SIZE; i_simd++){
67  double Bmag = B.magnitude(i_simd);
68  double Bp = B.poloidal_magnitude(i_simd);
69  double grad2 = gradpsi.r[i_simd]*gradpsi.r[i_simd] + gradpsi.z[i_simd]*gradpsi.z[i_simd];
70 
71  // get E-field in magnetic field-aligned basis
72  E.r[i_simd] = (fld.E.r[i_simd] *gradpsi.r[i_simd] + fld.E.z[i_simd] *gradpsi.z[i_simd])/grad2;
73  E.z[i_simd] = (fld.E.r[i_simd] *B.r[i_simd] + fld.E.z[i_simd] *B.z[i_simd] )/Bp;
74  E.phi[i_simd] = (E.z[i_simd]*Bp + fld.E.phi[i_simd]*B.phi[i_simd] )/Bmag; // parallel field
75 
76 # ifdef EXPLICIT_EM
77  // grad(A_h)
78  dAh.r[i_simd] = (fld.dAh.r[i_simd]*gradpsi.r[i_simd] + fld.dAh.z[i_simd]*gradpsi.z[i_simd])/grad2;
79  dAh.z[i_simd] = (fld.dAh.r[i_simd]*B.r[i_simd] + fld.dAh.z[i_simd]*B.z[i_simd] )/Bp;
80  dAh.phi[i_simd] = (dAh.z[i_simd] *Bp + fld.dAh.phi[i_simd]*B.phi[i_simd] )/Bmag;
81 
82  // These variables are for the mixed-variables formulation.
83  // Make sure that they are initialized to 0 if Hamiltonian
84  // formulation is used.
85  dAs.r[i_simd] = (fld.dAs.r[i_simd]*gradpsi.r[i_simd] + fld.dAs.z[i_simd]*gradpsi.z[i_simd])/grad2;
86  dAs.z[i_simd] = (fld.dAs.r[i_simd]*B.r[i_simd] + fld.dAs.z[i_simd]*B.z[i_simd] )/Bp;
87  dAs.phi[i_simd] = (dAs.z[i_simd] *Bp + fld.dAs.phi[i_simd]*B.phi[i_simd] )/Bmag;
88 # endif
89  }
90  }
91 
92  /* Convert stored electric/magnetic fields from field-aligned coordinates to cylindrical coordinates */
93  KOKKOS_INLINE_FUNCTION void convert_mag_to_cyl(const SimdVector &B, const SimdVector2D &gradpsi, LocalFields& fld) const{
94  for (int i_simd = 0; i_simd<SIMD_SIZE; i_simd++){
95  double Bmag = B.magnitude(i_simd);
96  double Bp = B.poloidal_magnitude(i_simd);
97  double Br_norm = B.r[i_simd] / Bp;
98  double Bz_norm = B.z[i_simd] / Bp;
99 
100  fld.E.r[i_simd] = E.r[i_simd]*gradpsi.r[i_simd]+E.z[i_simd]*Br_norm;
101  fld.E.z[i_simd] = E.r[i_simd]*gradpsi.z[i_simd]+E.z[i_simd]*Bz_norm;
102  fld.E.phi[i_simd] = (E.phi[i_simd]*Bmag - fld.E.r[i_simd]*B.r[i_simd] - fld.E.z[i_simd]*B.z[i_simd])/B.phi[i_simd];
103 
104 # ifdef EXPLICIT_EM
105  fld.dAh.r[i_simd] = dAh.r[i_simd]*gradpsi.r[i_simd]+dAh.z[i_simd]*Br_norm;
106  fld.dAh.z[i_simd] = dAh.r[i_simd]*gradpsi.z[i_simd]+dAh.z[i_simd]*Bz_norm;
107  fld.dAh.phi[i_simd] = (dAh.phi[i_simd]*Bmag - fld.dAh.r[i_simd]*B.r[i_simd] - fld.dAh.z[i_simd]*B.z[i_simd])/B.phi[i_simd];
108 
109  // The following are for the mixed-variable formulation
110  // Make sure that no all variables are initialized
111  fld.dAs.r[i_simd] = dAs.r[i_simd]*gradpsi.r[i_simd]+dAs.z[i_simd]*Br_norm;
112  fld.dAs.z[i_simd] = dAs.r[i_simd]*gradpsi.z[i_simd]+dAs.z[i_simd]*Bz_norm;
113  fld.dAs.phi[i_simd] = (dAs.phi[i_simd]*Bmag - fld.dAs.r[i_simd]*B.r[i_simd] - fld.dAs.z[i_simd]*B.z[i_simd])/B.phi[i_simd];
114 # endif
115  }
116  }
117 
118  public:
119 
120  /* Gather electric/magnetic fields and convert to field-aligned coordinates for storage */
121  template<class Device>
122  KOKKOS_INLINE_FUNCTION FieldAlignedLocalFields(const Grid<Device> &grid, const PushControls &push_controls, const Species<Device> &species,
124  // Declare arrays returned from B-field
125  SimdVector B;
126  SimdVector jacb[3];
127  Simd<double> psi;
128  SimdVector2D gradpsi;
129  SimdVector tdb;
130  Simd<bool> rz_outside;
131 
132  // obtain B-field information
133  magnetic_field.field(part.ph.x(), B, jacb, psi, gradpsi, tdb, rz_outside);
134  part.deactivate(rz_outside);
135 
136  // Get the current field following grid triangle and p weights
137  grid.get_grid_weights(magnetic_field, part.ph.v(), psi, grid_wts);
138 
139  // calculate rho weights
140  SimdGyroWeights<GyroKin> rho_wts(grid, magnetic_field, species, part.ph.v(), part.ct.mu);
141 
142  // get electric field
143  LocalFields fld;
144  gfpack.fields_at_point(push_controls, grid, B,gradpsi,grid_wts,rho_wts, fld);
145 
146  // Convert to field-aligned coordinates and store
147  convert_cyl_to_mag(B, gradpsi, fld);
148 
149  // These get stored without any conversion; Consolidate field_aligned_local_fields with local_fields objects to avoid this copy
150 #if defined(DELTAF_CONV) && !defined(EXPLICIT_EM)
151  E00.r=fld.E00.r; // Should E00 also be field-aligned?
152  E00.z=fld.E00.z;
153  E00.phi=fld.E00.phi;
154  ddpotdt = fld.ddpotdt;
155 #endif
156 #ifdef EXPLICIT_EM
157  Ah = fld.Ah;
158  As = fld.As;
159  Epar_em = fld.Epar_em;
160 #endif
161  }
162 
163  /* Return the stored field values for XGC1 ions */
164  template<class Device>
165  KOKKOS_INLINE_FUNCTION void fields_at_point(const GridFieldPack<Device, PhiInterpType::Planes> fields, const PushControls &push_controls, const Grid<Device> &grid,
166  const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights<Order::One, PhiInterpType::Planes>& grid_wts,
167  SimdGyroWeights<GyroKin>& rho_wts, LocalFields& fld) const{
168  convert_mag_to_cyl(B, gradpsi, fld);
169 
170  // These get stored without any conversion; Consolidate field_aligned_local_fields with local_fields objects to avoid this copy
171 #if defined(DELTAF_CONV) && !defined(EXPLICIT_EM)
172  fld.E00.r=E00.r; // Should E00 also be field-aligned?
173  fld.E00.z=E00.z;
174  fld.E00.phi=E00.phi;
175  fld.ddpotdt = ddpotdt;
176 #endif
177 #ifdef EXPLICIT_EM
178  fld.Ah = Ah;
179  fld.As = As;
180  fld.Epar_em = Epar_em;
181 #endif
182  }
183 };
184 
185 
186 #endif
Simd< double > r
Definition: simd.hpp:150
Definition: simd.hpp:149
Definition: grid_weights.hpp:73
Definition: grid_weights.hpp:47
Definition: push_controls.hpp:9
Definition: magnetic_field.hpp:12
Definition: grid.hpp:67
Definition: grid_field_pack.hpp:19
Definition: local_fields.hpp:7
KOKKOS_INLINE_FUNCTION double poloidal_magnitude(const int i_simd) const
Definition: simd.hpp:181
KOKKOS_INLINE_FUNCTION FieldAlignedLocalFields(const Grid< Device > &grid, const PushControls &push_controls, const Species< Device > &species, const MagneticField< Device > &magnetic_field, const GridFieldPack< Device, PIT > &gfpack, SimdParticles &part, SimdGridWeights< Order::One, PIT > &grid_wts)
Definition: field_aligned_local_fields.hpp:28
PhiInterpType
Definition: globals.hpp:95
Definition: gyro_radius.hpp:24
Definition: gyro_radius.hpp:84
KOKKOS_INLINE_FUNCTION SimdVector2D & x()
Definition: particles.hpp:27
KOKKOS_INLINE_FUNCTION SimdVector & v()
Definition: particles.hpp:39
KOKKOS_INLINE_FUNCTION void convert_cyl_to_mag(const SimdVector &B, const SimdVector2D &gradpsi, const LocalFields &fld)
Definition: field_aligned_local_fields.hpp:65
Simd< double > z
Definition: simd.hpp:141
Definition: globals.hpp:90
KOKKOS_INLINE_FUNCTION void fields_at_point(const GridFieldPack< Device, PhiInterpType::Planes > fields, const PushControls &push_controls, const Grid< Device > &grid, const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights< Order::One, PhiInterpType::Planes > &grid_wts, SimdGyroWeights< GyroKin > &rho_wts, LocalFields &fld) const
Definition: field_aligned_local_fields.hpp:165
SimdPhase ph
Definition: particles.hpp:59
Simd< double > phi
Definition: simd.hpp:152
Definition: particles.hpp:58
Simd< double > z
Definition: simd.hpp:151
KOKKOS_INLINE_FUNCTION FieldAlignedLocalFields(const Grid< Device > &grid, const PushControls &push_controls, const Species< Device > &species, const MagneticField< Device > &magnetic_field, const GridFieldPack< Device, PhiInterpType::Planes > &gfpack, SimdParticles &part, SimdGridWeights< Order::One, PhiInterpType::Planes > &grid_wts)
Definition: field_aligned_local_fields.hpp:122
Simd< double > r
Definition: simd.hpp:140
Definition: field_aligned_local_fields.hpp:21
Definition: magnetic_field.F90:1
SimdVector E
Definition: field_aligned_local_fields.hpp:47
Definition: simd.hpp:139
KOKKOS_INLINE_FUNCTION void field(const SimdVector2D &x, SimdVector &bvec, SimdVector(&jacb)[3], Simd< double > &psivec, SimdVector2D &gradpsi, SimdVector &tdb, Simd< bool > &rz_outside) const
Definition: magnetic_field.tpp:131
SimdConstants ct
Definition: particles.hpp:60
KOKKOS_INLINE_FUNCTION void convert_mag_to_cyl(const SimdVector &B, const SimdVector2D &gradpsi, LocalFields &fld) const
Definition: field_aligned_local_fields.hpp:93
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 &fld) const
Definition: grid_field_pack.hpp:174
Definition: species.hpp:75
KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField< Device > &magnetic_field, const SimdVector &v, const Simd< double > &psi, SimdVector2D &xff, SimdGridWeights< Order::One, PIT > &grid_wts) const
Definition: grid.tpp:705
KOKKOS_INLINE_FUNCTION double magnitude(const int i_simd) const
Definition: simd.hpp:176
Simd< double > mu
Definition: particles.hpp:52
KOKKOS_INLINE_FUNCTION void fields_at_point(const GridFieldPack< Device, PIT > fields, const PushControls &push_controls, const Grid< Device > &grid, const SimdVector &B, const SimdVector2D &gradpsi, const SimdGridWeights< Order::One, PIT > &grid_wts, SimdGyroWeights< KT > &rho_wts, LocalFields &fld) const
Definition: field_aligned_local_fields.hpp:33
SimdVector E
Definition: local_fields.hpp:8
KOKKOS_INLINE_FUNCTION void deactivate(const Simd< bool > &mask)
Definition: particles.hpp:76