Samples a Guiding Center Fast-ion distribution function at a given pass_grid index
Type | Intent | Optional | Attributes | Name | ||
---|---|---|---|---|---|---|
integer, | intent(in), | dimension(3) | :: | ind | pass_grid index |
|
type(LocalEMFields), | intent(out) | :: | fields | Electromagnetic fields at the guiding center |
||
real(kind=Float64), | intent(out) | :: | eb | Energy of the fast ion |
||
real(kind=Float64), | intent(out) | :: | ptch | Pitch of the fast ion |
||
real(kind=Float64), | intent(out) | :: | denf | Fast-ion density at guiding center |
||
integer, | intent(in), | optional | :: | output_coords | Indicates coordinate system of |
subroutine mc_fastion_pass_grid(ind,fields,eb,ptch,denf,output_coords)
!+ Samples a Guiding Center Fast-ion distribution function at a given [[libfida:pass_grid]] index
integer, dimension(3), intent(in) :: ind
!+ [[libfida:pass_grid]] index
type(LocalEMFields), intent(out) :: fields
!+ Electromagnetic fields at the guiding center
real(Float64), intent(out) :: eb
!+ Energy of the fast ion
real(Float64), intent(out) :: ptch
!+ Pitch of the fast ion
real(Float64), intent(out) :: denf
!+ Fast-ion density at guiding center
integer, intent(in), optional :: output_coords
!+ Indicates coordinate system of `fields`. Beam grid (0), machine (1) and cylindrical (2)
real(Float64), dimension(fbm%nenergy,fbm%npitch) :: fbeam
real(Float64), dimension(3) :: rg, rg_cyl
real(Float64), dimension(3) :: randomu3
real(Float64) :: rmin, rmax, zmin, phimin
integer, dimension(2,1) :: ep_ind
integer :: ocs
if(present(output_coords)) then
ocs = output_coords
else
ocs = 0
endif
denf=0.d0
call randu(randomu3)
rmin = pass_grid%r(ind(1))
rmax = rmin + pass_grid%dr
zmin = pass_grid%z(ind(2))
phimin = pass_grid%phi(ind(3))
! Sample uniformally in annulus
rg_cyl(1) = sqrt(randomu3(1)*(rmax**2 - rmin**2) + rmin**2)
rg_cyl(2) = zmin + randomu3(2)*pass_grid%dz
rg_cyl(3) = phimin + randomu3(3)*pass_grid%dphi
call cyl_to_uvw(rg_cyl, rg)
call get_fields(fields,pos=rg,input_coords=1,output_coords=ocs)
if(.not.fields%in_plasma) return
call get_distribution(fbeam,denf,coeffs=fields%b)
call randind(fbeam,ep_ind)
call randu(randomu3)
eb = fbm%energy(ep_ind(1,1)) + fbm%dE*(randomu3(1)-0.5)
ptch = fbm%pitch(ep_ind(2,1)) + fbm%dp*(randomu3(2)-0.5)
end subroutine mc_fastion_pass_grid