Reads in the fast-ion distribution function and stores the quantities in fbm
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| integer(kind=HID_T), | intent(inout) | :: | fid | HDF5 file ID |
||
| integer, | intent(out) | :: | error | Error code |
subroutine read_f(fid, error)
!+ Reads in the fast-ion distribution function and stores the quantities in [[libfida:fbm]]
integer(HID_T), intent(inout) :: fid
!+ HDF5 file ID
integer, intent(out) :: error
!+ Error code
integer(HSIZE_T), dimension(5) :: dims
real(Float64) :: deni_tot
integer :: ir,is
logical :: path_valid
if(inputs%verbose.ge.1) then
write(*,'(a)') '---- Fast-ion distribution settings ----'
endif
call h5ltread_dataset_int_scalar_f(fid,"/nenergy", fbm%nenergy, error)
call h5ltread_dataset_int_scalar_f(fid,"/npitch", fbm%npitch, error)
call h5ltread_dataset_int_scalar_f(fid,"/nr", fbm%nr, error)
call h5ltread_dataset_int_scalar_f(fid,"/nz", fbm%nz, error)
call h5ltpath_valid_f(fid, "/nphi", .True., path_valid, error)
if(path_valid) then
call h5ltread_dataset_int_scalar_f(fid,"/nphi", fbm%nphi, error)
else
fbm%nphi=1
endif
if(((fbm%nr.ne.inter_grid%nr).or.(fbm%nz.ne.inter_grid%nz)).or.(fbm%nphi.ne.inter_grid%nphi)) then
if(inputs%verbose.ge.0) then
write(*,'(a)') "READ_F: Distribution file has incompatable grid dimensions"
endif
stop
endif
if (fbm%nphi .eq. 1) then
allocate(fbm%energy(fbm%nenergy), fbm%pitch(fbm%npitch), fbm%r(fbm%nr), fbm%z(fbm%nz), fbm%phi(1))
allocate(fbm%denf(fbm%nr, fbm%nz,1))
allocate(fbm%f(fbm%nenergy, fbm%npitch, fbm%nr, fbm%nz,1))
else
allocate(fbm%energy(fbm%nenergy), fbm%pitch(fbm%npitch), fbm%r(fbm%nr), fbm%z(fbm%nz), fbm%phi(fbm%nphi))
allocate(fbm%denf(fbm%nr, fbm%nz, fbm%nphi))
allocate(fbm%f(fbm%nenergy, fbm%npitch, fbm%nr, fbm%nz, fbm%nphi))
endif
if (fbm%nphi .eq. 1) then
dims = [fbm%nenergy, fbm%npitch, fbm%nr, fbm%nz, 1]
call h5ltread_dataset_double_f(fid, "/denf",fbm%denf, dims(3:4), error)
call h5ltread_dataset_double_f(fid, "/f", fbm%f, dims(1:4), error)
else
dims = [fbm%nenergy, fbm%npitch, fbm%nr, fbm%nz, fbm%nphi]
call h5ltread_dataset_double_f(fid, "/denf",fbm%denf, dims(3:5), error)
call h5ltread_dataset_double_f(fid, "/f", fbm%f, dims(1:5), error)
endif
call h5ltread_dataset_double_f(fid, "/energy", fbm%energy, dims(1:1), error)
call h5ltread_dataset_double_f(fid, "/pitch", fbm%pitch, dims(2:2), error)
call h5ltread_dataset_double_f(fid, "/r", fbm%r, dims(3:3), error)
call h5ltread_dataset_double_f(fid, "/z", fbm%z, dims(4:4), error)
call h5ltpath_valid_f(fid, "/phi", .True., path_valid, error)
if(path_valid) then
call h5ltread_dataset_double_f(fid, "/phi", fbm%phi, dims(5:5), error)
else
fbm%phi=0.d0
endif
!call h5ltread_dataset_double_scalar_f(fid,"/A",fbm%A, error)
fbm%A = beam_mass
equil%plasma%denf = fbm%denf
fbm%dE = abs(fbm%energy(2) - fbm%energy(1))
fbm%dp = abs(fbm%pitch(2) - fbm%pitch(1))
fbm%dr = abs(fbm%r(2) - fbm%r(1))
fbm%dz = abs(fbm%z(2) - fbm%z(1))
if (fbm%nphi .eq. 1) then
fbm%dphi = 2*pi
else
fbm%dphi = abs(fbm%phi(2)-fbm%phi(1))
endif
fbm%emin = minval(fbm%energy,1)
fbm%emax = maxval(fbm%energy,1)
fbm%e_range = fbm%emax - fbm%emin
fbm%pmin = minval(fbm%pitch,1)
fbm%pmax = maxval(fbm%pitch,1)
fbm%p_range = fbm%pmax - fbm%pmin
fbm%rmin = minval(fbm%r,1)
fbm%rmax = maxval(fbm%r,1)
fbm%r_range = fbm%rmax - fbm%rmin
fbm%zmin = minval(fbm%z,1)
fbm%zmax = maxval(fbm%z,1)
fbm%z_range = fbm%zmax - fbm%zmin
fbm%phimin = minval(fbm%phi,1)
fbm%phimax = maxval(fbm%phi,1)
fbm%phi_range = fbm%phimax - fbm%phimin
deni_tot = 0.0
do ir=1,fbm%nr
fbm%n_tot = fbm%n_tot + fbm%dphi*fbm%dr*fbm%dz*sum(fbm%denf(ir,:,:))*fbm%r(ir)
do is=1,n_thermal
deni_tot = deni_tot + fbm%dphi*fbm%dr*fbm%dz*sum(equil%plasma(ir,:,:)%deni(is))*fbm%r(ir)
enddo
enddo
if(fbm%n_tot.ge.deni_tot) then
if(inputs%verbose.ge.0) then
write(*,'(a," (",ES10.3," >=",ES10.3,")")') &
"READ_F: The total of number of fast ions exceeded the total number of thermal ions.", &
fbm%n_tot, deni_tot
write(*,'(a)') "This is usually caused by zeff being incorrect."
endif
stop
endif
if(inputs%verbose.ge.1) then
if(fbm%nphi.gt.1) then
write(*,'(T2,"Distribution type: ",a)') "Non-axisymmetric Fast-ion Density Function F(energy,pitch,R,Z,Phi)"
else
write(*,'(T2,"Distribution type: ",a)') "Axisymmetric Fast-ion Density Function F(energy,pitch,R,Z)"
endif
write(*,'(T2,"Nenergy = ",i3)') fbm%nenergy
write(*,'(T2,"Npitch = ",i3)') fbm%npitch
write(*,'(T2,"Nr = ",i3)') fbm%nr
write(*,'(T2,"Nz = ",i3)') fbm%nz
if(fbm%nphi.gt.1) then
write(*,'(T2,"Nphi = ",i3)') fbm%nphi
endif
write(*,'(T2,"Energy range = [",f5.2,",",f6.2,"]")') fbm%emin,fbm%emax
write(*,'(T2,"Pitch range = [",f5.2,",",f5.2,"]")') fbm%pmin,fbm%pmax
write(*,'(T2,"R range = [",f6.2,",",f6.2,"]")') fbm%rmin,fbm%rmax
write(*,'(T2,"Z range = [",f7.2,",",f6.2,"]")') fbm%zmin,fbm%zmax
if(fbm%nphi.gt.1) then
write(*,'(T2,"Phi range = [",f5.2,",",f5.2,"]")') fbm%phimin,fbm%phimax
endif
write(*,'(T2,"Ntotal = ",ES10.3)') fbm%n_tot
write(*,*) ''
endif
end subroutine read_f