XGCa
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Macros Pages
Data Types | Public Member Functions | Public Attributes | Private Attributes | List of all members
resamp_module Module Reference

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 (grid, spall, psn, 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 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)
 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)
 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 setup_fullf (ptli, spec_type, den, flow, temp_ev, bvec, B, psi)
 Evaluate the background density, temperature, flow, magnetic field and psi at the input particle's position. More...
 
subroutine load_new_ptl2 (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 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 get_new_fullf_weight (grid, psn, part, psi, B, flow, den, temp_ev, spec, tri, bary, mark_den, err)
 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))\). 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 move_f0g_to_ptl (grid, psn, spall)
 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)
 
subroutine cce_fcoupling_resample (grid, psn, spall, istep)
 
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 ()
 

Public Attributes

integer resamp_rate = 2
 timesteps between subsequent resamples (multiplied by sml_f_source_period More...
 
real(8) resamp_var = 0.02D0
 threshold for relative standard deviation in bin for auto-resample More...
 
real(8) resamp_min_ratio = 0.2D0
 min ratio of (# of ptl)/(target # of ptl) in bin for auto-upsample More...
 
real(8) resamp_max_ratio = 2.0D0
 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-7
 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...
 
character(len=200) resamp_node_file
 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...
 
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
 

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 for normalization. More...
 
real(kind=8), dimension(:,:),
allocatable, private 
resamp_t_ev
 Reference temperature 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...
 
integer, parameter, private resamp_nongridmom = 8
 Number of bin moments to be conserved. 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, 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...
 
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...
 

Member Function/Subroutine Documentation

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.

Parameters
[in,out]binPhase-space bin data structure, type(resamp_bin_type)
[in]ptliData of a single particle, type(ptl_type)
[in]indIndex of that particle in the corresponding species data structure, integer

Here is the caller graph for this function:

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.

Parameters
[in,out]spParticle data structure for one species, type(species_type)
[in,out]binsPhase-space bin data structure, type(resamp_bin_type)
[in]gridXGC grid data structure, type(grid_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::cce_buffer_rcv_updated_bins ( integer, intent(in)  isp)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::cce_buffer_snd_updated_bins ( integer, intent(in)  isp)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::cce_fcoupling_resample ( type(grid_type), pointer  grid,
type(psn_type), pointer  psn,
type(species_type), dimension(:), pointer  spall,
integer, intent(in)  istep 
)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::cce_reset_particle_outside ( type(grid_type), intent(in)  grid,
type(species_type), intent(out)  sp 
)

Here is the call graph for this function:

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 
)

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.

Parameters
[in,out]binPhase-space bin data structure, type(resamp_bin_type)
[in]spParticle data, type(species_type)
[in]cce_sync_sendrcvfor the coupling
[in,out]resamp_auto_upsample_locCounter for upsampled bins
[in,out]resamp_auto_downsample_locCounter for downsampled bins
[in,out]resamp_var_resample_locCounter for variance-resampled bins

We must pass all resamp counters to be updated correctly with OMP.

Here is the caller graph for this function:

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)

Parameters
[in,out]binsBin data structure, type(resamp_bin_type)

Here is the caller graph for this function:

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.

Parameters
[in]gridXGC configuration space grid, type(grid_type)
[in]psnE-field/potential data structure, type(psn_type)
[in]ptliParticle data (phase+constants), type(ptl_type)
[in]ispParticle species index, integer
[out]velodX/dt and drho/dt, real(8)
[out]vradRadial particle velocity (dX/dt).grad(psi)

Here is the call graph for this function:

Here is the caller graph for this function:

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 
)

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\).

Parameters
[in]gridXGC grid data structure, type(grid_type)
[in]psnField data, type(psn_type)
[in,out]partInput particle data, type(ptl_type)
[in]psiPoloidal flux at particle position, real(8)
[in]BMagnetic field |B_0|, real(8)
[in]flowBackground mean parallel flow, real(8)
[in]denBackground density, real(8)
[in]temp_evBackground temperature
[in]specSpecies type of input particle, integer
[in]triTriangle in which the particle is, integer
[in]baryCorresponding barycentric coordinates, real(8)
[in]mark_denInverse of particle's phase space volume g=1/dV, real(8)
[out]errError flag, logical

Here is the call graph for this function:

