dcx Subroutine

public subroutine dcx()

Calculates Direct Charge Exchange (DCX) neutral density and spectra

Arguments

None

Calls

proc~~dcx~~CallsGraph proc~dcx dcx proc~mc_halo mc_halo proc~dcx->proc~mc_halo proc~get_nlaunch get_nlaunch proc~dcx->proc~get_nlaunch proc~ind2sub ind2sub proc~dcx->proc~ind2sub proc~get_plasma get_plasma proc~dcx->proc~get_plasma proc~track track proc~dcx->proc~track interface~parallel_sum parallel_sum proc~dcx->interface~parallel_sum proc~mc_halo->proc~get_plasma interface~randn randn proc~mc_halo->interface~randn interface~randu randu proc~mc_halo->interface~randu proc~rng_init rng_init proc~get_nlaunch->proc~rng_init interface~randind_cdf randind_cdf proc~get_nlaunch->interface~randind_cdf proc~cumsum cumsum proc~get_nlaunch->proc~cumsum proc~uvw_to_xyz uvw_to_xyz proc~get_plasma->proc~uvw_to_xyz proc~in_plasma in_plasma proc~get_plasma->proc~in_plasma proc~get_position get_position proc~get_plasma->proc~get_position proc~xyz_to_uvw xyz_to_uvw proc~get_plasma->proc~xyz_to_uvw proc~get_fields get_fields proc~track->proc~get_fields proc~doppler_stark doppler_stark proc~track->proc~doppler_stark proc~track->proc~in_plasma proc~get_indices get_indices proc~track->proc~get_indices proc~parallel_sum_i1 parallel_sum_i1 interface~parallel_sum->proc~parallel_sum_i1 proc~parallel_sum_i0 parallel_sum_i0 interface~parallel_sum->proc~parallel_sum_i0 proc~parallel_sum_d4 parallel_sum_d4 interface~parallel_sum->proc~parallel_sum_d4 proc~parallel_sum_d5 parallel_sum_d5 interface~parallel_sum->proc~parallel_sum_d5 proc~parallel_sum_d2 parallel_sum_d2 interface~parallel_sum->proc~parallel_sum_d2 proc~parallel_sum_d3 parallel_sum_d3 interface~parallel_sum->proc~parallel_sum_d3 proc~parallel_sum_d0 parallel_sum_d0 interface~parallel_sum->proc~parallel_sum_d0 proc~parallel_sum_d1 parallel_sum_d1 interface~parallel_sum->proc~parallel_sum_d1 proc~parallel_sum_i2 parallel_sum_i2 interface~parallel_sum->proc~parallel_sum_i2 proc~get_fields->proc~uvw_to_xyz proc~get_fields->proc~in_plasma proc~get_fields->proc~xyz_to_uvw proc~calc_perp_vectors calc_perp_vectors proc~get_fields->proc~calc_perp_vectors proc~rng_seed rng_seed proc~rng_init->proc~rng_seed proc~my_rank my_rank proc~rng_init->proc~my_rank mpi_allreduce mpi_allreduce proc~parallel_sum_i1->mpi_allreduce proc~parallel_sum_i0->mpi_allreduce proc~parallel_sum_d4->mpi_allreduce proc~in_plasma->proc~xyz_to_uvw interface~interpol_coeff interpol_coeff proc~in_plasma->interface~interpol_coeff proc~cyl_to_uvw cyl_to_uvw proc~in_plasma->proc~cyl_to_uvw proc~cyl_to_xyz cyl_to_xyz proc~get_position->proc~cyl_to_xyz proc~parallel_sum_d5->mpi_allreduce proc~parallel_sum_d2->mpi_allreduce proc~parallel_sum_d3->mpi_allreduce proc~parallel_sum_d0->mpi_allreduce proc~parallel_sum_d1->mpi_allreduce proc~parallel_sum_i2->mpi_allreduce proc~cyl_interpol3d_coeff cyl_interpol3D_coeff interface~interpol_coeff->proc~cyl_interpol3d_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~interpol2d_coeff interpol2D_coeff interface~interpol_coeff->proc~interpol2d_coeff proc~cyl_interpol3d_coeff_arr cyl_interpol3D_coeff_arr interface~interpol_coeff->proc~cyl_interpol3d_coeff_arr proc~interpol1d_coeff_arr interpol1D_coeff_arr interface~interpol_coeff->proc~interpol1d_coeff_arr proc~cyl_to_xyz->proc~uvw_to_xyz proc~cyl_to_xyz->proc~cyl_to_uvw proc~cyl_interpol3d_coeff->proc~interpol2d_coeff proc~interpol2d_coeff_arr->proc~interpol2d_coeff proc~cyl_interpol3d_coeff_arr->proc~cyl_interpol3d_coeff proc~cyl_interpol3d_coeff_arr->proc~interpol2d_coeff proc~interpol1d_coeff_arr->proc~interpol1d_coeff

