fida_mc Subroutine

  1. fidasim.f90
  2. libfida
  3. fida_mc

public subroutine fida_mc()

Calculate FIDA emission using a Monte Carlo Fast-ion distribution

Arguments

None

Calls

proc~~fida_mc~~CallsGraph proc~fida_mc fida_mc proc~get_fields get_fields proc~fida_mc->proc~get_fields proc~track track proc~fida_mc->proc~track proc~gyro_correction gyro_correction proc~fida_mc->proc~gyro_correction proc~randu randu proc~fida_mc->proc~randu proc~uvw_to_xyz uvw_to_xyz proc~fida_mc->proc~uvw_to_xyz proc~in_plasma in_plasma proc~get_fields->proc~in_plasma proc~calc_perp_vectors calc_perp_vectors proc~get_fields->proc~calc_perp_vectors proc~track->proc~in_plasma proc~get_indices get_indices proc~track->proc~get_indices proc~gyro_correction->proc~randu proc~gyro_step gyro_step proc~gyro_correction->proc~gyro_step proc~pitch_to_vec pitch_to_vec proc~gyro_correction->proc~pitch_to_vec proc~rng_uniform rng_uniform proc~randu->proc~rng_uniform omp_get_thread_num omp_get_thread_num proc~randu->omp_get_thread_num interface~interpol_coeff interpol_coeff proc~in_plasma->interface~interpol_coeff proc~xyz_to_uvw xyz_to_uvw proc~in_plasma->proc~xyz_to_uvw proc~cross_product cross_product proc~gyro_step->proc~cross_product proc~interpol2d_coeff interpol2D_coeff interface~interpol_coeff->proc~interpol2d_coeff proc~interpol1d_coeff interpol1D_coeff interface~interpol_coeff->proc~interpol1d_coeff proc~interpol2d_coeff_arr interpol2D_coeff_arr interface~interpol_coeff->proc~interpol2d_coeff_arr proc~interpol1d_coeff_arr interpol1D_coeff_arr interface~interpol_coeff->proc~interpol1d_coeff_arr proc~interpol2d_coeff_arr->proc~interpol2d_coeff proc~interpol1d_coeff_arr->proc~interpol1d_coeff

Called by

proc~~fida_mc~~CalledByGraph proc~fida_mc fida_mc program~fidasim fidasim program~fidasim->proc~fida_mc

Contents

Source Code

fida_mc

Source Code

subroutine fida_mc
    !+ Calculate FIDA emission using a Monte Carlo Fast-ion distribution
    integer :: iion,iphi,nphi
    type(FastIon) :: fast_ion
    type(LocalEMFields) :: fields
    type(LocalProfiles) :: plasma
    real(Float64) :: phi
    real(Float64), dimension(3) :: ri      !! start position
    real(Float64), dimension(3) :: vi      !! velocity of fast ions
    !! Determination of the CX probability
    real(Float64), dimension(nlevs) :: denn    !!  neutral dens (n=1-4)
    real(Float64), dimension(nlevs) :: rates    !! CX rates
    !! Collisiional radiative model along track
    real(Float64), dimension(nlevs) :: states  ! Density of n-states
    integer :: ncell
    type(ParticleTrack), dimension(beam_grid%ntrack) :: tracks
    logical :: los_intersect
    integer :: jj      !! counter along track
    real(Float64) :: photons !! photon flux
    integer, dimension(4) :: neut_types=[1,2,3,4]
    real(Float64), dimension(3) :: uvw, uvw_vi
    real(Float64)  :: s, c
    real(Float64)  :: maxcnt, inv_maxcnt, cnt
    real(Float64), dimension(1) :: randomu

    maxcnt=particles%nparticle
    inv_maxcnt = 100.d0/maxcnt
    nphi = ceiling(dble(inputs%n_fida)/particles%nparticle)
    if(inputs%verbose.ge.1) then
        write(*,'(T6,"# of markers: ",i9)') int(particles%nparticle*nphi,Int64)
    endif

    cnt=0.0
    !$OMP PARALLEL DO schedule(guided) private(iion,iphi,fast_ion,vi,ri,phi,tracks,s,c, &
    !$OMP& randomu,plasma,fields,uvw,uvw_vi,ncell,jj,rates,denn,los_intersect,states,photons)
    loop_over_fast_ions: do iion=1,particles%nparticle
        fast_ion = particles%fast_ion(iion)
        cnt=cnt+1
        if(fast_ion%vabs.eq.0) cycle loop_over_fast_ions
        if(.not.fast_ion%cross_grid) cycle loop_over_fast_ions
        phi_loop: do iphi=1,nphi
            !! Pick random torodial angle
            call randu(randomu)
            phi = fast_ion%phi_enter + fast_ion%delta_phi*randomu(1)
            s = sin(phi) ; c = cos(phi)

            !! Calculate position in machine coordinates
            uvw(1) = fast_ion%r*c
            uvw(2) = fast_ion%r*s
            uvw(3) = fast_ion%z

            !! Convert to beam grid coordinates
            call uvw_to_xyz(uvw, ri)

            if(inputs%dist_type.eq.2) then
                !! Get electomagnetic fields
                call get_fields(fields, pos=ri)

                !! Correct for gyro motion and get particle position and velocity
                call gyro_correction(fields, fast_ion%energy, fast_ion%pitch, ri, vi)
            else !! Full Orbit
                !! Calculate velocity vector
                uvw_vi(1) = c*fast_ion%vr - s*fast_ion%vt
                uvw_vi(2) = s*fast_ion%vr + c*fast_ion%vt
                uvw_vi(3) = fast_ion%vz
                vi = matmul(beam_grid%inv_basis,uvw_vi)
            endif

            !! Track particle through grid
            call track(ri, vi, tracks, ncell, los_intersect)
            if(.not.los_intersect) cycle phi_loop
            if(ncell.eq.0)cycle phi_loop

            !! Calculate CX probability
            call get_beam_cx_rate(tracks(1)%ind,ri,vi,beam_ion,neut_types,rates)
            if(sum(rates).le.0.)cycle phi_loop

            !! Weight CX rates by ion source density
            states=rates*fast_ion%weight/nphi

            !! Calculate the spectra produced in each cell along the path
            loop_along_track: do jj=1,ncell
                call get_plasma(plasma,pos=tracks(jj)%pos)

                call colrad(plasma,beam_ion, vi, tracks(jj)%time, states, denn, photons)

                call store_fida_photons(tracks(jj)%pos, vi, photons, fast_ion%class)
            enddo loop_along_track
        enddo phi_loop
        if (inputs%verbose.ge.2)then
          WRITE(*,'(f7.2,"% completed",a,$)') cnt*inv_maxcnt,char(13)
        endif
    enddo loop_over_fast_ions
    !$OMP END PARALLEL DO

end subroutine fida_mc