Here is the caller graph for this function:

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.

Parameters
[in]gridXGC configuration space grid, type(grid_type)

Here is the call graph for this function:

Here is the caller graph for this function:

logical function resamp_module::is_first_call ( )

Here is the caller graph for this function:

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.

Parameters
[in]gridXGC configuration space grid, type(grid_type)
[in,out]new_ptlParticle to store data in
[in]bin_idID of the bin into which particle is loaded
[out]itrTriangle into which the particle was loaded, integer
[out]pBarycentric coordinates of the loaded triangle, real(8)
[out]psiPol. flux at new particle's location, real(8)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::load_new_ptl2 ( 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 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.

Parameters
[in]gridXGC configuration space grid, type(grid_type)
[in,out]new_ptlParticle to store data in
[in]bin_idID of the bin into which particle is loaded
[out]itrTriangle into which the particle was loaded, integer
[out]pBarycentric coordinates of the loaded triangle, real(8)
[out]psiPol. flux at new particle's location, real(8)

Here is the call graph for this function:

subroutine resamp_module::move_f0g_to_ptl ( type(grid_type), intent(in)  grid,
type(psn_type), intent(in)  psn,
type(species_type), dimension(0:ptl_nsp_max), intent(inout)  spall 
)

Move phase space density from the 5D grid over to the particles for restarting simulatin with a different grid.

Parameters
[in]gridXGC grid data structure, type(grid_type)
[in]psnField data, type(psn_type)
[in,out]spallParticle data structure for all species, type(species_type)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::resamp_scale_phvol ( type(grid_type), intent(in)  grid,
type(species_type), intent(inout)  sp 
)

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::resample ( type (grid_type), intent(in)  grid,
type (species_type), dimension(0:ptl_nsp_max), intent(inout)  spall,
type (psn_type), intent(in)  psn,
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)  cce_sync_sendrcv,
integer, optional  tile_size 
)

Resample the marker particles in XGC to improve phase-space coverage and reduce marker noise.

Parameters
[in]gridXGC grid data structure, type(grid_type)
[in,out]spallParticle data, type(species_type)
[in]psnField data, type(psn_type)
[in]inode1Index of first mesh vertex in domain, integer
[in]inode2Index of last mesh vertex in domain, integer
[in]dsmuv_perp/mu resolution of XGC's total-f velocity grid, integer
[in]dvpv_para resolution of XGC's total-f v-grid, integer
[in]nmuNumber of v_perp/mu bins in total-f v-grid, integer
[in]nvpHalf the number of v_para bins in total-f v-grid, integer
[in]vp_maxv_para cutoff of total-f v-grid, real(8)
[in]mu_maxv_perp/mu cutoff of total-f v-grid, real(8)
[in]tile_sizeResampling combines tile_size x tile_size bins of the total-f v-grid to one macro bin, integer, optional
[in]cce_sendrcvTakes values 0(nothing), 1(sender) , 2(receiver)

Here is the call graph for this function:

Here is the caller graph for this function:

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 
)

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.

Parameters
[in,out]binPhase-space bin data structure, type(resamp_bin_type)
[in]gridXGC grid data structure, type(grid_type)
[in,out]spParticle data structure, type(species_type)
[in]psnField data, type(psn_type)
[in]cce_sync_sendrcv0:no f-coupling, 1 sender, 2 receiver
[in]start_wr_posWrite position in sp, integer
[in,out]resamp_failed_bins_locCounter for failed bins
[in,out]resamp_upsample_failed_locCounter for failed upsampled bins
[in,out]resamp_downsample_failed_locCounter for failed downsampled bins
[in,out]resamp_resample_failed_locCounter for failed var-resampled bins
[in,out]resamp_retried_bins_locCounter for retried bins
[in,out]resamp_retried_failed_locCounter for bins that failed on retry
[in,out]resamp_fullf_failed_locCounter for bins with failed full-f resampling
[in,out]resamp_phvol_added_locPhase space volume of particles added to empty bins

We must pass all resamp counters to be updated correctly with OMP.

Here is the call graph for this function:

Here is the caller graph for this function:

subroutine resamp_module::resample_finalize ( type(resamp_bin_type), dimension(:,:,:), allocatable  bins)

This routine cleans up the resampling module.