Called by

proc~~dcx~~CalledByGraph proc~dcx dcx program~fidasim fidasim program~fidasim->proc~dcx

Contents

Source Code

dcx

Source Code

subroutine dcx
    !+ Calculates Direct Charge Exchange (DCX) neutral density and spectra
    integer :: ic,i,j,k,ncell
    integer(Int64) :: idcx !! counter
    real(Float64), dimension(3) :: ri    !! start position
    real(Float64), dimension(3) :: vihalo
    integer,dimension(3) :: ind
    integer,dimension(3) :: neut_types = [1,2,3]
    !! Determination of the CX probability
    type(LocalProfiles) :: plasma
    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 :: ntrack
    type(ParticleTrack), dimension(beam_grid%ntrack) :: tracks  !! Particle tracks
    integer :: jj       !! counter along track
    real(Float64):: tot_denn, photons  !! photon flux
    integer, dimension(beam_grid%ngrid) :: cell_ind
    real(Float64), dimension(beam_grid%nx,beam_grid%ny,beam_grid%nz) :: papprox
    integer(Int32), dimension(beam_grid%nx,beam_grid%ny,beam_grid%nz) :: nlaunch
    real(Float64) :: fi_correction

    halo_iter_dens(dcx_type) = 0.d0
    papprox=0.d0
    tot_denn=0.d0
    do ic=1,beam_grid%ngrid
        call ind2sub(beam_grid%dims,ic,ind)
        i = ind(1) ; j = ind(2) ; k = ind(3)
        call get_plasma(plasma,ind=ind)
        if(.not.plasma%in_plasma) cycle
        tot_denn = sum(neut%full(:,i,j,k)) + &
                   sum(neut%half(:,i,j,k)) + &
                   sum(neut%third(:,i,j,k))
        papprox(i,j,k)= tot_denn*(plasma%denp-plasma%denf)
    enddo

    ncell = 0
    do ic=1,beam_grid%ngrid
        call ind2sub(beam_grid%dims,ic,ind)
        i = ind(1) ; j = ind(2) ; k = ind(3)
        if(papprox(i,j,k).gt.0.0) then
            ncell = ncell + 1
            cell_ind(ncell) = ic
        endif
    enddo

    call get_nlaunch(inputs%n_dcx,papprox,nlaunch)

    if(inputs%verbose.ge.1) then
       write(*,'(T6,"# of markers: ",i10)') sum(nlaunch)
    endif
    !$OMP PARALLEL DO schedule(dynamic,1) private(i,j,k,ic,idcx,ind,vihalo, &
    !$OMP& ri,tracks,ntrack,rates,denn,states,jj,photons,plasma,fi_correction)
    loop_over_cells: do ic = istart, ncell, istep
        call ind2sub(beam_grid%dims,cell_ind(ic),ind)
        i = ind(1) ; j = ind(2) ; k = ind(3)
        !! Loop over the markers
        loop_over_dcx: do idcx=1, nlaunch(i,j,k)
            !! Calculate ri,vhalo and track
            call mc_halo(ind,vihalo,ri)
            call track(ri,vihalo,tracks,ntrack)
            if(ntrack.eq.0) cycle loop_over_dcx

            !! Calculate CX probability
            call get_beam_cx_rate(tracks(1)%ind,ri,vihalo,thermal_ion,neut_types,rates)
            if(sum(rates).le.0.) cycle loop_over_dcx

            !! Solve collisional radiative model along track
            call get_plasma(plasma,pos=tracks(1)%pos)

            !! Weight CX rates by ion source density
            states = rates*plasma%denp
            fi_correction = max((plasma%denp - plasma%denf)/plasma%denp,0.d0)

            loop_along_track: do jj=1,ntrack
                call get_plasma(plasma,pos=tracks(jj)%pos)
                if(.not.plasma%in_plasma) exit loop_along_track
                call colrad(plasma,thermal_ion,vihalo,tracks(jj)%time,states,denn,photons)
                call store_neutrals(tracks(jj)%ind,dcx_type,denn/nlaunch(i,j,k),vihalo,plasma%in_plasma)

                if((photons.gt.0.d0).and.(inputs%calc_dcx.ge.1)) then
                    photons = fi_correction*photons !! Correct for including fast-ions in states
                    call store_bes_photons(tracks(jj)%pos,vihalo,photons/nlaunch(i,j,k),dcx_type)
                endif
            enddo loop_along_track
        enddo loop_over_dcx
    enddo loop_over_cells
    !$OMP END PARALLEL DO

#ifdef _MPI
    !! Combine densities
    call parallel_sum(neut%dcx)
    if(inputs%calc_dcx.ge.1) then
        call parallel_sum(spec%dcx)
    endif
    call parallel_sum(halo_iter_dens(dcx_type))
#endif

end subroutine dcx