XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
Public Member Functions | Public Attributes | List of all members
Grid< Device > Class Template Reference

#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
 

Constructor & Destructor Documentation

template<class Device>
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

Here is the call graph for this function:

template<class Device>
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 
)

Here is the call graph for this function:

template<class Device>
Grid< Device >::Grid ( NLReader::NamelistReader nlr,
const GridFiles grid_files 
)

< Simulate a wedge of 2pi/sml_wedge_n of the full torus

Here is the call graph for this function:

template<class Device>
Grid< Device >::Grid ( )
inline

Member Function Documentation

template<class Device>
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

Parameters
[in]magnetic_fieldThe magnetic field object (currently global so technically unnecessary)
[in]xVector of (r,z) coordinates
[in]phiVector of phi coordinates
[in]psi_inVector of psi coordinates
[out]xffVector of (r,z) coordinates mapped along B-field to midplane
[out]itrVector of grid triangles
[out]pVector of triangle vertex weights

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device>
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

Parameters
[in]magnetic_fieldThe magnetic field object (currently global so technically unnecessary)
[in]vVector of (r,z,phi) coordinates
[out]itrVector of grid triangles
[out]pVector of triangle vertex weights

Here is the call graph for this function:

template<class Device>
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

Parameters
[in]magnetic_fieldThe magnetic field object (currently global so technically unnecessary)
[in]vVector of (r,z,phi) coordinates
[in]psiVector of psi coordinates
[out]itrVector of grid triangles
[out]pVector of triangle vertex weights

Here is the call graph for this function:

template<class Device >
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

Parameters
[in]itrVector of grid triangles the particles are in
[out]wasnt_in_triangleVector of whether they are in a triangle

Here is the caller graph for this function:

template<class Device>
void Grid< Device >::draw_ascii_grid ( MagneticField< Device >  magnetic_field,
double  rmin,
double  rmax,
double  dr,
double  zmin,
double  zmax,
double  dz 
)
inline
template<class Device >
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.

Parameters
[in]just_left_the_gridVector of whether particles have left the grid
[in]xVector of (r,z) coordinates
[in]itrVector of grid triangles the particles were last in
[in]pVector of the vertex weighting in those triangles
[out]widxVector of wall point indices

Here is the caller graph for this function:

template<class Device>
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).

Parameters
[in]magnetic_fieldthe magnetic field object (currently global so technically unnecessary)
[in]nodegrid node to check (starts from 1)
Returns
bool true if inside, false if outside

Here is the caller graph for this function:

template<class Device >
KOKKOS_INLINE_FUNCTION void Grid< Device >::search_ptl_1d ( double  psi,
double &  wp,
int &  ip 
) const

Here is the caller graph for this function:

template<class Device >
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

Parameters
[in]xyVector of (r,z) coordinates
[out]itrVector of triangle locations
[out]poutVector of triangle vertex weights

Here is the call graph for this function:

Here is the caller graph for this function:

template<class Device >
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

Parameters
[in]xVector of (r,z) coordinates
[in]old_itrVector of old triangle location
[out]itrVector of new triangle location (same as old if unchanged, -1 if search is needed)
[out]pVector of p weights if in same triangle
template<class Device >
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.

Parameters
[in]xVector of (r,z) coordinates
[in]itrVector of triangles
[out]pVector of the weights of the triangle vertices
template<class Device>
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

Here is the caller graph for this function:

template<class Device >
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

Parameters
[in,out]phi_modVector of phi coordinates

Here is the caller graph for this function:

Member Data Documentation

template<class Device>
Kokkos::View<int*, Kokkos::LayoutRight,Device> Grid< Device >::basis

A basis for the guesses?

template<class Device>
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::bfield

Magnetic field magnitude at nodes.

template<class Device>
double Grid< Device >::delta_phi

Distance between planes.

template<class Device>
double Grid< Device >::dpsi00
template<class Device>
double Grid< Device >::drho
template<class Device>
double Grid< Device >::eps_flux_surface

The minimum difference in psi that counts as being on two different flux surfaces.

template<class Device>
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_count
template<class Device>
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::guess_list
template<class Device>
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_list_1d
template<class Device>
int Grid< Device >::guess_list_len
template<class Device>
double Grid< Device >::guess_max1
template<class Device>
double Grid< Device >::guess_max2
template<class Device>
double Grid< Device >::guess_min1
template<class Device>
double Grid< Device >::guess_min2
template<class Device>
Kokkos::View<int**,Kokkos::LayoutRight,Device> Grid< Device >::guess_xtable
template<class Device>
Kokkos::View<RZPair*,Kokkos::LayoutRight,Device> Grid< Device >::gx
template<class Device>
double Grid< Device >::inv_delta_phi

1/delta_phi

template<class Device>
double Grid< Device >::inv_guess_d1
template<class Device>
double Grid< Device >::inv_guess_d2
template<class Device>
int Grid< Device >::iphi_offset
template<class Device>
Kokkos::View<VertexMap*,Kokkos::LayoutRight,Device> Grid< Device >::mapping
template<class Device>
double Grid< Device >::maxval_psi_surf2
template<class Device>
double Grid< Device >::minval_psi_surf2
template<class Device>
int Grid< Device >::nguess
template<class Device>
int Grid< Device >::nnode

Number of grid nodes.

template<class Device>
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::node_to_wall
template<class Device>
Kokkos::View<Vertex*,Kokkos::LayoutRight,Device> Grid< Device >::nodes
template<class Device>
int Grid< Device >::nplanes

Number of planes.

template<class Device>
int Grid< Device >::npsi00
template<class Device>
int Grid< Device >::npsi_surf2
template<class Device>
int Grid< Device >::nrho
template<class Device>
int Grid< Device >::ntriangle

Number of grid triangles.

template<class Device>
int Grid< Device >::nwall
template<class Device>
double Grid< Device >::phimax
template<class Device>
double Grid< Device >::phimin
template<class Device>
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::psi

An array of psi coordinates.

template<class Device>
double Grid< Device >::psi00max
template<class Device>
double Grid< Device >::psi00min
template<class Device>
Kokkos::View<double*,Kokkos::LayoutRight,Device> Grid< Device >::psi_surf2
template<class Device>
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::rgn
template<class Device>
double Grid< Device >::rhomax
template<class Device>
Kokkos::View<int*,Kokkos::LayoutRight,Device> Grid< Device >::wall_nodes
template<class Device>
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)


The documentation for this class was generated from the following files: