XGCa
|
Data Types | |
type | resamp_bin_type |
Public Member Functions | |
subroutine | set_cce_buffer_id (id, inode1, inode2) |
subroutine | set_cce_buffer_and_overlap_id (id1, id2, inode1, inode2) |
subroutine | resample_init (grid, inode1, inode2, dsmu, dvp, nmu, nvp, vp_max, mu_max, tile_size, bins, cce_sync_sendrcv) |
This routine initializes the resampling module and sets some global parameters such as v-grid size, etc. More... | |
subroutine | resample_finalize (bins) |
This routine cleans up the resampling module. More... | |
subroutine | resample_spall_setup (inode1, inode2, dsmu, dvp, nmu, nvp, vp_max, mu_max, cce_sync_sendrcv, tile_size) |
Resample the marker particles in XGC to improve phase-space coverage and reduce marker noise. More... | |
subroutine | set_eq_bins (eq_bin, ineq_bin) |
subroutine | use_default_eq_bins () |
subroutine | resample_one_sp (isp, cce_sync_sendrcv) |
Resample the marker particles of one species in XGC to improve phase-space coverage and reduce marker noise. More... | |
subroutine | transpose_from_aos_wrapper (isp) |
subroutine | resample_spall_finish () |
subroutine | build_bins (sp, bins, grid, cce_sync) |
Loops over particles, puts them in the appropriate bin, computes partial mean This routine sorts particles in the appropriate phase-space bin and computes the total charge in the bin. More... | |
subroutine | add_to_bin (bin, ptli, ind) |
Adds the data of one particle to a phase-space bin. More... | |
subroutine | destroy_bins (bins) |
Cleans up bin memory without deallocating the array of bins itself (can be kept for next species) More... | |
subroutine | update_cmat (grid, psn, nconstraint, nptl, vert, patch, bary, tri, psi, B, bvec, spec, pindex, part, Cmat, op_mode) |
This subroutine sets up the matrix C (constraint matrix) with the coefficients needed to compute the conserved quantities charge, parallel momentum, magnetic moment, parallel and perp. kinetic energy and the grid charge on the vertices connected to the bin. Those quantities are then obtained by multiplying the matrix C with the delta-f (–> w1) weight vector from the right. This matrix forms the basis for the QP optimization program that finds the delta-f weights for the resampled particle set. More... | |
subroutine | check_bin (bin, sp, cce_sync_sendrcv, resamp_auto_upsample_loc, resamp_auto_downsample_loc, resamp_var_resample_loc, outside_bin) |
This routine checks whether a bin needs to be resampled and determines the type of resampling (up-/down-/re-sampling). The chosen resampling option is stored in the bin data structure: 0) No resampling needed 1) Resampling due to high weight variance 2) Upsampling due to low particle count 3) Downsampling due to high particle count. More... | |
subroutine | resample_bin (bin, grid, sp, psn, cce_sync_sendrcv, start_wr_pos, resamp_failed_bins_loc, resamp_upsample_failed_loc, resamp_downsample_failed_loc, resamp_resample_failed_loc, resamp_retried_bins_loc, resamp_retried_failed_loc, resamp_fullf_failed_loc, resamp_phvol_added_loc, resamp_var_inc_upsamp_loc, resamp_var_inc_downsamp_loc, dump_problem_bins, resamp_var_inc_loc, resamp_up_hist_loc, resamp_down_hist_loc) |
The actual process of resampling the particle set in a bin happens in this routine. A bin is passed as input, the pre-determined type of resampling is executed and the weights of the new particle set are computed with a QP optimization. If the optimization is successful, the old particle set in the input variable "sp" is replaced with the new particle set. More... | |
subroutine | load_new_ptl (grid, new_ptl, bin_id, itr, p, psi) |
This routine loads new particles into a bin with uniform marker density in configuration space. Particles are drawn from a rectacngle enclosing the bin volume and are accepted if inside the Voronoi volume of the bin. A more efficient loading can be achieved by pre-selecting the triangle in which to load. This functionality will be provided in a separate function. More... | |
subroutine | load_shift_p_ptl (grid, new_ptl, itr, p, psi, my_node, max_shift) |
subroutine | setup_fullf (ptli, spec_type, den, flow, temp_ev, bvec, B, psi) |
Placeholder routine for velocity space shift, may not be implemented. More... | |
subroutine | get_new_fullf_weight (grid, psn, part, psi, B, flow, den, temp_ev, spec, tri, bary, mark_den, err) |
This routine loads new particles into a bin with uniform marker density in configuration space. Particles are drawn from the bin volume. First, a triangle is chosen using importance sampling with the triangle area as criterion. Then a random position in the triangle is drawn. The drawn position is accepted if inside the Voronoi volume of the bin. More... | |
subroutine | setup_opt (nconstraint, part_num, Cmat, constraint, weights, n_eq_bin, n_ineq_bin, eq_bin, ineq_bin, grid_ineq_on, ineq_tol, positivity, ierr) |
This routine sets up the input for the QP optimization that computes the delta-f weights of the new particles. There are two types of constraints, bin and grid constraints. Bin constraints correspond to moments of the distribution function evaluated with nearest neighbor interpolation on the configuration space grid. Grid constraints correspond to moments evaluated with linear (barycentric) interpolation on the conf. space grid. Bin constraints can be enforced as equality and inequality constraints: use idx_constraint and idx_constraint2 to choose which constraint is an equality and which is an inequality constraint. Grid constraints are enforced as inequality constraints. Since grid constraints make the optimization problem much harder due to their large number, using them is optional (grid_ineq_on) A positivity constraint for the weights can be activated with the variable "positivity". This is used for full-f optimization. More... | |
subroutine | restart_resample_set_params () |
Move phase space density from the 5D grid over to the particles for restarting simulatin with a different grid. More... | |
subroutine | restart_resample_restore_params () |
Move phase space density from the 5D grid over to the particles for restarting simulatin with a different grid. More... | |
subroutine | get_marker_velo (grid, psn, ptli, isp, velo, vrad) |
Evaluates the marker particle equations of motion. The results are the time derivatives of the marker position and parallel velocity (or rho, to be precise), as well as the radial velocity v.grad(psi). The results can be used in the moment preserving QP weight optimization. More... | |
subroutine | get_sampling_ratio (grid) |
This routine upsamples the marker population in spall for a new configuration space grid. The new grid should be denser than the current grid, i.e. require upsampling, not downsampling. This is due to particle data locality. Upsampling to a denser mesh is done on the coarse mesh, while downsampling from a fine mesh is better done with particle data arranged for the coarse mesh. The routine reads in the vertex positions of the new mesh and runs a search for each of them to associate them with the vertices of the current mesh. The ratio of new mesh vertices per current mesh vertex is the upsampling ratio. If the upsampling ratio is smaller than unity, that vertex is skipped. More... | |
subroutine | resamp_scale_phvol (grid, sp) |
Loops over particles and scales their phase space volume to account for particles added to empty bins. More... | |
subroutine | f0_upsampling_set_params () |
set parameters for f0 upsampling More... | |
subroutine | f0_upsampling_restore_params () |
restore parameters after f0 upsampling More... | |
integer function | cce_do_resampling (idbuf) |
integer function | cce_buffer_sndrcvs_accessor (idbuf) |
integer function | get_cce_buffer_nb () |
subroutine | set_cce_buffers (istep, idbuf) |
subroutine | cce_reset_particle_outside (grid, sp) |
subroutine | cce_buffer_snd_updated_bins (isp) |
subroutine | cce_buffer_rcv_updated_bins (isp) |
logical function | is_first_call () |
subroutine | nrm_std_dev (npart, mean, weights, var) |
subroutine | output_bin (npnew, nconst, sp, bin, Cmat, constraint, bin_new_ptl, new_tri, new_psi, ierr, dump_problem_bins, weights) |
Routine for outputting a flagged resampling bin. This is typically due to either the failure of the quadratic program OR a high variance result, and is controlled by the input flag resamp_output_problem_bins. More... | |
subroutine | get_v_part (smu_n, vp_n, mu, rho, temp, B, spec) |
This subroutine returns particle velocities given normalized v-grid coordinates. More... | |
subroutine | poisson_disc (initial, corners, bin_part, bin_target, new_set, k) |
logical function | in_neighborhood (point, indices, grid, new_set, cellsize, radius_squared) |
logical function | in_bin (point, corners) |
subroutine | get_v_n (smu_n, vp_n, B, temp_ev, part, spec) |
This subroutine returns normalized v-grid coordinates for a particle. More... | |
subroutine | get_refer_v_n (smu_n, vp_n, smu_ref, vp_ref, bin_id, err) |
This subroutine returns reference velocity coordinates in the resampling bin. More... | |
real(kind=8) function, dimension(2, 2) | local_lin_int (smu_ref, vp_ref) |
This function returns the 2x2 local element matrix of bilinear nodal interpolation coefficients corresponding to the particles reference coordinates in the cell. More... | |
subroutine | find_bin (smu_n, vp_n, imu, ivp) |
This subroutine returns the v-grid bin a particle is in given its normalized coordinates, including high velocity bins. More... | |
subroutine | refer_to_v_n (smu_ref, vp_ref, smu_n, vp_n, bin_id) |
This subroutine maps reference coordinates to normalized coordinates. More... | |
subroutine | particle_not_in_vbin (smu_ref, vp_ref, smu_n, vp_n, bin, part, num_part, index, rtemp, B, b_grid) |
This subroutine is error handling for when a particle(usually generated) does not like in the correct velocity space bin. More... | |
Public Attributes | |
integer | resamp_rate = 2 |
timesteps between subsequent resamples, placedholder, in practice ~ sml_f_source_period More... | |
real(8) | resamp_var = 1D4 |
threshold for relative standard deviation in bin for auto-resample More... | |
real(8) | resamp_min_ratio = 0.5D0 |
min ratio of (# of ptl)/(target # of ptl) in bin for auto-upsample More... | |
real(8) | resamp_max_ratio = 1.5D0 |
max ratio of (# of ptl)/(target # of ptl) in bin for auto-downsample More... | |
integer | resamp_max_target = 4 |
Overrides the number of constraints in determining the target # of ptl of a bin. More... | |
integer | resamp_tile_size = 2 |
Bin size on the velocity space grid in cells (not vertices) (input parameter) More... | |
integer | resamp_tile_size_now = 2 |
Used to override resamp_tile_size. More... | |
real(8) | resamp_ineq_tol = 1D-4 |
Threshold for relative error in the inequality constraints in the QP optimization. More... | |
logical | resamp_restart = .false. |
Perform resampling before dumping the final restart file. More... | |
logical | resamp_retry = .false. |
Retry QP optimization for failed bins with relaxed inequality constraints. More... | |
real(8) | resamp_ineq_tol_max = 1D-3 |
Maximal threshold for relative error in inequality constraints for retried bins. More... | |
real(8) | resamp_highv_max = 10D0 |
energy cutoff of the high velocity bins v_para>f0_vp_max and v_perp>f0_smu_max More... | |
real(8) | resamp_highv_max_ratio = 4D0 |
Downsampling threshold for high-velocity bins. More... | |
real(8) | resamp_var_limit = 3D0 |
Increase in relative bin variance for flagging for possible rejection. More... | |
logical | resamp_discard_var_bins = .false. |
Discard resampled bins that increase the variance by factor of resamp_var_limit. More... | |
logical | resamp_keep_upsamples = .false. |
Retain upsampling results with high variance, for filling for pseudoinverse. Only relevant if resamp_discard_var_bins is .true. More... | |
logical | resamp_keep_downsamples = .false. |
Retain downsampling results with high variance, mainly for preventing buildup. Only relevant if resamp_discard_var_bins is .true. More... | |
character(len=200) | resamp_node_file = 'dum.node' |
File containing the vertex positions of the mesh for which to resample. More... | |
integer | resamp_nphi_new =1 |
Number of poloidal planes in simulation with new mesh. More... | |
logical | resamp_restart_read =.false. |
Whether to read a restart file written from a simulation with different grid. More... | |
logical | resamp_fill_empty =.false. |
Whether to fill empty bins. More... | |
logical | resamp_fullf_on =.false. |
Whether to resample the full-f weights in addition to delta-f weights. More... | |
logical | resamp_grid_ineq_on =.false. |
Switch for using inequality constraints for the grid charge for resampling. More... | |
logical | resamp_distribute_evenly_subbins =.false. |
Whether to fill/remove evenly in 1x1 velocity cells in the bin when resamp_tile_size > 1 and upsampling/downsampling. More... | |
logical | resamp_fill_empty_subbins =.false. |
Whether to fill all empty 1x1 velocity cells in the bin and there are already particles in the bin. More... | |
logical | resamp_fill_empty_subbins_skip_full_bins |
If resamp_fill_empty_subbins=.true., skip a bin if it is already filled enough to do the pseudo-inverse interpolation. More... | |
integer | resamp_fill_empty_subbins_target =1 |
If resamp_fill_empty_subbins=.true., minimum target of sub bins. More... | |
integer | resamp_fill_empty_subbins_corner_cell_target =1 |
If resamp_fill_empty_subbins=.true., minimum target of corner cells. More... | |
logical | resamp_output_problem_bins =.false. |
Switch to output failed or high-variance bins as to .bp files. More... | |
character(16) | resamp_up_choice ='random' |
option arg for upsampled new particle selection: 'random','copy','poisson' currently added random: random velocity coordinates for new particles in bin copy: add new particles as copies of old particles with largest absolute w0*w1 weight poision: add new particles as poisson disc samples generated from the old particle set. More... | |
character(16) | resamp_down_choice ='weight' |
option arg for downsampled new particle selection: 'random','weight','weight+replace','volume' weight: keep particles with largest w0*w1 absolute weight, deterministic weight+replace: keep particles with largest w0*w1 absolute weight, allow multiple copies of same particle if split weight is still larger than next unsplit particle volume: keep particles with largest phase space volume, deterministic More... | |
real(kind=8) | resamp_max_shift = 1D-1 |
maximum shift in local coordinates for 'copy', 'weight+replace' More... | |
integer, dimension(:,:), allocatable | resamp_patches |
!< node numbers of vertices contained in triangle patch of voronoi cell More... | |
integer, dimension(:), allocatable | resamp_patch_size |
total number of patch vertices in Voronoi cell More... | |
real(kind=8), dimension(:,:), allocatable | resamp_patch_rzdims |
R-Z boundaries of the Voronoi-triangle patch. More... | |
integer | cce_buffer_id |
Index of currently resampled buffer. More... | |
integer | cce_buffer_first_node |
integer | cce_buffer_last_node |
node interval of current buffer More... | |
integer | cce_grid_nnode |
Number of node on the grid. More... | |
integer, parameter | cce_nbconstraint =5 |
real(kind=8), dimension(:,:,:,:,:,:), allocatable | fbins |
integer | cce_inode1 |
integer | cce_inode2 |
logical | cce_skip_nodes |
type(ptl_type), dimension(:), allocatable, target | local_particles |
type(resamp_bin_type), dimension(:,:,:), allocatable | resamp_bins |
logical | tmp_resamp_fill_empty |
logical | tmp_resamp_retry |
logical | tmp_resamp_fullf_on |
logical | tmp_resamp_grid_ineq_on |
logical | tmp_resamp_fill_empty_subbins |
logical | tmp_resamp_distribute_evenly_subbins |
logical | tmp_resamp_fill_empty_subbins_skip_full_bins |
real(8) | tmp_resamp_var |
real(8) | tmp_resamp_min_ratio |
real(8) | tmp_resamp_max_ratio |
integer | tmp_resamp_tile_size |
integer | tmp_resamp_max_target |
real(8) | tmp_resamp_ineq_tol |
real(8) | tmp_resamp_ineq_tol_max |
real(8) | tmp_resamp_highv_max |
real(8) | tmp_resamp_highv_max_ratio |
character(len=200) | tmp_resamp_node_file |
integer | tmp_resamp_nphi_new |
integer | tmp_resamp_fill_empty_subbins_target |
integer | tmp_resamp_fill_empty_subbins_corner_cell_target |
character(len=16) | tmp_up_choice |
character(len=16) | tmp_down_choice |
Private Attributes | |
real(kind=8), private | resamp_bin_max_fac |
threshold for auto downsample More... | |
real(kind=8), private | resamp_inv_bin_num |
Inverse of the number of macro bins in v-space. More... | |
integer, private | resamp_nmu_in |
Number of cells on original vperp/mu grid (=f0_nmu) More... | |
integer, private | resamp_nvp_in |
Number of cells on original v_para grid (=f0_nvp) More... | |
real(kind=8), private | resamp_dsmu |
v_perp/mu spacing of original v-grid (=f0_dsmu) More... | |
real(kind=8), private | resamp_dvp |
v_para spacing of original v-grid (=f0_dvp) More... | |
real(kind=8), private | resamp_mu_max |
v_perp/mu cutoff of v-grid (=f0_mu_max) More... | |
real(kind=8), private | resamp_vp_max |
v_para cutoff of v-grid (=f0_vp_max) More... | |
integer, private | resamp_musize |
Number of sampling bins in v_perp/mu direction. More... | |
integer, private | resamp_vpsize |
Number of sampling bins in v_para direction. More... | |
integer, private | resamp_inode1 |
Index of the first node of local patch of the configuration space mesh (=f0_inode1) More... | |
integer, private | resamp_inode2 |
Index of the last node of local patch of the configuration space mesh (=f0_inode2) More... | |
real(kind=8), dimension(:,:), allocatable, private | resamp_den |
Reference density per species for normalization. More... | |
real(kind=8), dimension(:,:), allocatable, private | resamp_t_ev |
Reference temperature per species for normalization. More... | |
real(kind=8), dimension(:,:), allocatable, private | resamp_inv_ph_vol |
Inverse phase-space volume element. More... | |
real(kind=8), dimension(:), allocatable, private | resamp_samp_ratio |
Upsampling ratio for restarts with finer grid. More... | |
logical, private | resamp_restart_mode =.false. |
real(kind=8), dimension(:), allocatable, private | resamp_phvol_sum |
Sum of phase space volume of all particles on a node. More... | |
real(kind=8), dimension(:), allocatable, private | resamp_phvol_added |
Sum of phase space volume of particles added in empty bins. More... | |
real(kind=8), dimension(:,:), allocatable, private | p_save |
Holds barycentric coordinates for particles of current resampled species. More... | |
real(kind=8), dimension(:), allocatable, private | psi_save |
Caches psi interpolation results for multiple use during resampling. More... | |
integer, parameter, private | resamp_hist_size = 10 |
Size of variance histogram for resampling statistics, bin width = 1 normalized std deviation. More... | |
integer, dimension(resamp_hist_size), private | resamp_up_hist |
Histogram of upsampling results. More... | |
integer, dimension(resamp_hist_size), private | resamp_down_hist |
Histogram of downsampling results. More... | |
integer, parameter, private | resamp_nongridmom = 8 |
Maximum Number of bin moments to be conserved. More... | |
integer, dimension(resamp_nongridmom), private | resamp_eq_bin |
which constraints are active as equality constraints More... | |
integer, dimension(resamp_nongridmom), private | resamp_ineq_bin |
which constraints are active as inequality constraints More... | |
integer, dimension(:), allocatable, private | resamp_off_vgrid |
Counter for number of particles that are not on the v-grid. More... | |
integer, dimension(:), allocatable, private | resamp_on_vgrid |
Counter for number of particles that are on the v-grid. More... | |
integer, dimension(:), allocatable, private | tr_save |
integer, private | resamp_empty_bins |
Number of empty sampling bins. More... | |
integer, private | resamp_resampled_bins |
Number of bins that were picked for re-/up-/down-sampling. More... | |
integer, private | resamp_failed_bins |
Number of bins for which re-/up-/down-sampling failed. More... | |
integer, private | resamp_var_resample |
Number of variance resampled bins. More... | |
integer, private | resamp_resample_failed |
Number of bins for which variance resampling failed. More... | |
integer, private | resamp_auto_upsample |
Number of upsampled bins. More... | |
integer, private | resamp_upsample_failed |
Number of bins for which upsampling failed. More... | |
integer, private | resamp_auto_downsample |
Number of downsampled bins. More... | |
integer, private | resamp_downsample_failed |
Number of bins for which downsampling failed. More... | |
integer, private | resamp_retried_bins |
Number of bins retried after failure. More... | |
integer, private | resamp_retried_failed |
Number of bins failed again after retry. More... | |
integer, private | resamp_fullf_failed |
Number of full-f failed bins. More... | |
integer, private | resamp_wrong_particles |
Number of particles in wwrong bins. More... | |
integer, private | resamp_var_inc_upsamp |
Number of bins that increased the variance tenfold while upsampling. More... | |
integer, private | resamp_var_inc_downsamp |
Number of bins that increased the variance tenfold while downsampling. More... | |
integer, private | resamp_var_inc |
Number of bins that increased the variance when resampling, RARE, usually result of ill-conditioning of matrix) More... | |
real(kind=8), private | resamp_smu_fac |
Factor used in the context of tapering off the marker density at high velocity. More... | |
real(kind=8), private | resamp_vp_fac |
Factor used in the context of tapering off the marker density at high velocity. More... | |
integer, private | resamp_err_count1 |
Counter for w0=0 errors in get_new_fullf_weight. More... | |
integer, private | resamp_err_count2 |
Counter for w0=0 errors in resample_bin. More... | |
subroutine resamp_module::add_to_bin | ( | type(resamp_bin_type), intent(inout) | bin, |
type(ptl_type), intent(in) | ptli, | ||
integer, intent(in) | ind | ||
) |
Adds the data of one particle to a phase-space bin.
[in,out] | bin | Phase-space bin data structure, type(resamp_bin_type) |
[in] | ptli | Data of a single particle, type(ptl_type) |
[in] | ind | Index of that particle in the corresponding species data structure, integer |
subroutine resamp_module::build_bins | ( | type(species_type), intent(inout) | sp, |
type(resamp_bin_type), dimension(resamp_inode1:resamp_inode2,resamp_musize+1,0:resamp_vpsize+1), intent(inout) | bins, | ||
type(grid_type), intent(in) | grid, | ||
logical, intent(in) | cce_sync | ||
) |
Loops over particles, puts them in the appropriate bin, computes partial mean This routine sorts particles in the appropriate phase-space bin and computes the total charge in the bin.
[in,out] | sp | Particle data structure for one species, type(species_type) |
[in,out] | bins | Phase-space bin data structure, type(resamp_bin_type) |
[in] | grid | XGC grid data structure, type(grid_type) |
subroutine resamp_module::cce_buffer_rcv_updated_bins | ( | integer, intent(in) | isp | ) |
subroutine resamp_module::cce_buffer_snd_updated_bins | ( | integer, intent(in) | isp | ) |
integer function resamp_module::cce_buffer_sndrcvs_accessor | ( | integer, intent(in) | idbuf | ) |
integer function resamp_module::cce_do_resampling | ( | integer, intent(in) | idbuf | ) |
subroutine resamp_module::cce_reset_particle_outside | ( | type(grid_type), intent(in) | grid, |
type(species_type), intent(out) | sp | ||
) |
subroutine resamp_module::check_bin | ( | type(resamp_bin_type), intent(inout) | bin, |
type(species_type), intent(in) | sp, | ||
integer, intent(in) | cce_sync_sendrcv, | ||
integer, intent(inout) | resamp_auto_upsample_loc, | ||
integer, intent(inout) | resamp_auto_downsample_loc, | ||
integer, intent(inout) | resamp_var_resample_loc, | ||
logical, intent(in) | outside_bin | ||
) |
This routine checks whether a bin needs to be resampled and determines the type of resampling (up-/down-/re-sampling). The chosen resampling option is stored in the bin data structure: 0) No resampling needed 1) Resampling due to high weight variance 2) Upsampling due to low particle count 3) Downsampling due to high particle count.
[in,out] | bin | Phase-space bin data structure, type(resamp_bin_type) |
[in] | sp | Particle data, type(species_type) |
[in] | cce_sync_sendrcv | for the coupling |
[in,out] | resamp_auto_upsample_loc | Counter for upsampled bins |
[in,out] | resamp_auto_downsample_loc | Counter for downsampled bins |
[in,out] | resamp_var_resample_loc | Counter for variance-resampled bins |
We must pass all resamp counters to be updated correctly with OMP.
subroutine resamp_module::destroy_bins | ( | type(resamp_bin_type), dimension(resamp_inode1:resamp_inode2,1:resamp_musize+1,0:resamp_vpsize+1), intent(inout) | bins | ) |
Cleans up bin memory without deallocating the array of bins itself (can be kept for next species)
[in,out] | bins | Bin data structure, type(resamp_bin_type) |
subroutine resamp_module::f0_upsampling_restore_params | ( | ) |
restore parameters after f0 upsampling
subroutine resamp_module::f0_upsampling_set_params | ( | ) |
set parameters for f0 upsampling
subroutine resamp_module::find_bin | ( | real (kind=8), intent(in) | smu_n, |
real (kind=8), intent(in) | vp_n, | ||
integer, intent(out) | imu, | ||
integer, intent(out) | ivp | ||
) |
This subroutine returns the v-grid bin a particle is in given its normalized coordinates, including high velocity bins.
[in] | smu_n | normalized square root mu coordinate |
[in] | vp_n | normalized parallel velocity coordinate |
[out] | imu | integer bin mu/vperp index |
[out] | ivp | integer bin parallel velocity index |
integer function resamp_module::get_cce_buffer_nb | ( | ) |
subroutine resamp_module::get_marker_velo | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(in) | psn, | ||
type(ptl_type), intent(in) | ptli, | ||
integer, intent(in) | isp, | ||
real (kind=8), dimension(4), intent(out) | velo, | ||
real (kind=8), intent(out) | vrad | ||
) |
Evaluates the marker particle equations of motion. The results are the time derivatives of the marker position and parallel velocity (or rho, to be precise), as well as the radial velocity v.grad(psi). The results can be used in the moment preserving QP weight optimization.
[in] | grid | XGC configuration space grid, type(grid_type) |
[in] | psn | E-field/potential data structure, type(psn_type) |
[in] | ptli | Particle data (phase+constants), type(ptl_type) |
[in] | isp | Particle species index, integer |
[out] | velo | dX/dt and drho/dt, real(8) |
[out] | vrad | Radial particle velocity (dX/dt).grad(psi) |
subroutine resamp_module::get_new_fullf_weight | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(in) | psn, | ||
type(ptl_type), intent(inout) | part, | ||
real (8), intent(in) | psi, | ||
real (8), intent(in) | B, | ||
real (8), intent(in) | flow, | ||
real (8), intent(in) | den, | ||
real (8), intent(in) | temp_ev, | ||
integer, intent(in) | spec, | ||
integer, intent(in) | tri, | ||
real (8), dimension(3), intent(in) | bary, | ||
real (8), intent(in) | mark_den, | ||
logical, intent(out) | err | ||
) |
This routine loads new particles into a bin with uniform marker density in configuration space. Particles are drawn from the bin volume. First, a triangle is chosen using importance sampling with the triangle area as criterion. Then a random position in the triangle is drawn. The drawn position is accepted if inside the Voronoi volume of the bin.
[in] | grid | XGC configuration space grid, type(grid_type) |
[in,out] | new_ptl | Particle to store data in |
[in] | bin_id | ID of the bin into which particle is loaded |
[out] | itr | Triangle into which the particle was loaded, integer |
[out] | p | Barycentric coordinates of the loaded triangle, real(8) |
[out] | psi | Pol. flux at new particle's location, real(8) Compute the new full-f weight of the input particle \(w_0 = f/g\). Note that either \(w_0\) or \(f\) can be chosen freely, while the particle's phasespace volume \(\delta V = 1/g\) is determined by Liousville's theorem of phase space volume conservation. The distribution function \(f\) can in principle be anything here as long as the \(w_2\) weight is set accordingly. The goal is to choose \(w_0\) and \(w_2\) such that the weight evolution equation does not have to be changed. This is from Seung-Hoe Ku's notes: \(f_p = w_1 w_0 g\) \(\frac{\mathrm{D}f_p}{\mathrm{D}t} = w_0 g \frac{\mathrm{d}w_1}{\mathrm{d}t}\) \(=-\frac{\mathrm{D}f_0}{\mathrm{D}t} + S\) \(w_1(t_2)-w_1(t_1) = -\frac{f_0(t_2)-f_0(t_1)}{w_0 g}\) \(w_2(t_2)-w_2(t=0) = -\frac{f_0(t_2)-f_0(t=0)}{w_0 g}\) Easiest choice: \(w_0 g = f_0(t=0),\quad w_2(t=0) = 0\) Then: \(1-w_2 = \frac{f_0(t)}{w_0 g} = \frac{f_0(t)}{f_0(t=0)}\) and \(w_1(t_2)-w_1(t_1) = -\frac{f_0(t_2)}{f_0(t=0)} + (1-w_2(t_1))\) |
Over time \(f_0(t) = f_M + f_g(t)\) can become negative, so \(w_0=f_0(t)/g\) with \(w_2(t)=0\) would be negative. But we want positive full-f weights. Therefore, we choose \(w_0=\max(f_0(t_{resample}),\alpha f_M)/g > 0\). This yields \(1-w_2(t) = 1-w_2(t_{resample})-\frac{f_0(t_{resample})}{w_0 g} + \frac{f_0(t)}{w_0 g}\) In order to get \(1-w_2 = \frac{f_0(t)}{w_0 g}\), we need to set \(w_2(t_{resample})\) such that \(1- w_2(t_{resample})-\frac{f_0(t_{resample})}{w_0 g}=0\).
[in] | grid | XGC grid data structure, type(grid_type) |
[in] | psn | Field data, type(psn_type) |
[in,out] | part | Input particle data, type(ptl_type) |
[in] | psi | Poloidal flux at particle position, real(8) |
[in] | B | Magnetic field |B_0|, real(8) |
[in] | flow | Background mean parallel flow, real(8) |
[in] | den | Background density, real(8) |
[in] | temp_ev | Background temperature |
[in] | spec | Species type of input particle, integer |
[in] | tri | Triangle in which the particle is, integer |
[in] | bary | Corresponding barycentric coordinates, real(8) |
[in] | mark_den | Inverse of particle's phase space volume g=1/dV, real(8) |
[out] | err | Error flag, logical |
subroutine resamp_module::get_refer_v_n | ( | real (kind=8), intent(in) | smu_n, |
real (kind=8), intent(in) | vp_n, | ||
real (kind=8), intent(out) | smu_ref, | ||
real (kind=8), intent(out) | vp_ref, | ||
integer, dimension(3), intent(in) | bin_id, | ||
logical, intent(out) | err | ||
) |
This subroutine returns reference velocity coordinates in the resampling bin.
[in] | smu_n | normalized square root mu coordinate |
[in] | vp_n | normalized parallel velocity coordinate |
[out] | smu_ref | reference (bin) square root mu coordinate |
[out] | vp_ref | reference (bin) parallel velocity coordinate |
[in] | bin_id | integer array of 3D bin locators |
[out] | err | logical, returns .true. if ref. coords. not in [0,1]x[0,1] |
subroutine resamp_module::get_sampling_ratio | ( | type(grid_type), intent(in) | grid | ) |
This routine upsamples the marker population in spall for a new configuration space grid. The new grid should be denser than the current grid, i.e. require upsampling, not downsampling. This is due to particle data locality. Upsampling to a denser mesh is done on the coarse mesh, while downsampling from a fine mesh is better done with particle data arranged for the coarse mesh. The routine reads in the vertex positions of the new mesh and runs a search for each of them to associate them with the vertices of the current mesh. The ratio of new mesh vertices per current mesh vertex is the upsampling ratio. If the upsampling ratio is smaller than unity, that vertex is skipped.
[in] | grid | XGC configuration space grid, type(grid_type) |
subroutine resamp_module::get_v_n | ( | real (kind=8), intent(out) | smu_n, |
real (kind=8), intent(out) | vp_n, | ||
real (kind=8), intent(in) | B, | ||
real (kind=8), intent(in) | temp_ev, | ||
type (ptl_type), intent(in) | part, | ||
integer, intent(in) | spec | ||
) |
This subroutine returns normalized v-grid coordinates for a particle.
[in] | B | magnetic field magnitude |
[in] | temp_ev | particle/reference temperature |
[in] | part | type ptl_type, current particle |
[in] | spec | species index, integer |
[out] | smu_n | normalized square root mu coordinate |
[out] | vp_n | normalized parallel velocity coordinate |
[in] | grid_b | nodal magnetic field component |
subroutine resamp_module::get_v_part | ( | real (kind=8), intent(in) | smu_n, |
real (kind=8), intent(in) | vp_n, | ||
real (kind=8), intent(out) | mu, | ||
real (kind=8), intent(out) | rho, | ||
real (kind=8), intent(in) | temp, | ||
real (kind=8), intent(in) | B, | ||
integer, intent(in) | spec | ||
) |
This subroutine returns particle velocities given normalized v-grid coordinates.
[in] | smu_n | normalized square root mu coordinate |
[in] | vp_n | normalized parallel velocity coordinate |
[out] | mu | particle mu |
[out] | rho | particle rho |
[in] | temp | reference temperate |
[in] | B | magnetic field magnitude |
[in] | spec | species index, integer |
[in] | nod_b4 | nodal magnetic field component |
logical function resamp_module::in_bin | ( | real (kind=8), dimension(2), intent(in) | point, |
real (kind=8), dimension(2,2), intent(in) | corners | ||
) |
logical function resamp_module::in_neighborhood | ( | real (kind=8), dimension(2), intent(in) | point, |
integer, dimension(2), intent(in) | indices, | ||
integer, dimension(:,:), intent(in) | grid, | ||
real (kind=8), dimension(:,:), intent(in) | new_set, | ||
real (kind=8), intent(in) | cellsize, | ||
real (kind=8), intent(in) | radius_squared | ||
) |
logical function resamp_module::is_first_call | ( | ) |
subroutine resamp_module::load_new_ptl | ( | type(grid_type), intent(in) | grid, |
type(ptl_type), intent(inout) | new_ptl, | ||
integer, dimension(3), intent(in) | bin_id, | ||
integer, intent(out) | itr, | ||
real (kind=8), dimension(3), intent(out) | p, | ||
real (kind=8), intent(out) | psi | ||
) |
This routine loads new particles into a bin with uniform marker density in configuration space. Particles are drawn from a rectacngle enclosing the bin volume and are accepted if inside the Voronoi volume of the bin. A more efficient loading can be achieved by pre-selecting the triangle in which to load. This functionality will be provided in a separate function.
[in] | grid | XGC configuration space grid, type(grid_type) |
[in,out] | new_ptl | Particle to store data in |
[in] | bin_id | ID of the bin into which particle is loaded |
[out] | itr | Triangle into which the particle was loaded, integer |
[out] | p | Barycentric coordinates of the loaded triangle, real(8) |
[out] | psi | Pol. flux at new particle's location, real(8) |
subroutine resamp_module::load_shift_p_ptl | ( | type(grid_type), intent(in) | grid, |
type(ptl_type), intent(inout) | new_ptl, | ||
integer, intent(inout) | itr, | ||
real (kind=8), dimension(3), intent(out) | p, | ||
real (kind=8), intent(out) | psi, | ||
integer, intent(in) | my_node, | ||
real (kind=8), intent(in) | max_shift | ||
) |
real (kind=8) function, dimension(2,2) resamp_module::local_lin_int | ( | real (kind=8), intent(in) | smu_ref, |
real (kind=8), intent(in) | vp_ref | ||
) |
This function returns the 2x2 local element matrix of bilinear nodal interpolation coefficients corresponding to the particles reference coordinates in the cell.
[in] | smu_ref | real reference coordinate of square-root mu variable |
[in] | vp_ref | real reference coordinate of parallel velocity |
subroutine resamp_module::nrm_std_dev | ( | integer, intent(in) | npart, |
real (8), intent(in) | mean, | ||
real (8), dimension(:), intent(in) | weights, | ||
real (8), intent(out) | var | ||
) |
subroutine resamp_module::output_bin | ( | integer, intent(in) | npnew, |
integer, intent(in) | nconst, | ||
type (species_type), intent(in) | sp, | ||
type (resamp_bin_type), intent(in) | bin, | ||
real (kind=8), dimension(nconst,npnew), intent(in) | Cmat, | ||
real (kind=8), dimension(nconst), intent(in) | constraint, | ||
type (ptl_type), dimension(npnew), intent(in) | bin_new_ptl, | ||
integer, dimension(:), intent(in) | new_tri, | ||
real (kind=8), dimension(:), intent(in) | new_psi, | ||
integer, intent(in) | ierr, | ||
logical, intent(inout) | dump_problem_bins, | ||
real (kind=8), dimension(:), intent(in), optional | weights | ||
) |
Routine for outputting a flagged resampling bin. This is typically due to either the failure of the quadratic program OR a high variance result, and is controlled by the input flag resamp_output_problem_bins.
[in] | npnew | integer, number of target particles for resampling |
[in] | nconst | integer, number of total constraints for resampling, for outputing Cmat |
[in] | sp | species_type, current species being resampled |
[in] | bin | resamp_bin_type, current bin being resampled and output |
[in] | bin_new_ptl | ptl_type, current set of candidate particle locations after resampling |
[in] | Cmat | real array, constraint matrix (including inequalities) for quadratic program |
[in] | constraint | constraint vector for QP |
[in] | ierr | integer, error flag |
[in] | new_tri | integer vector of triangles for candidate particles in bin |
[in] | new_psi | real 8 vector of canditate particles psi values |
[in,out] | dump_problem_bins | thread-local value to prevent massive IO from occuring |
in,optional] | weights if QP produced solution with unaccetpable variance, new weights |
subroutine resamp_module::particle_not_in_vbin | ( | real (kind=8), intent(in) | smu_ref, |
real (kind=8), intent(in) | vp_ref, | ||
real (kind=8), intent(in) | smu_n, | ||
real (kind=8), intent(in) | vp_n, | ||
type (resamp_bin_type), intent(in) | bin, | ||
type (ptl_type), intent(in) | part, | ||
integer, intent(in) | num_part, | ||
integer, intent(in) | index, | ||
real (kind=8), intent(in) | rtemp, | ||
real (kind=8), intent(in) | B, | ||
real (kind=8), intent(in) | b_grid | ||
) |
This subroutine is error handling for when a particle(usually generated) does not like in the correct velocity space bin.
[in] | smu_ref | computed reference smu coordinate |
[in] | vp_ref | computed refernce para velocity coordinate |
[in] | smu_n | normalized mu velocity coordinate |
[in] | vp_n | normalized parallel velocity coordinate |
[in] | bin | resamp_bin_type, current resampling bin |
[in] | part | ptl_type, current particle |
[in] | num_part | integer, target number of particles |
[in] | index | integer, which candidate particle in not_in_bin |
[in] | rtemp | particles reference temperate |
[in] | B | magnetic field magnitude at particle location |
[in] | b_grid | real, needed component of grid magnetic field for normalization |
subroutine resamp_module::poisson_disc | ( | real (kind=8), dimension(:,:), intent(in) | initial, |
real (kind=8), dimension(:,:), intent(in) | corners, | ||
integer, intent(in) | bin_part, | ||
integer, intent(in) | bin_target, | ||
real (kind=8), dimension(:,:), intent(inout) | new_set, | ||
integer, intent(in) | k | ||
) |
subroutine resamp_module::refer_to_v_n | ( | real (kind=8), intent(in) | smu_ref, |
real (kind=8), intent(in) | vp_ref, | ||
real (kind=8), intent(out) | smu_n, | ||
real (kind=8), intent(out) | vp_n, | ||
integer, dimension(3), intent(in) | bin_id | ||
) |
This subroutine maps reference coordinates to normalized coordinates.
[out] | smu_n | normalized square root mu coordinate |
[out] | vp_n | normalized parallel velocity coordinate |
[in] | smu_ref | reference (bin) square root mu coordinate |
[in] | vp_ref | reference (bin) parallel velocity coordinate |
[in] | bin_id | integer array of 3D bin locators |
subroutine resamp_module::resamp_scale_phvol | ( | type(grid_type), intent(in) | grid, |
type(species_type), intent(inout) | sp | ||
) |
Loops over particles and scales their phase space volume to account for particles added to empty bins.
[in] | grid | XGC grid data structure, type(grid_type) |
[in,out] | sp | Particle data structure for one species, type(species_type) |
subroutine resamp_module::resample_bin | ( | type (resamp_bin_type), intent(inout) | bin, |
type (grid_type), intent(in) | grid, | ||
type (species_type), intent(inout) | sp, | ||
type (psn_type), intent(in) | psn, | ||
integer, intent(in) | cce_sync_sendrcv, | ||
integer, intent(in) | start_wr_pos, | ||
integer, intent(inout) | resamp_failed_bins_loc, | ||
integer, intent(inout) | resamp_upsample_failed_loc, | ||
integer, intent(inout) | resamp_downsample_failed_loc, | ||
integer, intent(inout) | resamp_resample_failed_loc, | ||
integer, intent(inout) | resamp_retried_bins_loc, | ||
integer, intent(inout) | resamp_retried_failed_loc, | ||
integer, intent(inout) | resamp_fullf_failed_loc, | ||
real (8), intent(inout) | resamp_phvol_added_loc, | ||
integer, intent(inout) | resamp_var_inc_upsamp_loc, | ||
integer, intent(inout) | resamp_var_inc_downsamp_loc, | ||
logical, intent(inout) | dump_problem_bins, | ||
integer, intent(inout) | resamp_var_inc_loc, | ||
integer, dimension(resamp_hist_size), intent(inout) | resamp_up_hist_loc, | ||
integer, dimension(resamp_hist_size), intent(inout) | resamp_down_hist_loc | ||
) |
The actual process of resampling the particle set in a bin happens in this routine. A bin is passed as input, the pre-determined type of resampling is executed and the weights of the new particle set are computed with a QP optimization. If the optimization is successful, the old particle set in the input variable "sp" is replaced with the new particle set.
[in,out] | bin | Phase-space bin data structure, type(resamp_bin_type) |
[in] | grid | XGC grid data structure, type(grid_type) |
[in,out] | sp | Particle data structure, type(species_type) |
[in] | psn | Field data, type(psn_type) |
[in] | cce_sync_sendrcv | 0:no f-coupling, 1 sender, 2 receiver |
[in] | start_wr_pos | Write position in sp, integer |
[in,out] | resamp_failed_bins_loc | Counter for failed bins |
[in,out] | resamp_upsample_failed_loc | Counter for failed upsampled bins |
[in,out] | resamp_downsample_failed_loc | Counter for failed downsampled bins |
[in,out] | resamp_resample_failed_loc | Counter for failed var-resampled bins |
[in,out] | resamp_retried_bins_loc | Counter for retried bins |
[in,out] | resamp_retried_failed_loc | Counter for bins that failed on retry |
[in,out] | resamp_fullf_failed_loc | Counter for bins with failed full-f resampling |
[in,out] | resamp_phvol_added_loc | Phase space volume of particles added to empty bins |
[in,out] | resamp_var_inc_upsamp_loc | Counter for bins with large relative variance after upsampling |
[in,out] | resamp_var_inc_downsamp_loc | Counter for bins with large relative variance after downsampling |
subroutine resamp_module::resample_finalize | ( | type(resamp_bin_type), dimension(:,:,:), allocatable | bins | ) |
This routine cleans up the resampling module.
[in,out] | bins | Resampling bin data structure, type(resamp_bin_type) |
subroutine resamp_module::resample_init | ( | type(grid_type), intent(in) | grid, |
integer, intent(in) | inode1, | ||
integer, intent(in) | inode2, | ||
real (kind=8), intent(in) | dsmu, | ||
real (kind=8), intent(in) | dvp, | ||
integer, intent(in) | nmu, | ||
integer, intent(in) | nvp, | ||
real (kind=8), intent(in) | vp_max, | ||
real (kind=8), intent(in) | mu_max, | ||
integer, intent(in) | tile_size, | ||
type(resamp_bin_type), dimension(:,:,:), allocatable | bins, | ||
integer, intent(in) | cce_sync_sendrcv | ||
) |
This routine initializes the resampling module and sets some global parameters such as v-grid size, etc.
[in] | grid | XGC mesh data structure, type(grid_type) |
[in] | inode1 | Index of first vertex in the mesh domain, integer |
[in] | inode2 | Index of last vertex in the mesh domain, integer |
[in] | dsmu | v_perp/mu resolution of XGC's total-f velocity grid, integer |
[in] | dvp | v_para resolution of XGC's total-f v-grid, integer |
[in] | nmu | Number of v_perp/mu bins in total-f v-grid, integer |
[in] | nvp | Half the number of v_para bins in total-f v-grid, integer |
[in] | vp_max | v_para cutoff of total-f v-grid, real(8) |
[in] | mu_max | v_perp/mu cutoff of total-f v-grid, real(8) |
[in] | tile_size | Resampling combines tile_size x tile_size bins of the total-f v-grid to one macro bin, integer |
[in,out] | bins | Resampling phase space bins, type(resamp_bin_type) |
subroutine resamp_module::resample_one_sp | ( | integer | isp, |
integer | cce_sync_sendrcv | ||
) |
Resample the marker particles of one species in XGC to improve phase-space coverage and reduce marker noise.
[in] | isp | The species to resample |
[in] | cce_sendrcv | Takes values 0(nothing), 1(sender) , 2(receiver) |
subroutine resamp_module::resample_spall_finish | ( | ) |
subroutine resamp_module::resample_spall_setup | ( | integer | inode1, |
integer | inode2, | ||
real (kind=8) | dsmu, | ||
real (kind=8) | dvp, | ||
integer | nmu, | ||
integer | nvp, | ||
real (kind=8) | vp_max, | ||
real (kind=8) | mu_max, | ||
integer | cce_sync_sendrcv, | ||
integer | tile_size | ||
) |
Resample the marker particles in XGC to improve phase-space coverage and reduce marker noise.
[in] | inode1 | Index of first mesh vertex in domain, integer |
[in] | inode2 | Index of last mesh vertex in domain, integer |
[in] | dsmu | v_perp/mu resolution of XGC's total-f velocity grid, integer |
[in] | dvp | v_para resolution of XGC's total-f v-grid, integer |
[in] | nmu | Number of v_perp/mu bins in total-f v-grid, integer |
[in] | nvp | Half the number of v_para bins in total-f v-grid, integer |
[in] | vp_max | v_para cutoff of total-f v-grid, real(8) |
[in] | mu_max | v_perp/mu cutoff of total-f v-grid, real(8) |
[in] | tile_size | Resampling combines tile_size x tile_size bins of the total-f v-grid to one macro bin, integer, optional |
[in] | cce_sendrcv | Takes values 0(nothing), 1(sender) , 2(receiver) |
subroutine resamp_module::restart_resample_restore_params | ( | ) |
Move phase space density from the 5D grid over to the particles for restarting simulatin with a different grid.
subroutine resamp_module::restart_resample_set_params | ( | ) |
Move phase space density from the 5D grid over to the particles for restarting simulatin with a different grid.
subroutine resamp_module::set_cce_buffer_and_overlap_id | ( | integer, intent(in) | id1, |
integer, intent(in) | id2, | ||
integer, intent(in) | inode1, | ||
integer, intent(in) | inode2 | ||
) |
subroutine resamp_module::set_cce_buffer_id | ( | integer, intent(in) | id, |
integer, intent(in) | inode1, | ||
integer, intent(in) | inode2 | ||
) |
subroutine resamp_module::set_cce_buffers | ( | integer, intent(in) | istep, |
integer, intent(in) | idbuf | ||
) |
subroutine resamp_module::set_eq_bins | ( | integer, dimension(resamp_nongridmom), intent(in) | eq_bin, |
integer, dimension(resamp_nongridmom), intent(in) | ineq_bin | ||
) |
subroutine resamp_module::setup_fullf | ( | type(ptl_type), intent(in) | ptli, |
integer, intent(in) | spec_type, | ||
real (8), intent(out) | den, | ||
real (8), intent(out) | flow, | ||
real (8), intent(out) | temp_ev, | ||
real (8), dimension(3), intent(out) | bvec, | ||
real (8), intent(out) | B, | ||
real (8), intent(inout) | psi | ||
) |
Placeholder routine for velocity space shift, may not be implemented.
Evaluate the background density, temperature, flow, magnetic field and psi at the input particle's position
[in] | ptli | Input particle data, type(ptl_type) |
in[ | spec_type Species type of input particle, integer | |
[out] | den | Background density, real(8) |
[out] | flow | Background toroidal flow velocity, real(8) |
[out] | temp_ev | Background temperature, real(8) |
[out] | bvec | Magnetic field vector B, real(8) |
[out] | B | Magnetic field |B_0|, real(8) |
[out] | psi | Poloidal magnetic flux, real(8) |
subroutine resamp_module::setup_opt | ( | integer, intent(in) | nconstraint, |
integer, intent(in) | part_num, | ||
real (8), dimension(nconstraint,part_num), intent(in) | Cmat, | ||
real (8), dimension(nconstraint), intent(inout) | constraint, | ||
real (8), dimension(part_num), intent(inout) | weights, | ||
integer, intent(in) | n_eq_bin, | ||
integer, intent(in) | n_ineq_bin, | ||
integer, dimension(resamp_nongridmom), intent(in) | eq_bin, | ||
integer, dimension(resamp_nongridmom), intent(in) | ineq_bin, | ||
logical, intent(in) | grid_ineq_on, | ||
real (8), intent(in) | ineq_tol, | ||
logical, intent(in) | positivity, | ||
integer, intent(out) | ierr | ||
) |
This routine sets up the input for the QP optimization that computes the delta-f weights of the new particles. There are two types of constraints, bin and grid constraints. Bin constraints correspond to moments of the distribution function evaluated with nearest neighbor interpolation on the configuration space grid. Grid constraints correspond to moments evaluated with linear (barycentric) interpolation on the conf. space grid. Bin constraints can be enforced as equality and inequality constraints: use idx_constraint and idx_constraint2 to choose which constraint is an equality and which is an inequality constraint. Grid constraints are enforced as inequality constraints. Since grid constraints make the optimization problem much harder due to their large number, using them is optional (grid_ineq_on) A positivity constraint for the weights can be activated with the variable "positivity". This is used for full-f optimization.
[in] | nconstraint | Number of constraints, integer |
[in] | part_num | Number of particles, integer |
[in] | Cmat | Constraint matrix, N_constraint x N_ptl, real(8) |
[in] | constraint | Values of the constraints, real(8) |
[in,out] | weights | Delta-f weights w1, N_ptl, real(8) |
[in] | grid_ineq_on | Whether to conserve the grid charge or not, logical |
[in] | ineq_tol | Error tolerance for the inequality constraints, real(8) |
[in] | positivity | Whether to require weight>0, logical |
[in,out] | ierr | Error flag, integer |
subroutine resamp_module::transpose_from_aos_wrapper | ( | integer | isp | ) |
subroutine resamp_module::update_cmat | ( | type(grid_type), intent(in) | grid, |
type(psn_type), intent(in) | psn, | ||
integer, intent(in) | nconstraint, | ||
integer, intent(in) | nptl, | ||
integer, intent(in) | vert, | ||
integer, dimension(1:vert), intent(in) | patch, | ||
real (8), dimension(3), intent(in) | bary, | ||
integer, intent(in) | tri, | ||
real (8), intent(in) | psi, | ||
real (8), intent(in) | B, | ||
real (8), dimension(3), intent(in) | bvec, | ||
integer, intent(in) | spec, | ||
integer, intent(in) | pindex, | ||
type(ptl_type), intent(in) | part, | ||
real (8), dimension(nconstraint,nptl), intent(out) | Cmat, | ||
integer, intent(in) | op_mode | ||
) |
This subroutine sets up the matrix C (constraint matrix) with the coefficients needed to compute the conserved quantities charge, parallel momentum, magnetic moment, parallel and perp. kinetic energy and the grid charge on the vertices connected to the bin. Those quantities are then obtained by multiplying the matrix C with the delta-f (–> w1) weight vector from the right. This matrix forms the basis for the QP optimization program that finds the delta-f weights for the resampled particle set.
[in] | grid | XGC grid data structure, type(grid_type) |
[in] | psn | XGC E-field/potential data, type(psn_type) |
[in] | nconstraint | Number of equality+inequality constraints, integer |
[in] | nptl | Number of particles in bin, integer |
[in] | vert | Number of vertices connected to this bin (patch size), integer |
[in] | patch | Indices of the vertices connected to this bin, integer |
[in] | bary | Barycentric coordinates of the particle, real(8) |
[in] | tri | Triangle in which the particle is, integer |
[in] | psi | Pol. flux at particle position, real(8) |
[in] | B | Local magnetic field |B|, real(8) |
[in] | bvec | Magnetic field vector B, real(8) |
[in] | spec | Species index, integer |
[in] | pindex | Position of the particle in the particle set, integer |
[in] | part | Data of one particle, type(ptl_type) |
[in,out] | Cmat | The matrix C from the (see documentation), real(8) |
[in] | op_mode | Switch for optimizing full-f weights w0 (0)or delta-f weights (1), integer |
subroutine resamp_module::use_default_eq_bins | ( | ) |
integer resamp_module::cce_buffer_first_node |
integer resamp_module::cce_buffer_id |
Index of currently resampled buffer.
integer resamp_module::cce_buffer_last_node |
node interval of current buffer
integer resamp_module::cce_grid_nnode |
Number of node on the grid.
integer resamp_module::cce_inode1 |
integer resamp_module::cce_inode2 |
integer, parameter resamp_module::cce_nbconstraint =5 |
logical resamp_module::cce_skip_nodes |
real (kind=8), dimension(:,:,:,:,:,:), allocatable resamp_module::fbins |
type(ptl_type), dimension(:), allocatable, target resamp_module::local_particles |
|
private |
Holds barycentric coordinates for particles of current resampled species.
|
private |
Caches psi interpolation results for multiple use during resampling.
|
private |
Number of downsampled bins.
|
private |
Number of upsampled bins.
|
private |
threshold for auto downsample
type (resamp_bin_type), dimension(:,:,:), allocatable resamp_module::resamp_bins |
|
private |
Reference density per species for normalization.
logical resamp_module::resamp_discard_var_bins = .false. |
Discard resampled bins that increase the variance by factor of resamp_var_limit.
logical resamp_module::resamp_distribute_evenly_subbins =.false. |
Whether to fill/remove evenly in 1x1 velocity cells in the bin when resamp_tile_size > 1 and upsampling/downsampling.
character(16) resamp_module::resamp_down_choice ='weight' |
option arg for downsampled new particle selection: 'random','weight','weight+replace','volume' weight: keep particles with largest w0*w1 absolute weight, deterministic weight+replace: keep particles with largest w0*w1 absolute weight, allow multiple copies of same particle if split weight is still larger than next unsplit particle volume: keep particles with largest phase space volume, deterministic
|
private |
Histogram of downsampling results.
|
private |
Number of bins for which downsampling failed.
|
private |
v_perp/mu spacing of original v-grid (=f0_dsmu)
|
private |
v_para spacing of original v-grid (=f0_dvp)
|
private |
Number of empty sampling bins.
|
private |
which constraints are active as equality constraints
|
private |
Counter for w0=0 errors in get_new_fullf_weight.
|
private |
Counter for w0=0 errors in resample_bin.
|
private |
Number of bins for which re-/up-/down-sampling failed.
logical resamp_module::resamp_fill_empty =.false. |
Whether to fill empty bins.
logical resamp_module::resamp_fill_empty_subbins =.false. |
Whether to fill all empty 1x1 velocity cells in the bin and there are already particles in the bin.
integer resamp_module::resamp_fill_empty_subbins_corner_cell_target =1 |
If resamp_fill_empty_subbins=.true., minimum target of corner cells.
logical resamp_module::resamp_fill_empty_subbins_skip_full_bins |
If resamp_fill_empty_subbins=.true., skip a bin if it is already filled enough to do the pseudo-inverse interpolation.
integer resamp_module::resamp_fill_empty_subbins_target =1 |
If resamp_fill_empty_subbins=.true., minimum target of sub bins.
|
private |
Number of full-f failed bins.
logical resamp_module::resamp_fullf_on =.false. |
Whether to resample the full-f weights in addition to delta-f weights.
logical resamp_module::resamp_grid_ineq_on =.false. |
Switch for using inequality constraints for the grid charge for resampling.
real (8) resamp_module::resamp_highv_max = 10D0 |
energy cutoff of the high velocity bins v_para>f0_vp_max and v_perp>f0_smu_max
real (8) resamp_module::resamp_highv_max_ratio = 4D0 |
Downsampling threshold for high-velocity bins.
|
private |
Size of variance histogram for resampling statistics, bin width = 1 normalized std deviation.
|
private |
which constraints are active as inequality constraints
real (8) resamp_module::resamp_ineq_tol = 1D-4 |
Threshold for relative error in the inequality constraints in the QP optimization.
real (8) resamp_module::resamp_ineq_tol_max = 1D-3 |
Maximal threshold for relative error in inequality constraints for retried bins.
|
private |
Index of the first node of local patch of the configuration space mesh (=f0_inode1)
|
private |
Index of the last node of local patch of the configuration space mesh (=f0_inode2)
|
private |
Inverse of the number of macro bins in v-space.
|
private |
Inverse phase-space volume element.
logical resamp_module::resamp_keep_downsamples = .false. |
Retain downsampling results with high variance, mainly for preventing buildup. Only relevant if resamp_discard_var_bins is .true.
logical resamp_module::resamp_keep_upsamples = .false. |
Retain upsampling results with high variance, for filling for pseudoinverse. Only relevant if resamp_discard_var_bins is .true.
real (8) resamp_module::resamp_max_ratio = 1.5D0 |
max ratio of (# of ptl)/(target # of ptl) in bin for auto-downsample
real (kind=8) resamp_module::resamp_max_shift = 1D-1 |
maximum shift in local coordinates for 'copy', 'weight+replace'
integer resamp_module::resamp_max_target = 4 |
Overrides the number of constraints in determining the target # of ptl of a bin.
real (8) resamp_module::resamp_min_ratio = 0.5D0 |
min ratio of (# of ptl)/(target # of ptl) in bin for auto-upsample
|
private |
v_perp/mu cutoff of v-grid (=f0_mu_max)
|
private |
Number of sampling bins in v_perp/mu direction.
|
private |
Number of cells on original vperp/mu grid (=f0_nmu)
character (len=200) resamp_module::resamp_node_file = 'dum.node' |
File containing the vertex positions of the mesh for which to resample.
|
private |
Maximum Number of bin moments to be conserved.
integer resamp_module::resamp_nphi_new =1 |
Number of poloidal planes in simulation with new mesh.
|
private |
Number of cells on original v_para grid (=f0_nvp)
|
private |
Counter for number of particles that are not on the v-grid.
|
private |
Counter for number of particles that are on the v-grid.
logical resamp_module::resamp_output_problem_bins =.false. |
Switch to output failed or high-variance bins as to .bp files.
real (kind=8), dimension(:,:), allocatable resamp_module::resamp_patch_rzdims |
R-Z boundaries of the Voronoi-triangle patch.
integer, dimension(:), allocatable resamp_module::resamp_patch_size |
total number of patch vertices in Voronoi cell
integer, dimension(:,:), allocatable resamp_module::resamp_patches |
!< node numbers of vertices contained in triangle patch of voronoi cell
|
private |
Sum of phase space volume of particles added in empty bins.
|
private |
Sum of phase space volume of all particles on a node.
integer resamp_module::resamp_rate = 2 |
timesteps between subsequent resamples, placedholder, in practice ~ sml_f_source_period
|
private |
Number of bins for which variance resampling failed.
|
private |
Number of bins that were picked for re-/up-/down-sampling.
logical resamp_module::resamp_restart = .false. |
Perform resampling before dumping the final restart file.
|
private |
logical resamp_module::resamp_restart_read =.false. |
Whether to read a restart file written from a simulation with different grid.
|
private |
Number of bins retried after failure.
|
private |
Number of bins failed again after retry.
logical resamp_module::resamp_retry = .false. |
Retry QP optimization for failed bins with relaxed inequality constraints.
|
private |
Upsampling ratio for restarts with finer grid.
|
private |
Factor used in the context of tapering off the marker density at high velocity.
|
private |
Reference temperature per species for normalization.
integer resamp_module::resamp_tile_size = 2 |
Bin size on the velocity space grid in cells (not vertices) (input parameter)
integer resamp_module::resamp_tile_size_now = 2 |
Used to override resamp_tile_size.
character(16) resamp_module::resamp_up_choice ='random' |
option arg for upsampled new particle selection: 'random','copy','poisson' currently added random: random velocity coordinates for new particles in bin copy: add new particles as copies of old particles with largest absolute w0*w1 weight poision: add new particles as poisson disc samples generated from the old particle set.
|
private |
Histogram of upsampling results.
|
private |
Number of bins for which upsampling failed.
real (8) resamp_module::resamp_var = 1D4 |
threshold for relative standard deviation in bin for auto-resample
|
private |
Number of bins that increased the variance when resampling, RARE, usually result of ill-conditioning of matrix)
|
private |
Number of bins that increased the variance tenfold while downsampling.
|
private |
Number of bins that increased the variance tenfold while upsampling.
real (8) resamp_module::resamp_var_limit = 3D0 |
Increase in relative bin variance for flagging for possible rejection.
|
private |
Number of variance resampled bins.
|
private |
Factor used in the context of tapering off the marker density at high velocity.
|
private |
v_para cutoff of v-grid (=f0_vp_max)
|
private |
Number of sampling bins in v_para direction.
|
private |
Number of particles in wwrong bins.
character (len=16) resamp_module::tmp_down_choice |
logical resamp_module::tmp_resamp_distribute_evenly_subbins |
logical resamp_module::tmp_resamp_fill_empty |
logical resamp_module::tmp_resamp_fill_empty_subbins |
integer resamp_module::tmp_resamp_fill_empty_subbins_corner_cell_target |
logical resamp_module::tmp_resamp_fill_empty_subbins_skip_full_bins |
integer resamp_module::tmp_resamp_fill_empty_subbins_target |
logical resamp_module::tmp_resamp_fullf_on |
logical resamp_module::tmp_resamp_grid_ineq_on |
real (8) resamp_module::tmp_resamp_highv_max |
real (8) resamp_module::tmp_resamp_highv_max_ratio |
real (8) resamp_module::tmp_resamp_ineq_tol |
real (8) resamp_module::tmp_resamp_ineq_tol_max |
real (8) resamp_module::tmp_resamp_max_ratio |
integer resamp_module::tmp_resamp_max_target |
real (8) resamp_module::tmp_resamp_min_ratio |
character (len=200) resamp_module::tmp_resamp_node_file |
integer resamp_module::tmp_resamp_nphi_new |
logical resamp_module::tmp_resamp_retry |
integer resamp_module::tmp_resamp_tile_size |
real (8) resamp_module::tmp_resamp_var |
character (len=16) resamp_module::tmp_up_choice |
|
private |