XGCa
|
#include <grid.hpp>
Public Member Functions | |
Grid (const MagneticField< Device > &magnetic_field, int ntriangle, int nnode, int nplanes, double delta_phi, int nguess, int guess_list_len, double guess_min1, double guess_min2, double guess_max1, double guess_max2, double inv_guess_d1, double inv_guess_d2, int npsi_surf2, int iphi_offset, int npsi00, double psi00min, double psi00max, double dpsi00, int nwall, double wedge_angle, double eq_x_psi, double *psi, int *basis, Vertex *nodes, VertexMap *mapping, int *guess_list, int *guess_xtable, int *guess_count, int *guess_list_1d, double *psi_surf2, double minval_psi_surf2, double maxval_psi_surf2, int *wall_nodes, int *node_to_wall, int *rgn, RZPair *gx, int nrho, double rhomax, double phimax) | |
Grid (GridFiles &grid_files, const MagneticField< Device > &magnetic_field, int nplanes_in, int wedge_n_in, int iphi_offset_in, int guess_table_size, int npsi00_in, int nrho_in, double rhomax_in) | |
Grid (NLReader::NamelistReader &nlr, const GridFiles &grid_files) | |
Grid () | |
KOKKOS_INLINE_FUNCTION void | t_coeff (const SimdVector2D &x, const Simd< int > &itr, SimdGridVec &p) const |
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 |
KOKKOS_INLINE_FUNCTION void | search_tr2 (const SimdVector2D &xy, Simd< int > &itr, SimdGridVec &pout) const |
KOKKOS_INLINE_FUNCTION void | search_tr_check_guess (const SimdVector2D &x, const Simd< int > &old_itr, Simd< int > &itr, SimdGridVec &p) const |
KOKKOS_INLINE_FUNCTION void | search_ptl_1d (double psi, double &wp, int &ip) const |
KOKKOS_INLINE_FUNCTION void | charge_search_index (const MagneticField< Device > &magnetic_field, const SimdVector2D &x, const Simd< double > &phi, const Simd< double > &psi, SimdVector2D &xff, Simd< int > &itr, SimdGridVec &p) const |
KOKKOS_INLINE_FUNCTION void | charge_search_index (const MagneticField< Device > &magnetic_field, const SimdVector &v, Simd< int > &itr, SimdGridVec &p) const |
KOKKOS_INLINE_FUNCTION void | charge_search_index (const MagneticField< Device > &magnetic_field, const SimdVector &v, const Simd< double > &psi, Simd< int > &itr, SimdGridVec &p) const |
KOKKOS_INLINE_FUNCTION void | get_wall_index (const Simd< bool > &just_left_the_grid, const SimdVector2D &x, const Simd< int > &itr, const SimdGridVec &p, Simd< int > &widx) const |
KOKKOS_INLINE_FUNCTION void | wedge_modulo_phi (Simd< double > &phi_mod) const |
KOKKOS_INLINE_FUNCTION void | check_triangle (const Simd< int > &itr, Simd< bool > &wasnt_in_triangle) const |
KOKKOS_INLINE_FUNCTION bool | node_is_inside_psi_range (const MagneticField< Device > &magnetic_field, const int node) const |
void | draw_ascii_grid (MagneticField< Device > magnetic_field, double rmin, double rmax, double dr, double zmin, double zmax, double dz) |
Public Attributes | |
int | ntriangle |
Number of grid triangles. More... | |
int | nnode |
Number of grid nodes. More... | |
int | nplanes |
Number of planes. More... | |
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 modeled) More... | |
double | delta_phi |
Distance between planes. More... | |
double | inv_delta_phi |
1/delta_phi More... | |
Kokkos::View< double *, Kokkos::LayoutRight, Device > | psi |
An array of psi coordinates. More... | |
Kokkos::View< double *, Kokkos::LayoutRight, Device > | bfield |
Magnetic field magnitude at nodes. More... | |
Kokkos::View< int *, Kokkos::LayoutRight, Device > | basis |
A basis for the guesses? More... | |
Kokkos::View< Vertex *, Kokkos::LayoutRight, Device > | nodes |
Kokkos::View< VertexMap *, Kokkos::LayoutRight, Device > | mapping |
int | nguess |
int | guess_list_len |
double | guess_min1 |
double | guess_min2 |
double | guess_max1 |
double | guess_max2 |
double | inv_guess_d1 |
double | inv_guess_d2 |
Kokkos::View< int *, Kokkos::LayoutRight, Device > | guess_list |
Kokkos::View< int **, Kokkos::LayoutRight, Device > | guess_xtable |
Kokkos::View< int **, Kokkos::LayoutRight, Device > | guess_count |
Kokkos::View< int **, Kokkos::LayoutRight, Device > | guess_list_1d |
int | npsi_surf2 |
double | minval_psi_surf2 |
double | maxval_psi_surf2 |
Kokkos::View< double *, Kokkos::LayoutRight, Device > | psi_surf2 |
int | iphi_offset |
int | npsi00 |
double | psi00min |
double | psi00max |
double | dpsi00 |
Kokkos::View< int *, Kokkos::LayoutRight, Device > | rgn |
Kokkos::View< RZPair *, Kokkos::LayoutRight, Device > | gx |
int | nwall |
Kokkos::View< int *, Kokkos::LayoutRight, Device > | wall_nodes |
Kokkos::View< int *, Kokkos::LayoutRight, Device > | node_to_wall |
double | eps_flux_surface |
The minimum difference in psi that counts as being on two different flux surfaces. More... | |
int | nrho |
double | rhomax |
double | drho |
double | phimin |
double | phimax |
Grid< Device >::Grid | ( | const MagneticField< Device > & | magnetic_field, |
int | ntriangle_in, | ||
int | nnode_in, | ||
int | nplanes_in, | ||
double | delta_phi_in, | ||
int | nguess_in, | ||
int | guess_list_len_in, | ||
double | guess_min1_in, | ||
double | guess_min2_in, | ||
double | guess_max1_in, | ||
double | guess_max2_in, | ||
double | inv_guess_d1_in, | ||
double | inv_guess_d2_in, | ||
int | npsi_surf2_in, | ||
int | iphi_offset_in, | ||
int | npsi00_in, | ||
double | psi00min_in, | ||
double | psi00max_in, | ||
double | dpsi00_in, | ||
int | nwall_in, | ||
double | wedge_angle_in, | ||
double | eq_x_psi, | ||
double * | psi_in, | ||
int * | basis_in, | ||
Vertex * | nodes_in, | ||
VertexMap * | mapping_in, | ||
int * | guess_list_in, | ||
int * | guess_xtable_in, | ||
int * | guess_count_in, | ||
int * | guess_list_1d_in, | ||
double * | psi_surf2_in, | ||
double | minval_psi_surf2_in, | ||
double | maxval_psi_surf2_in, | ||
int * | wall_nodes_in, | ||
int * | node_to_wall_in, | ||
int * | rgn_in, | ||
RZPair * | gx_in, | ||
int | nrho_in, | ||
double | rhomax_in, | ||
double | phimax_in | ||
) |
Constructor for grid
Grid< Device >::Grid | ( | GridFiles & | grid_files, |
const MagneticField< Device > & | magnetic_field, | ||
int | nplanes_in, | ||
int | wedge_n_in, | ||
int | iphi_offset_in, | ||
int | guess_table_size, | ||
int | npsi00_in, | ||
int | nrho_in, | ||
double | rhomax_in | ||
) |
Grid< Device >::Grid | ( | NLReader::NamelistReader & | nlr, |
const GridFiles & | grid_files | ||
) |
< Simulate a wedge of 2pi/sml_wedge_n of the full torus
KOKKOS_INLINE_FUNCTION void Grid< Device >::charge_search_index | ( | const MagneticField< Device > & | magnetic_field, |
const SimdVector2D & | x, | ||
const Simd< double > & | phi, | ||
const Simd< double > & | psi_in, | ||
SimdVector2D & | xff, | ||
Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
Get the grid triangle and vertex weights for a vector of locations. Name is historical and should change
[in] | magnetic_field | The magnetic field object (currently global so technically unnecessary) |
[in] | x | Vector of (r,z) coordinates |
[in] | phi | Vector of phi coordinates |
[in] | psi_in | Vector of psi coordinates |
[out] | xff | Vector of (r,z) coordinates mapped along B-field to midplane |
[out] | itr | Vector of grid triangles |
[out] | p | Vector of triangle vertex weights |
KOKKOS_INLINE_FUNCTION void Grid< Device >::charge_search_index | ( | const MagneticField< Device > & | magnetic_field, |
const SimdVector & | v, | ||
Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
Get the grid triangle and vertex weights for a vector of locations. Name is historical and should change
[in] | magnetic_field | The magnetic field object (currently global so technically unnecessary) |
[in] | v | Vector of (r,z,phi) coordinates |
[out] | itr | Vector of grid triangles |
[out] | p | Vector of triangle vertex weights |
KOKKOS_INLINE_FUNCTION void Grid< Device >::charge_search_index | ( | const MagneticField< Device > & | magnetic_field, |
const SimdVector & | v, | ||
const Simd< double > & | psi, | ||
Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
Get the grid triangle and vertex weights for a vector of locations. Name is historical and should change
[in] | magnetic_field | The magnetic field object (currently global so technically unnecessary) |
[in] | v | Vector of (r,z,phi) coordinates |
[in] | psi | Vector of psi coordinates |
[out] | itr | Vector of grid triangles |
[out] | p | Vector of triangle vertex weights |
KOKKOS_INLINE_FUNCTION void Grid< Device >::check_triangle | ( | const Simd< int > & | itr, |
Simd< bool > & | wasnt_in_triangle | ||
) | const |
Check if vector of particles are in triangle
[in] | itr | Vector of grid triangles the particles are in |
[out] | wasnt_in_triangle | Vector of whether they are in a triangle |
|
inline |
KOKKOS_INLINE_FUNCTION void Grid< Device >::get_wall_index | ( | const Simd< bool > & | just_left_the_grid, |
const SimdVector2D & | x, | ||
const Simd< int > & | itr, | ||
const SimdGridVec & | p, | ||
Simd< int > & | widx | ||
) | const |
Find nearest wall point for a vector of locations.
[in] | just_left_the_grid | Vector of whether particles have left the grid |
[in] | x | Vector of (r,z) coordinates |
[in] | itr | Vector of grid triangles the particles were last in |
[in] | p | Vector of the vertex weighting in those triangles |
[out] | widx | Vector of wall point indices |
KOKKOS_INLINE_FUNCTION bool Grid< Device >::node_is_inside_psi_range | ( | const MagneticField< Device > & | magnetic_field, |
const int | node | ||
) | const |
Checks if a grid node lies within the simulation domain (radially).
[in] | magnetic_field | the magnetic field object (currently global so technically unnecessary) |
[in] | node | grid node to check (starts from 1) |
KOKKOS_INLINE_FUNCTION void Grid< Device >::search_ptl_1d | ( | double | psi, |
double & | wp, | ||
int & | ip | ||
) | const |
KOKKOS_INLINE_FUNCTION void Grid< Device >::search_tr2 | ( | const SimdVector2D & | xy, |
Simd< int > & | itr, | ||
SimdGridVec & | pout | ||
) | const |
Determine triangle location and triangle vertex weights of vector of locations, first checking the particles have left their previous triangle, and if so doing a search
[in] | xy | Vector of (r,z) coordinates |
[out] | itr | Vector of triangle locations |
[out] | pout | Vector of triangle vertex weights |
KOKKOS_INLINE_FUNCTION void Grid< Device >::search_tr_check_guess | ( | const SimdVector2D & | x, |
const Simd< int > & | old_itr, | ||
Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
A check to see if vector of particles has stayed in the same triangle since last time. If so just recalculate p weights
[in] | x | Vector of (r,z) coordinates |
[in] | old_itr | Vector of old triangle location |
[out] | itr | Vector of new triangle location (same as old if unchanged, -1 if search is needed) |
[out] | p | Vector of p weights if in same triangle |
KOKKOS_INLINE_FUNCTION void Grid< Device >::t_coeff | ( | const SimdVector2D & | x, |
const Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
Calculate a vector of barycentric coordinates p from real-space coordinates x in triangle itr.
[in] | x | Vector of (r,z) coordinates |
[in] | itr | Vector of triangles |
[out] | p | Vector of the weights of the triangle vertices |
KOKKOS_INLINE_FUNCTION void Grid< Device >::t_coeff_mod | ( | const MagneticField< Device > & | magnetic_field, |
const SimdVector2D & | xy, | ||
const Simd< double > & | psiin, | ||
const Simd< int > & | itr, | ||
SimdGridVec & | p | ||
) | const |
KOKKOS_INLINE_FUNCTION void Grid< Device >::wedge_modulo_phi | ( | Simd< double > & | phi_mod | ) | const |
Modulo phi in place: Keep phi between 0 and wedge_angle. Best algorithm as long as the particles haven't left the wedge multiple times over
[in,out] | phi_mod | Vector of phi coordinates |
Kokkos::View<int*, Kokkos::LayoutRight,Device> Grid< Device >::basis |
A basis for the guesses?
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::bfield |
Magnetic field magnitude at nodes.
double Grid< Device >::delta_phi |
Distance between planes.
double Grid< Device >::dpsi00 |
double Grid< Device >::drho |
double Grid< Device >::eps_flux_surface |
The minimum difference in psi that counts as being on two different flux surfaces.
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_count |
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::guess_list |
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_list_1d |
int Grid< Device >::guess_list_len |
double Grid< Device >::guess_max1 |
double Grid< Device >::guess_max2 |
double Grid< Device >::guess_min1 |
double Grid< Device >::guess_min2 |
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_xtable |
double Grid< Device >::inv_delta_phi |
1/delta_phi
double Grid< Device >::inv_guess_d1 |
double Grid< Device >::inv_guess_d2 |
int Grid< Device >::iphi_offset |
double Grid< Device >::maxval_psi_surf2 |
double Grid< Device >::minval_psi_surf2 |
int Grid< Device >::nguess |
int Grid< Device >::nnode |
Number of grid nodes.
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::node_to_wall |
int Grid< Device >::nplanes |
Number of planes.
int Grid< Device >::npsi00 |
int Grid< Device >::npsi_surf2 |
int Grid< Device >::nrho |
int Grid< Device >::ntriangle |
Number of grid triangles.
int Grid< Device >::nwall |
double Grid< Device >::phimax |
double Grid< Device >::phimin |
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::psi |
An array of psi coordinates.
double Grid< Device >::psi00max |
double Grid< Device >::psi00min |
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::psi_surf2 |
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::rgn |
double Grid< Device >::rhomax |
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::wall_nodes |
double Grid< Device >::wedge_angle |
The size of the wedge (the model is periodic in phi, a angle of e.g. pi means half the tokamak is modeled)