Parameters
[in,out]binsResampling bin data structure, type(resamp_bin_type)

Here is the caller graph for this function:

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.

Parameters
[in]gridXGC mesh data structure, type(grid_type)
[in]inode1Index of first vertex in the mesh domain, integer
[in]inode2Index of last vertex in the mesh domain, integer
[in]dsmuv_perp/mu resolution of XGC's total-f velocity grid, integer
[in]dvpv_para resolution of XGC's total-f v-grid, integer
[in]nmuNumber of v_perp/mu bins in total-f v-grid, integer
[in]nvpHalf the number of v_para bins in total-f v-grid, integer
[in]vp_maxv_para cutoff of total-f v-grid, real(8)
[in]mu_maxv_perp/mu cutoff of total-f v-grid, real(8)
[in]tile_sizeResampling combines tile_size x tile_size bins of the total-f v-grid to one macro bin, integer
[in,out]binsResampling phase space bins, type(resamp_bin_type)

Here is the call graph for this function:

Here is the caller graph for this function:

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 
)

Here is the caller graph for this function:

subroutine resamp_module::set_cce_buffer_id ( integer, intent(in)  id,
integer, intent(in)  inode1,
integer, intent(in)  inode2 
)

Here is the caller graph for this function:

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(out)  psi 
)

Evaluate the background density, temperature, flow, magnetic field and psi at the input particle's position.

