XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
grid.hpp
Go to the documentation of this file.
1 #ifndef GRID_HPP
2 #define GRID_HPP
3 
4 #include "magnetic_field.hpp"
5 #include "grid_structs.hpp"
6 #include "grid_weights.hpp"
7 #include "grid_files.hpp"
8 #include "guess_table.hpp"
9 #include "uniform_range.hpp"
10 #include "guess_list_1d.hpp"
11 #include "gradient_matrices.hpp"
12 
13 extern "C" void init_fortran_vols_and_areas(double* tr_area_h, double* tr_vol_h,
14 #ifdef XGC1
15  double* node_vol_ff_h,
16 #endif
17  double* node_vol_nearest_h, double* node_vol_h, double* node_area_h);
18 
19 // Volume- and Area-related Host Views. Should find a different location for these so that these unused views aren't sent to GPU
21  View<double*,CLayout,HostType> tr_area_h;
22  View<double*,CLayout,HostType> tr_vol_h;
23  View<double**,CLayout,HostType> node_vol_ff_h;
24  View<double*,CLayout,HostType> node_vol_nearest_h;
25  View<double*,CLayout,HostType> node_vol_h;
26  View<double*,CLayout,HostType> node_area_h;
27 
29 
30  VolumesAndAreas(int ntriangle_in, int nnode_in)
31  : tr_area_h(NoInit("tr_area"), ntriangle_in),
32  tr_vol_h(NoInit("tr_vol"), ntriangle_in),
33 #ifdef XGC1
34  node_vol_ff_h(NoInit("node_vol_ff"),2, nnode_in),
35 #endif
36  node_vol_nearest_h(NoInit("node_vol_nearest"),nnode_in),
37  node_vol_h(NoInit("node_vol"),nnode_in),
38  node_area_h(NoInit("node_area"),nnode_in)
39  {}
40 
42 #ifndef NO_FORTRAN_MODULES
44 #ifdef XGC1
45  node_vol_ff_h.data(),
46 #endif
47  node_vol_nearest_h.data(), node_vol_h.data(), node_area_h.data());
48 #endif
49  }
50 };
51 
52 template<class Device>
53 struct Projection{
54  View<int**,CLayout,Device> tr;
55  View<double***,CLayout,Device> p;
56  View<double**,CLayout,Device> dx;
57 
59  Projection(int nnode)
60  : tr(NoInit("tr"),2,nnode),
61  p(NoInit("p"),2,nnode,3),
62  dx(NoInit("dx"),2,nnode){}
63 };
64 
65 // Grid class
66 template<class Device>
67 class Grid {
68  // This line enables private access between the host and device implementations of this class
69  // i.e. Makes the copy from host to device more straightforward
70  template<class Device2> friend class Grid;
71 
72  public:
73 
74  // Constructor from namelist
76 
77  // Default constructor
78  Grid(){}
79 
80  // Create a mirror with a different device type
81  template<class Device2>
82  inline Grid<Device2> mirror() const{
83  Grid<Device2> m;
84 
85  m.ntriangle = ntriangle;
86  m.nnode = nnode;
87  m.nplanes = nplanes;
89  m.delta_phi = delta_phi;
91 
92  m.psi = my_mirror_view(psi, Device2());
93  mirror_copy(m.psi, psi);
94  m.bfield = my_mirror_view(bfield, Device2());
96  m.basis = my_mirror_view(basis, Device2());
98  m.nodes = my_mirror_view(nodes, Device2());
100  m.mapping = my_mirror_view(mapping, Device2());
102 
103  m.guess = guess.template mirror<Device2>();
104 
108  m.psi_surf2 = my_mirror_view(psi_surf2, Device2());
110 
111  m.psi00 = psi00;
112  m.psi_guess = psi_guess.template mirror<Device2>();
113 
114  m.rgn = my_mirror_view(rgn, Device2());
115  mirror_copy(m.rgn, rgn);
116  m.gx = my_mirror_view(gx, Device2());
117  mirror_copy(m.gx, gx);
118 
119  m.nwall = nwall;
120  m.wall_nodes = my_mirror_view(wall_nodes, Device2());
122  m.node_to_wall = my_mirror_view(node_to_wall, Device2());
124 
126 
127  m.volumes_and_areas = volumes_and_areas; // Doesn't require mirror view since everything is on host
128 
131 
132 #ifdef NO_FORTRAN_MODULES
133  m.gradient_matrices_h = gradient_matrices_h;
134 #endif
135 
136  return m;
137  }
138 
139  KOKKOS_INLINE_FUNCTION void t_coeff(const SimdVector2D &x, const Simd<int>& itr, SimdGridVec &p) const;
140 
141  KOKKOS_INLINE_FUNCTION void t_coeff_mod(const MagneticField<Device> &magnetic_field, const double r, const double z, const double psiin, const int itr, double (&p)[3] ) const;
142 
143  KOKKOS_INLINE_FUNCTION void t_coeff_mod(const MagneticField<Device> &magnetic_field, const SimdVector2D &xy, const Simd<double>& psiin, const Simd<int>& itr, SimdGridVec &p ) const;
144 
145  KOKKOS_INLINE_FUNCTION void search_tr2_no_precheck(const double r, const double z, int& itr, double (&p)[3]) const;
146 
147  KOKKOS_INLINE_FUNCTION void search_tr2( const SimdVector2D &xy, Simd<int>& itr, SimdGridVec &pout ) const;
148 
149  KOKKOS_INLINE_FUNCTION void search_tr_check_guess(const SimdVector2D &x, const Simd<int>& old_itr,Simd<int>& itr,SimdGridVec &p ) const;
150 
151  KOKKOS_INLINE_FUNCTION void psi_search(double psi,double &wp,int &ip) const;
152 
153  KOKKOS_INLINE_FUNCTION void follow_field_to_midplane(const MagneticField<Device> &magnetic_field, const SimdVector2D &x, const Simd<double>& phi, SimdVector2D &xff) const;
154 
155  KOKKOS_INLINE_FUNCTION void get_grid_weights_ff(const MagneticField<Device> &magnetic_field, const SimdVector &v, const Simd<double>& psi_in,SimdVector2D &xff,SimdGridWeights<Order::One, PhiInterpType::Planes>& grid_wts) const;
156 
157  KOKKOS_INLINE_FUNCTION void get_grid_weights_no_ff(const MagneticField<Device> &magnetic_field, const SimdVector2D &x, const Simd<double>& psi_in,SimdGridWeights<Order::One, PhiInterpType::None>& grid_wts) const;
158 
159  template<PhiInterpType PIT>
160  KOKKOS_INLINE_FUNCTION void nearest_node(const SimdGridWeights<Order::One, PIT>& grid_wts, SimdGridWeights<Order::Zero, PIT>& grid_wts0) const;
161 
162  template<PhiInterpType PIT>
163  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;
164 
165  template<PhiInterpType PIT>
166  KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField<Device> &magnetic_field, const SimdVector &v, SimdVector2D &xff, SimdGridWeights<Order::One, PIT>& grid_wts) const;
167 
168  template<PhiInterpType PIT>
169  KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField<Device> &magnetic_field, const SimdVector &v, const Simd<double>& psi_in, SimdGridWeights<Order::Zero, PIT>& grid_wts) const;
170 
171  template<PhiInterpType PIT>
172  KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField<Device> &magnetic_field, const SimdVector &v, const Simd<double>& psi_in, SimdGridWeights<Order::One, PIT>& grid_wts) const;
173 
174  template<PhiInterpType PIT>
175  KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField<Device> &magnetic_field, const SimdVector &v, SimdGridWeights<Order::Zero, PIT>& grid_wts) const;
176 
177  template<PhiInterpType PIT>
178  KOKKOS_INLINE_FUNCTION void get_grid_weights(const MagneticField<Device> &magnetic_field, const SimdVector &v, SimdGridWeights<Order::One, PIT>& grid_wts) const;
179 
180  KOKKOS_INLINE_FUNCTION void get_wall_index(const Simd<bool>& just_left_the_grid,const SimdVector2D &x, const SimdGridWeights<Order::One, PIT_GLOBAL>& grid_wts, Simd<int>& widx) const;
181 
182  KOKKOS_INLINE_FUNCTION void wedge_modulo_phi(Simd<double>& phi_mod) const;
183 
184  KOKKOS_INLINE_FUNCTION double wedge_modulo_phi(double phi) const;
185 
186  KOKKOS_INLINE_FUNCTION bool node_is_inside_psi_range(const MagneticField<Device> &magnetic_field, const int node) const;
187 
188  KOKKOS_INLINE_FUNCTION int get_plane_index(double phi) const;
189 
190  KOKKOS_INLINE_FUNCTION int get_node_index(int triangle_index, int tri_vertex_index) const;
191 
192  KOKKOS_INLINE_FUNCTION double node_area_to_volume(const MagneticField<Device> &magnetic_field, double area, int node_index) const;
193 
194  KOKKOS_INLINE_FUNCTION void get_triangle_area_and_volume(const MagneticField<Device> &magnetic_field, int i, double& area, double& volume) const;
195 
196  KOKKOS_INLINE_FUNCTION double get_r_center_of_mass(int itr) const;
197 
198  KOKKOS_INLINE_FUNCTION bool node_is_in_region_1_or_2_no_wall(const int inode) const;
199 
200  KOKKOS_INLINE_FUNCTION bool node_is_in_private_region_no_wall(const int inode) const;
201 
202  KOKKOS_INLINE_FUNCTION bool node_is_in_included_region(const int inode, const bool exclude_private_region) const;
203 
204  KOKKOS_INLINE_FUNCTION bool node_is_in_region_1_or_2(const MagneticField<Device> &magnetic_field, const int inode) const;
205 
206  KOKKOS_INLINE_FUNCTION int uses_rz_basis(const int inode) const;
207 
208  KOKKOS_INLINE_FUNCTION double get_r(const int inode) const;
209 
210  KOKKOS_INLINE_FUNCTION double get_dist2_from_node(const int inode, const double r, const double z) const;
211 
212  KOKKOS_INLINE_FUNCTION double get_dist_from_node(const int inode, const double r, const double z) const;
213 
214  KOKKOS_INLINE_FUNCTION void get_rz_coordinates(const int inode, double& r, double& z) const;
215 
216  KOKKOS_INLINE_FUNCTION void get_rz_coordinates(const Simd<int>& grid_inds, SimdVector2D& x) const;
217 
218  KOKKOS_INLINE_FUNCTION void get_rz_coordinates(const int itr, const double (&p)[3], double& r, double& z) const;
219 
220  KOKKOS_INLINE_FUNCTION void get_rz_coordinates(const Simd<int>& itr, const SimdGridVec& p, SimdVector2D& x) const;
221 
222  KOKKOS_INLINE_FUNCTION int get_nearest_node(const Simd<int>& itr, const SimdGridVec& p, int i_simd) const;
223 
224  KOKKOS_INLINE_FUNCTION int get_nearest_node(const int itr, const double (&p)[3]) const;
225 
226  KOKKOS_INLINE_FUNCTION double get_nearest_midplane(double phi) const;
227 
228  KOKKOS_INLINE_FUNCTION void set_gradient_mat_triangle(const View<int*,CLayout,DeviceType>& num_t_node,
229  const View<int**,CLayout,DeviceType>& tr_node, const View<double*,CLayout,DeviceType>& tr_area,
230  const View<double***, CLayout, DeviceType>& unit_vecs,
231  const int i, const GradientMatrices<Device>& gradient_matrices) const;
232 
233  KOKKOS_INLINE_FUNCTION int get_characteristic_length(const View<int*,CLayout,DeviceType>& num_t_node,
234  const View<int**,CLayout,DeviceType>& tr_node, const int i,
235  double& dist_triangle, double& dist_psi, double& dist_theta) const;
236 
237  KOKKOS_INLINE_FUNCTION void set_grad_matrix_from_psi_theta(bool is_psi_dir, int itr_pos, int itr_neg, double (&p_pos)[3], double (&p_neg)[3], double dl_pos, double dl_neg, int inode, const Matrix<Device>& matrix) const;
238 
239  void draw_ascii_grid(MagneticField<Device> magnetic_field, double rmin, double rmax, double dr, double zmin, double zmax, double dz){
240  printf("\n\n******** PLOT OF GRID ******* \n");
241  Kokkos::View<int*,HostType> save_vert("save_vert", 1+((rmax-rmin)/dr));
242  for (double zplot = zmax; zplot>=zmin; zplot-=dz){
243  printf("\n");
244  int save_num = 0;
245  int rsave_ind = -1;
246  for (double rplot = rmin; rplot<=rmax; rplot+=dr){
247  rsave_ind++;
248  SimdVector v;
249  v.r = rplot;
250  v.z = zplot;
251  v.phi = 0.001;
252 
253  // Get itr (triangle index)
254  SimdGridVec p;
255  Simd<int> itr = -1;
256  get_grid_weights(magnetic_field, v, itr, p);
257  if(itr[0]!=save_num || itr[0]!=save_vert[rsave_ind]){
258  if(itr[0]<10 && itr[0]>=0) printf("%d ",itr[0]);
259  else printf("%d ",itr[0]);
260  save_num = itr[0];
261  save_vert[rsave_ind] = itr[0];
262  } else {printf(" ");}
263  }
264  }
265  }
266 
267  int ntriangle;
268  int nnode;
269  int nplanes;
270  double wedge_angle;
271  double delta_phi;
272  double inv_delta_phi;
273 
274  Kokkos::View<double*,Kokkos::LayoutRight,Device> psi;
275  Kokkos::View<double*,Kokkos::LayoutRight,Device> bfield;
276 
277  private:
278  Kokkos::View<int*, Kokkos::LayoutRight,Device> basis;
279  Kokkos::View<Vertex*,Kokkos::LayoutRight,Device> nodes;
280  Kokkos::View<VertexMap*,Kokkos::LayoutRight,Device> mapping;
282 
283  public:
287  Kokkos::View<double*,Kokkos::LayoutRight,Device> psi_surf2;
289 
290  private:
292 
293  Kokkos::View<int*,Kokkos::LayoutRight,Device> rgn;
294  Kokkos::View<RZPair*,Kokkos::LayoutRight,Device> gx;
295 
296  public:
297  int nwall;
298 
299  private:
300  Kokkos::View<int*,Kokkos::LayoutRight,Device> wall_nodes;
301  Kokkos::View<int*,Kokkos::LayoutRight,Device> node_to_wall;
302 
304 
305  public:
306 
307  // Volume- and Area-related Host Views. Should find a different location for these so that these unused views aren't sent to GPU
309 
310  // Projections - reside in host memory for now
313 
314 
315 #ifdef NO_FORTRAN_MODULES
316  GradientMatrices<HostType> gradient_matrices_h;
317 #endif
318 };
319 
320 #include "grid.tpp"
321 #endif
Projection< HostType > one_plane_ff
Definition: grid.hpp:312
KOKKOS_INLINE_FUNCTION void t_coeff(const SimdVector2D &x, const Simd< int > &itr, SimdGridVec &p) const
Definition: grid.tpp:494
double inv_delta_phi
1/delta_phi
Definition: grid.hpp:272
VolumesAndAreas volumes_and_areas
Definition: grid.hpp:308
Simd< double > r
Definition: simd.hpp:150
Definition: guess_list_1d.hpp:7
KOKKOS_INLINE_FUNCTION void get_grid_weights_ff(const MagneticField< Device > &magnetic_field, const SimdVector &v, const Simd< double > &psi_in, SimdVector2D &xff, SimdGridWeights< Order::One, PhiInterpType::Planes > &grid_wts) const
Definition: grid.tpp:658
Definition: simd.hpp:149
GuessList1D< Device > psi_guess
Definition: grid.hpp:291
Definition: grid.hpp:53
Projection()
Definition: grid.hpp:58
KOKKOS_INLINE_FUNCTION double get_r_center_of_mass(int itr) const
Definition: grid.tpp:942
Grid()
Definition: grid.hpp:78
Definition: gradient_matrices.hpp:9
void mirror_copy(T1 &view_dest, const T2 &view_src)
Definition: my_mirror_view.hpp:122
void set_fortran_ptrs()
Definition: grid.hpp:41
Definition: grid_weights.hpp:77
Kokkos::View< double *, Kokkos::LayoutRight, Device > bfield
Magnetic field magnitude at nodes.
Definition: grid.hpp:275
int ntriangle
Number of grid triangles.
Definition: grid.hpp:267
KOKKOS_INLINE_FUNCTION double get_dist_from_node(const int inode, const double r, const double z) const
Definition: grid.tpp:1049
double minval_psi_surf2
Definition: grid.hpp:285
KOKKOS_INLINE_FUNCTION bool node_is_in_private_region_no_wall(const int inode) const
Definition: grid.tpp:969
Definition: grid_weights.hpp:51
KOKKOS_INLINE_FUNCTION int get_characteristic_length(const View< int *, CLayout, DeviceType > &num_t_node, const View< int **, CLayout, DeviceType > &tr_node, const int i, double &dist_triangle, double &dist_psi, double &dist_theta) const
Definition: grid.tpp:1236
Definition: NamelistReader.hpp:193
Definition: magnetic_field.hpp:12
Grid< Device2 > mirror() const
Definition: grid.hpp:82
KOKKOS_INLINE_FUNCTION bool node_is_in_included_region(const int inode, const bool exclude_private_region) const
Definition: grid.tpp:982
View< double *, CLayout, HostType > node_area_h
Definition: grid.hpp:26
Definition: grid_files.hpp:17
Definition: grid.hpp:67
KOKKOS_INLINE_FUNCTION void get_triangle_area_and_volume(const MagneticField< Device > &magnetic_field, int i, double &area, double &volume) const
Definition: grid.tpp:916
KOKKOS_INLINE_FUNCTION void t_coeff_mod(const MagneticField< Device > &magnetic_field, const double r, const double z, const double psiin, const int itr, double(&p)[3]) const
Definition: grid.tpp:373
int nplanes
Number of planes.
Definition: grid.hpp:269
KOKKOS_INLINE_FUNCTION void get_wall_index(const Simd< bool > &just_left_the_grid, const SimdVector2D &x, const SimdGridWeights< Order::One, PIT_GLOBAL > &grid_wts, Simd< int > &widx) const
Definition: grid.tpp:797
Projection< HostType > half_plane_ff
Definition: grid.hpp:311
KOKKOS_INLINE_FUNCTION double get_nearest_midplane(double phi) const
Definition: grid.tpp:1159
KOKKOS_INLINE_FUNCTION int uses_rz_basis(const int inode) const
Definition: grid.tpp:1008
Kokkos::View< double *, Kokkos::LayoutRight, Device > psi_surf2
Definition: grid.hpp:287
double wedge_angle
The size of the wedge (the model is periodic in phi, a angle of e.g. pi means half the tokamak is mod...
Definition: grid.hpp:270
int nwall
Definition: grid.hpp:297
Definition: grid.hpp:20
View< double ***, CLayout, Device > p
Definition: grid.hpp:55
KOKKOS_INLINE_FUNCTION void search_tr_check_guess(const SimdVector2D &x, const Simd< int > &old_itr, Simd< int > &itr, SimdGridVec &p) const
Definition: grid.tpp:512
void init_fortran_vols_and_areas(double *tr_area_h, double *tr_vol_h, double *node_vol_nearest_h, double *node_vol_h, double *node_area_h)
Kokkos::View< double *, Kokkos::LayoutRight, Device > psi
An array of psi coordinates.
Definition: grid.hpp:274
KOKKOS_INLINE_FUNCTION void get_rz_coordinates(const int inode, double &r, double &z) const
Definition: grid.tpp:1074
Definition: guess_table.hpp:8
GuessTable< Device > guess
Definition: grid.hpp:281
Kokkos::View< int *, Kokkos::LayoutRight, Device > rgn
Definition: grid.hpp:293
View< double *, CLayout, HostType > tr_area_h
Definition: grid.hpp:21
Kokkos::View< int *, Kokkos::LayoutRight, Device > basis
0 for Psi-theta, 1 for R-Z ; should be enum (and bool/char?)
Definition: grid.hpp:278
KOKKOS_INLINE_FUNCTION void get_grid_weights_no_ff(const MagneticField< Device > &magnetic_field, const SimdVector2D &x, const Simd< double > &psi_in, SimdGridWeights< Order::One, PhiInterpType::None > &grid_wts) const
Definition: grid.tpp:637
int npsi_surf2
Definition: grid.hpp:284
KOKKOS_INLINE_FUNCTION int get_nearest_node(const Simd< int > &itr, const SimdGridVec &p, int i_simd) const
Definition: grid.tpp:1131
KOKKOS_INLINE_FUNCTION int get_node_index(int triangle_index, int tri_vertex_index) const
Definition: grid.tpp:885
KOKKOS_INLINE_FUNCTION double get_r(const int inode) const
Definition: grid.tpp:1020
KOKKOS_INLINE_FUNCTION bool node_is_inside_psi_range(const MagneticField< Device > &magnetic_field, const int node) const
Definition: grid.tpp:859
VolumesAndAreas(int ntriangle_in, int nnode_in)
Definition: grid.hpp:30
KOKKOS_INLINE_FUNCTION void nearest_node(const SimdGridWeights< Order::One, PIT > &grid_wts, SimdGridWeights< Order::Zero, PIT > &grid_wts0) const
Definition: grid.tpp:729
KOKKOS_INLINE_FUNCTION int get_plane_index(double phi) const
Definition: grid.tpp:872
KOKKOS_INLINE_FUNCTION void search_tr2_no_precheck(const double r, const double z, int &itr, double(&p)[3]) const
Definition: grid.tpp:247
KOKKOS_INLINE_FUNCTION void wedge_modulo_phi(Simd< double > &phi_mod) const
Definition: grid.tpp:767
Simd< double > phi
Definition: simd.hpp:152
View< double *, CLayout, HostType > tr_vol_h
Definition: grid.hpp:22
KOKKOS_INLINE_FUNCTION void search_tr2(const SimdVector2D &xy, Simd< int > &itr, SimdGridVec &pout) const
Definition: grid.tpp:298
Definition: matrix.hpp:13
KOKKOS_INLINE_FUNCTION void follow_field_to_midplane(const MagneticField< Device > &magnetic_field, const SimdVector2D &x, const Simd< double > &phi, SimdVector2D &xff) const
Definition: grid.tpp:603
Simd< double > z
Definition: simd.hpp:151
KOKKOS_INLINE_FUNCTION double node_area_to_volume(const MagneticField< Device > &magnetic_field, double area, int node_index) const
Definition: grid.tpp:901
Definition: grid_structs.hpp:7
KOKKOS_INLINE_FUNCTION bool node_is_in_region_1_or_2_no_wall(const int inode) const
Definition: grid.tpp:957
View< T *, CLayout, Device > my_mirror_view(const View< T *, CLayout, Device > &view, Device nd)
Definition: my_mirror_view.hpp:14
Kokkos::View< int *, Kokkos::LayoutRight, Device > node_to_wall
Definition: grid.hpp:301
Kokkos::View< int *, Kokkos::LayoutRight, Device > wall_nodes
Definition: grid.hpp:300
KOKKOS_INLINE_FUNCTION double get_dist2_from_node(const int inode, const double r, const double z) const
Definition: grid.tpp:1035
Definition: grid_weights.hpp:56
View< double *, CLayout, HostType > node_vol_h
Definition: grid.hpp:25
UniformRange psi00
Definition: grid.hpp:288
View< int **, CLayout, Device > tr
Definition: grid.hpp:54
Projection(int nnode)
Definition: grid.hpp:59
Definition: magnetic_field.F90:1
void draw_ascii_grid(MagneticField< Device > magnetic_field, double rmin, double rmax, double dr, double zmin, double zmax, double dz)
Definition: grid.hpp:239
Definition: simd.hpp:139
Definition: uniform_range.hpp:6
KOKKOS_INLINE_FUNCTION void set_gradient_mat_triangle(const View< int *, CLayout, DeviceType > &num_t_node, const View< int **, CLayout, DeviceType > &tr_node, const View< double *, CLayout, DeviceType > &tr_area, const View< double ***, CLayout, DeviceType > &unit_vecs, const int i, const GradientMatrices< Device > &gradient_matrices) const
Definition: grid.tpp:1174
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:682
double maxval_psi_surf2
Definition: grid.hpp:286
int nnode
Number of grid nodes.
Definition: grid.hpp:268
double delta_phi
Distance between planes.
Definition: grid.hpp:271
double eps_flux_surface
The minimum difference in psi that counts as being on two different flux surfaces.
Definition: grid.hpp:303
Kokkos::View< VertexMap *, Kokkos::LayoutRight, Device > mapping
Definition: grid.hpp:280
KOKKOS_INLINE_FUNCTION void psi_search(double psi, double &wp, int &ip) const
Definition: grid.tpp:556
View< double **, CLayout, Device > dx
Definition: grid.hpp:56
KOKKOS_INLINE_FUNCTION void set_grad_matrix_from_psi_theta(bool is_psi_dir, int itr_pos, int itr_neg, double(&p_pos)[3], double(&p_neg)[3], double dl_pos, double dl_neg, int inode, const Matrix< Device > &matrix) const
Definition: grid.tpp:1361
Kokkos::ViewAllocateWithoutInitializing NoInit
Definition: space_settings.hpp:68
View< double *, CLayout, HostType > node_vol_nearest_h
Definition: grid.hpp:24
VolumesAndAreas()
Definition: grid.hpp:28
KOKKOS_INLINE_FUNCTION bool node_is_in_region_1_or_2(const MagneticField< Device > &magnetic_field, const int inode) const
Definition: grid.tpp:996
Kokkos::View< Vertex *, Kokkos::LayoutRight, Device > nodes
Definition: grid.hpp:279
Kokkos::View< RZPair *, Kokkos::LayoutRight, Device > gx
Definition: grid.hpp:294
View< double **, CLayout, HostType > node_vol_ff_h
Definition: grid.hpp:23