XGCa
|
Data Types | |
interface | gradparx2_cpp |
type | grid_type |
Public Member Functions | |
subroutine | set_fortran_flx_aif_format (old_flx_aif_format_int) |
subroutine | init_fortran_grid (grid_nnode_in, grid_ntriangle_in, grid_nwall_in, grid_x_in, grid_rgn_in, grid_nd_in, grid_mapping_in, grid_wall_nodes_in, grid_node_to_wall_in, basis_in) |
subroutine | init_fortran_grid_surf (n_xpt, i_x1, i_x2, nsurf1, nsurf2, nsurf3a, nsurf3b, isurf_sep1, isurf_sep2, npsi_surf, surf_maxlen, num_non_aligned, surf_len, surf_idx, non_aligned_vert, non_aligned_nsurf, non_aligned_surf_idx) |
subroutine | init_secondary_grid_arrays (grid_psi_in) |
subroutine | init_guess_table_from_cpp (guess_n1, guess_list_len, guess_n_in, guess_min_in, guess_max_in, inv_guess_d_in, guess_list_in, guess_xtable_in, guess_count_in) |
subroutine | init_fortran_vols_and_areas (tr_area_h, tr_vol_h, |
subroutine | init_fortran_psi_grid_props (psi00_n, psi00_val_min, psi00_val_max, psi00_dval, npsi_surf2, psi_surf2, psi_guess_list) |
subroutine | alloc_fortran_gradient_mat_unit_vec (m, n, width) |
subroutine | init_grid (grid) |
Initializes the grid data structure from input files containing i) The (R,Z) coordinates of all vertices including whether they are wall vertices (.node file), ii) the connectivity of the vertices defining the triangles of the unstructured mesh (.ele file), iii) the flux-surface connectivity, i.e. which vertices are on the same flux-surface, and for vertices not aligned on flux-surfaces between which flux-surfaces those vertices lie. More... | |
subroutine | init_triangle (grid) |
subroutine | t_coeff (grid, itr, x, p) |
subroutine | t_coeff_mod (grid, xy, psiin, itr, p) |
subroutine | search_tr_with_guess (grid, x, init, itr, p, count) |
subroutine | calc_gen_theta_psi (grid) |
Calculates the average psi on the complete flux-surfaces and the generalized poloidal angle (=normalized poloidal arc length) More... | |
subroutine | search_ptl_1d (grid, psi, wp, ip) |
1D binary search on a non-uniform psi grid with initial guess based on uniform psi-grid (Assumes that flux-surfaces on the non-uniform grid are in ascending order) More... | |
subroutine | calc_mag_drift_flux_avg (grid) |
Computes the flux-surface average of 1/B^3 Bxgrad(B) and 1/B^2 curl(B). Those are zero analytically, but not numerically. This numerical cancellation error affects the calculation of radial fluxes a lot –> can be removed from the calculation. More... | |
subroutine | init_wall (grid) |
subroutine | sendright_recvleft (grid, input) |
Sends the iphi=1 plane to the toroidal domain to the right and receives the iphi=0 plane from the toroidal domain to the left. It is used for exchanging field information that is known only on the iphi=1 plane on all processes. More... | |
subroutine | sendright_recvleft_sendleft_recvright (grid, input, output) |
Sends the (R,Z) plane data to the toroidal domain to the right and receives the (R,Z) plane data from the toroidal domain to the left and vice versa to collect the data from 3 adjacent plances It is typically used for exchanging field information that is known only on the iphi=1 plane on all processes. More... | |
subroutine | send_recv_potential (dpot, recvr, sendl, nnode) |
Collect (R,Z) plane data on 4 adjacent planes iphi=(-1:2) iphi=1 is assumed to be local data on the MPI rank. More... | |
subroutine | gradparx (grid, tr, p, dx, bd, input, output) |
Calculates the parallel gradient like for the right plane (iphi=1). The input is supposed to contain the data of four adjacent planes (0:2) More... | |
subroutine | gradparx2 (grid, tr, p, dx, bd, input, output) |
Calculates the parallel gradient like GradParX, but for the left and the right plane at the same time. The input is supposed to contain the data of four adjacent planes (-1:2 !!!) The routine calculates the parallel gradient for both planes in the sector (0 and 1) More... | |
integer function | get_nsurfs_for_avg () |
Public Attributes | |
integer, parameter | grid_rgn_wall = 100 |
real(8), parameter | grid_search_tr2_psi_null = -1D10 |
integer, parameter | grid_min_surf_size = 5 |
type(grid_type), target | grid_global |
subroutine grid_class::alloc_fortran_gradient_mat_unit_vec | ( | integer(c_int) | m, |
integer(c_int) | n, | ||
integer(c_int) | width | ||
) |
subroutine grid_class::calc_gen_theta_psi | ( | type(grid_type), intent(inout) | grid | ) |
Calculates the average psi on the complete flux-surfaces and the generalized poloidal angle (=normalized poloidal arc length)
grid | (inout) type(grid_type), grid information |
subroutine grid_class::calc_mag_drift_flux_avg | ( | type(grid_type), intent(inout) | grid | ) |
Computes the flux-surface average of 1/B^3 Bxgrad(B) and 1/B^2 curl(B). Those are zero analytically, but not numerically. This numerical cancellation error affects the calculation of radial fluxes a lot –> can be removed from the calculation.
[in,out] | grid | grid data, type(grid_type) |
integer function grid_class::get_nsurfs_for_avg | ( | ) |
subroutine grid_class::gradparx | ( | type(grid_type), intent(in) | grid, |
integer, dimension(grid%nnode,0:1), intent(in) | tr, | ||
real (kind=8), dimension(3,grid%nnode,0:1), intent(in) | p, | ||
real (kind=8), dimension(grid%nnode,0:1), intent(in) | dx, | ||
type(boundary2_type), intent(in) | bd, | ||
real (kind=8), dimension(grid%nnode,3), intent(in) | input, | ||
real (kind=8), dimension(grid%nnode), intent(out) | output | ||
) |
Calculates the parallel gradient like for the right plane (iphi=1). The input is supposed to contain the data of four adjacent planes (0:2)
[in] | grid | XGC grid data structure, type(grid_type) |
[in] | tr | parallel triangle connectivity, integer |
[in] | p | interpolation weights corresponding to tr, real(8) |
[in] | dx | distances along the field line, real(8) |
[in] | bd | boundary conditions, type(boundary2_type) |
[in] | input | input data for all planes from 0 to 2, real(8) |
[out] | output | parallel derivative on planes 0 and 1, real(8) |
subroutine grid_class::gradparx2 | ( | type(grid_type), intent(in) | grid, |
integer, dimension(grid%nnode,0:1), intent(in) | tr, | ||
real (kind=8), dimension(3,grid%nnode,0:1), intent(in) | p, | ||
real (kind=8), dimension(grid%nnode,0:1), intent(in) | dx, | ||
type(boundary2_type), intent(in) | bd, | ||
real (kind=8), dimension(grid%nnode,-1:2), intent(in) | input, | ||
real (kind=8), dimension(grid%nnode,0:1), intent(out) | output | ||
) |
Calculates the parallel gradient like GradParX, but for the left and the right plane at the same time. The input is supposed to contain the data of four adjacent planes (-1:2 !!!) The routine calculates the parallel gradient for both planes in the sector (0 and 1)
[in] | grid | XGC grid data structure, type(grid_type) |
[in] | tr | parallel triangle connectivity, integer |
[in] | p | interpolation weights corresponding to tr, real(8) |
[in] | dx | distances along the field line, real(8) |
[in] | bd | boundary conditions, type(boundary2_type) |
[in] | input | input data for all planes from -1 to 2, real(8) |
[out] | output | parallel derivative on planes 0 and 1, real(8) |
subroutine grid_class::init_fortran_grid | ( | integer(c_int), intent(in) | grid_nnode_in, |
integer(c_int), intent(in) | grid_ntriangle_in, | ||
integer(c_int), intent(in) | grid_nwall_in, | ||
real (kind=8), dimension(2,grid_nnode_in), target | grid_x_in, | ||
integer(c_int), dimension(grid_nnode_in), target | grid_rgn_in, | ||
integer(c_int), dimension(3,grid_ntriangle_in), target | grid_nd_in, | ||
real (kind=8), dimension(2,3,grid_ntriangle_in), target | grid_mapping_in, | ||
integer(c_int), dimension(grid_nwall_in), target | grid_wall_nodes_in, | ||
integer(c_int), dimension(grid_nnode_in), target | grid_node_to_wall_in, | ||
integer(c_int), dimension(grid_nnode_in), target | basis_in | ||
) |
subroutine grid_class::init_fortran_grid_surf | ( | integer(c_int), intent(in) | n_xpt, |
integer(c_int), intent(in) | i_x1, | ||
integer(c_int), intent(in) | i_x2, | ||
integer(c_int), intent(in) | nsurf1, | ||
integer(c_int), intent(in) | nsurf2, | ||
integer(c_int), intent(in) | nsurf3a, | ||
integer(c_int), intent(in) | nsurf3b, | ||
integer(c_int), intent(in) | isurf_sep1, | ||
integer(c_int), intent(in) | isurf_sep2, | ||
integer(c_int), intent(in) | npsi_surf, | ||
integer(c_int), intent(in) | surf_maxlen, | ||
integer(c_int), intent(in) | num_non_aligned, | ||
integer(c_int), dimension(npsi_surf), target | surf_len, | ||
integer(c_int), dimension(surf_maxlen, npsi_surf), target | surf_idx, | ||
integer(c_int), dimension(num_non_aligned), target | non_aligned_vert, | ||
integer(c_int), dimension(num_non_aligned), target | non_aligned_nsurf, | ||
integer(c_int), dimension(4,num_non_aligned), target | non_aligned_surf_idx | ||
) |
subroutine grid_class::init_fortran_psi_grid_props | ( | integer(c_int) | psi00_n, |
real(c_double) | psi00_val_min, | ||
real(c_double) | psi00_val_max, | ||
real(c_double) | psi00_dval, | ||
integer(c_int) | npsi_surf2, | ||
real(c_double), dimension(npsi_surf2), target | psi_surf2, | ||
integer(c_int), dimension(psi00_n, 2), target | psi_guess_list | ||
) |
subroutine grid_class::init_fortran_vols_and_areas | ( | real(c_double), dimension(grid_global%ntriangle), target | tr_area_h, |
real(c_double), dimension(grid_global%ntriangle), target | tr_vol_h | ||
) |
subroutine grid_class::init_grid | ( | type(grid_type), intent(inout) | grid | ) |
Initializes the grid data structure from input files containing i) The (R,Z) coordinates of all vertices including whether they are wall vertices (.node file), ii) the connectivity of the vertices defining the triangles of the unstructured mesh (.ele file), iii) the flux-surface connectivity, i.e. which vertices are on the same flux-surface, and for vertices not aligned on flux-surfaces between which flux-surfaces those vertices lie.
[in,out] | grid | memory for grid data, type(grid_type) |
[in] | filename3 | File containing the flux-surface connectivity |
subroutine grid_class::init_guess_table_from_cpp | ( | integer(c_int) | guess_n1, |
integer(c_int) | guess_list_len, | ||
integer(c_int), dimension(2) | guess_n_in, | ||
real (kind=8), dimension(2) | guess_min_in, | ||
real (kind=8), dimension(2) | guess_max_in, | ||
real (kind=8), dimension(2) | inv_guess_d_in, | ||
integer(c_int), dimension(guess_list_len), target | guess_list_in, | ||
integer(c_int), dimension(guess_n1,guess_n1), target | guess_xtable_in, | ||
integer(c_int), dimension(guess_n1,guess_n1), target | guess_count_in | ||
) |
subroutine grid_class::init_secondary_grid_arrays | ( | real (kind=8), dimension(grid_global%nnode), target | grid_psi_in | ) |
subroutine grid_class::init_triangle | ( | type(grid_type) | grid | ) |
subroutine grid_class::init_wall | ( | type(grid_type) | grid | ) |
subroutine grid_class::search_ptl_1d | ( | type(grid_type), intent(in) | grid, |
real (kind=8), intent(in) | psi, | ||
real (kind=8), intent(out) | wp, | ||
integer, intent(out) | ip | ||
) |
1D binary search on a non-uniform psi grid with initial guess based on uniform psi-grid (Assumes that flux-surfaces on the non-uniform grid are in ascending order)
[in] | grid | (grid_type), grid-data |
[in] | psi | (real (8)), psi-value to be located |
[out] | wp | (real (8)), linear interpolation weight for the lower bound |
[out] | ip | (integer), index of the lower boundary surface |
subroutine grid_class::search_tr_with_guess | ( | type(grid_type), target | grid, |
real (kind=8), dimension(2), intent(in) | x, | ||
integer, intent(in) | init, | ||
integer, intent(out) | itr, | ||
real (kind=8), dimension(3), intent(out) | p, | ||
integer | count | ||
) |
subroutine grid_class::send_recv_potential | ( | real (kind=8), dimension(nnode,-1:2) | dpot, |
real (kind=8), dimension(nnode) | recvr, | ||
real (kind=8), dimension(nnode) | sendl, | ||
integer, intent(in) | nnode | ||
) |
Collect (R,Z) plane data on 4 adjacent planes iphi=(-1:2) iphi=1 is assumed to be local data on the MPI rank.
subroutine grid_class::sendright_recvleft | ( | type(grid_type), intent(in) | grid, |
real (kind=8), dimension(grid%nnode,0:1), intent(inout) | input | ||
) |
Sends the iphi=1 plane to the toroidal domain to the right and receives the iphi=0 plane from the toroidal domain to the left. It is used for exchanging field information that is known only on the iphi=1 plane on all processes.
[in] | grid | XGC grid data structure, type(grid_type) |
[in,out] | input | Data that is to be exchanged, real(8) |
subroutine grid_class::sendright_recvleft_sendleft_recvright | ( | type(grid_type), intent(in) | grid, |
real (kind=8), dimension(grid%nnode), intent(in) | input, | ||
real (kind=8), dimension(grid%nnode,3), intent(out) | output | ||
) |
Sends the (R,Z) plane data to the toroidal domain to the right and receives the (R,Z) plane data from the toroidal domain to the left and vice versa to collect the data from 3 adjacent plances It is typically used for exchanging field information that is known only on the iphi=1 plane on all processes.
[in] | grid | grid information,type(grid_type) |
[in] | input | The scalar that is to be exchanged, real(8) |
[out] | output | The collected data of 3 planes, real(8) |
subroutine grid_class::set_fortran_flx_aif_format | ( | integer(c_int), intent(in) | old_flx_aif_format_int | ) |
subroutine grid_class::t_coeff | ( | type(grid_type), intent(in) | grid, |
integer, intent(in) | itr, | ||
real (kind=8), dimension(2), intent(in) | x, | ||
real (kind=8), dimension(3), intent(out) | p | ||
) |
subroutine grid_class::t_coeff_mod | ( | type (grid_type), intent(in) | grid, |
real (kind=8), dimension(2), intent(in) | xy, | ||
real (kind=8), intent(in) | psiin, | ||
integer, intent(in) | itr, | ||
real (kind=8), dimension(3), intent(out) | p | ||
) |
type(grid_type), target grid_class::grid_global |
integer, parameter grid_class::grid_min_surf_size = 5 |
integer, parameter grid_class::grid_rgn_wall = 100 |
real (8), parameter grid_class::grid_search_tr2_psi_null = -1D10 |