Parameters
[in]ptliInput particle data, type(ptl_type)
in[spec_type Species type of input particle, integer
[out]denBackground density, real(8)
[out]flowBackground toroidal flow velocity, real(8)
[out]temp_evBackground temperature, real(8)
[out]bvecMagnetic field vector B, real(8)
[out]BMagnetic field |B_0|, real(8)
[out]psiPoloidal magnetic flux, real(8)

Here is the call graph for this function:

Here is the caller graph for this function:

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.

Parameters
[in]nconstraintNumber of constraints, integer
[in]part_numNumber of particles, integer
[in]CmatConstraint matrix, N_constraint x N_ptl, real(8)
[in]constraintValues of the constraints, real(8)
[in,out]weightsDelta-f weights w1, N_ptl, real(8)
[in]grid_ineq_onWhether to conserve the grid charge or not, logical
[in]ineq_tolError tolerance for the inequality constraints, real(8)
[in]positivityWhether to require weight>0, logical
[in,out]ierrError flag, integer

Here is the call graph for this function:

Here is the caller graph for this function:

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.

Parameters
[in]gridXGC grid data structure, type(grid_type)
[in]psnXGC E-field/potential data, type(psn_type)
[in]nconstraintNumber of equality+inequality constraints, integer
[in]nptlNumber of particles in bin, integer
[in]vertNumber of vertices connected to this bin (patch size), integer
[in]patchIndices of the vertices connected to this bin, integer
[in]baryBarycentric coordinates of the particle, real(8)
[in]triTriangle in which the particle is, integer
[in]psiPol. flux at particle position, real(8)
[in]BLocal magnetic field |B|, real(8)
[in]bvecMagnetic field vector B, real(8)
[in]specSpecies index, integer
[in]pindexPosition of the particle in the particle set, integer
[in]partData of one particle, type(ptl_type)
[in,out]CmatThe matrix C from the (see documentation), real(8)
[in]op_modeSwitch for optimizing full-f weights w0 (0)or delta-f weights (1), integer

Here is the call graph for this function:

Here is the caller graph for this function:

Member Data Documentation

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
integer, private resamp_module::resamp_auto_downsample
private

Number of downsampled bins.

integer, private resamp_module::resamp_auto_upsample
private

Number of upsampled bins.

real (kind=8), private resamp_module::resamp_bin_max_fac
private

threshold for auto downsample

real (kind=8), dimension(:), allocatable, private resamp_module::resamp_den
private

Reference density for normalization.

integer, private resamp_module::resamp_downsample_failed
private

Number of bins for which downsampling failed.

real (kind=8), private resamp_module::resamp_dsmu
private

v_perp/mu spacing of original v-grid (=f0_dsmu)

real (kind=8), private resamp_module::resamp_dvp
private

v_para spacing of original v-grid (=f0_dvp)

integer, private resamp_module::resamp_empty_bins
private

Number of empty sampling bins.

integer, private resamp_module::resamp_failed_bins
private

Number of bins for which re-/up-/down-sampling failed.

logical resamp_module::resamp_fill_empty =.false.

Whether to fill empty bins.

integer, private resamp_module::resamp_fullf_failed
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.

real (8) resamp_module::resamp_ineq_tol = 1D-7

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.

integer, private resamp_module::resamp_inode1
private

Index of the first node of local patch of the configuration space mesh (=f0_inode1)

integer, private resamp_module::resamp_inode2
private

Index of the last node of local patch of the configuration space mesh (=f0_inode2)

real (kind=8), private resamp_module::resamp_inv_bin_num
private

Inverse of the number of macro bins in v-space.

real (kind=8), dimension(:,:), allocatable, private resamp_module::resamp_inv_ph_vol
private

Inverse phase-space volume element.

real (8) resamp_module::resamp_max_ratio = 2.0D0

max ratio of (# of ptl)/(target # of ptl) in bin for auto-downsample

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.2D0

min ratio of (# of ptl)/(target # of ptl) in bin for auto-upsample

real (kind=8), private resamp_module::resamp_mu_max
private

v_perp/mu cutoff of v-grid (=f0_mu_max)

integer, private resamp_module::resamp_musize
private

Number of sampling bins in v_perp/mu direction.

integer, private resamp_module::resamp_nmu_in
private

Number of cells on original vperp/mu grid (=f0_nmu)

character (len=200) resamp_module::resamp_node_file

File containing the vertex positions of the mesh for which to resample.

integer, parameter, private resamp_module::resamp_nongridmom = 8
private

Number of bin moments to be conserved.

integer resamp_module::resamp_nphi_new =1

Number of poloidal planes in simulation with new mesh.

integer, private resamp_module::resamp_nvp_in
private

Number of cells on original v_para grid (=f0_nvp)

integer, dimension(:), allocatable, private resamp_module::resamp_off_vgrid
private

Counter for number of particles that are not on the v-grid.

integer, dimension(:), allocatable, private resamp_module::resamp_on_vgrid
private

Counter for number of particles that are on the v-grid.

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

real (kind=8), dimension(:), allocatable, private resamp_module::resamp_phvol_added
private

Sum of phase space volume of particles added in empty bins.

real (kind=8), dimension(:), allocatable, private resamp_module::resamp_phvol_sum
private

Sum of phase space volume of all particles on a node.

integer resamp_module::resamp_rate = 2

timesteps between subsequent resamples (multiplied by sml_f_source_period

integer, private resamp_module::resamp_resample_failed
private

Number of bins for which variance resampling failed.

integer, private resamp_module::resamp_resampled_bins
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.

logical, private resamp_module::resamp_restart_mode =.false.
private
logical resamp_module::resamp_restart_read =.false.

Whether to read a restart file written from a simulation with different grid.

integer, private resamp_module::resamp_retried_bins
private

Number of bins retried after failure.

integer, private resamp_module::resamp_retried_failed
private

Number of bins failed again after retry.

logical resamp_module::resamp_retry = .false.

Retry QP optimization for failed bins with relaxed inequality constraints.

real (kind=8), dimension(:), allocatable, private resamp_module::resamp_samp_ratio
private

Upsampling ratio for restarts with finer grid.

real (kind=8), private resamp_module::resamp_smu_fac
private

Factor used in the context of tapering off the marker density at high velocity.

real (kind=8), dimension(:,:), allocatable, private resamp_module::resamp_t_ev
private

Reference temperature 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.

integer, private resamp_module::resamp_upsample_failed
private

Number of bins for which upsampling failed.

real (8) resamp_module::resamp_var = 0.02D0

threshold for relative standard deviation in bin for auto-resample

integer, private resamp_module::resamp_var_resample
private

Number of variance resampled bins.

real (kind=8), private resamp_module::resamp_vp_fac
private

Factor used in the context of tapering off the marker density at high velocity.

real (kind=8), private resamp_module::resamp_vp_max
private

v_para cutoff of v-grid (=f0_vp_max)

integer, private resamp_module::resamp_vpsize
private

Number of sampling bins in v_para direction.

integer, private resamp_module::resamp_wrong_particles
private

Number of particles in wwrong bins.


The documentation for this module was generated from the following file: