neutron_f Subroutine

public subroutine neutron_f()

Calculate neutron emission rate using a fast-ion distribution function F(E,p,r,z)

Arguments

None

Calls

proc~~neutron_f~~CallsGraph proc~neutron_f neutron_f proc~get_fields get_fields proc~neutron_f->proc~get_fields proc~store_neutrons store_neutrons proc~neutron_f->proc~store_neutrons proc~get_plasma get_plasma proc~neutron_f->proc~get_plasma proc~get_neutron_rate get_neutron_rate proc~neutron_f->proc~get_neutron_rate proc~uvw_to_xyz uvw_to_xyz proc~neutron_f->proc~uvw_to_xyz proc~gyro_correction gyro_correction proc~neutron_f->proc~gyro_correction proc~write_neutrons write_neutrons proc~neutron_f->proc~write_neutrons 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~get_plasma->proc~in_plasma interface~interpol_coeff interpol_coeff proc~get_neutron_rate->interface~interpol_coeff proc~randu randu 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 h5fcreate_f h5fcreate_f proc~write_neutrons->h5fcreate_f h5ltmake_dataset_double_f h5ltmake_dataset_double_f proc~write_neutrons->h5ltmake_dataset_double_f h5open_f h5open_f proc~write_neutrons->h5open_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_neutrons->h5ltset_attribute_string_f h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_neutrons->h5ltmake_dataset_int_f h5close_f h5close_f proc~write_neutrons->h5close_f interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_neutrons->interface~h5ltmake_compressed_dataset_double_f h5fclose_f h5fclose_f proc~write_neutrons->h5fclose_f 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~interpol1d_coeff_arr interpol1D_coeff_arr interface~interpol_coeff->proc~interpol1d_coeff_arr proc~rng_uniform rng_uniform proc~randu->proc~rng_uniform omp_get_thread_num omp_get_thread_num proc~randu->omp_get_thread_num 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~h5ltmake_compressed_dataset_double_f_6 h5ltmake_compressed_dataset_double_f_6 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_6 proc~h5ltmake_compressed_dataset_double_f_5 h5ltmake_compressed_dataset_double_f_5 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_5 proc~h5ltmake_compressed_dataset_double_f_4 h5ltmake_compressed_dataset_double_f_4 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_4 proc~h5ltmake_compressed_dataset_double_f_3 h5ltmake_compressed_dataset_double_f_3 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_3 proc~h5ltmake_compressed_dataset_double_f_2 h5ltmake_compressed_dataset_double_f_2 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_2 proc~h5ltmake_compressed_dataset_double_f_1 h5ltmake_compressed_dataset_double_f_1 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_1 proc~h5ltmake_compressed_dataset_double_f_7 h5ltmake_compressed_dataset_double_f_7 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_7 proc~interpol2d_coeff_arr->proc~interpol2d_coeff proc~h5ltmake_compressed_dataset_double_f_6->h5ltmake_dataset_double_f h5dwrite_f h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_6->h5dwrite_f h5dclose_f h5dclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5dclose_f h5sclose_f h5sclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5sclose_f h5pset_chunk_f h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_chunk_f h5pcreate_f h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5pcreate_f proc~chunk_size chunk_size proc~h5ltmake_compressed_dataset_double_f_6->proc~chunk_size h5pset_deflate_f h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_deflate_f h5pset_shuffle_f h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_shuffle_f h5screate_simple_f h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_6->h5screate_simple_f h5pclose_f h5pclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5pclose_f h5dcreate_f h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_5->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_5->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_5->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_5->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_5->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_5->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_5->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_5->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_4->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_4->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_4->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_4->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_4->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_4->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_4->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_3->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_3->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_3->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_3->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_3->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_3->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_3->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_3->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_2->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_2->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_2->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_2->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_2->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_2->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_2->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_2->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_1->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_1->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_1->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_1->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_1->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_1->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_1->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_1->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_7->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_7->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_7->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_7->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_7->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_7->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_7->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_7->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5dcreate_f proc~interpol1d_coeff_arr->proc~interpol1d_coeff

Called by

proc~~neutron_f~~CalledByGraph proc~neutron_f neutron_f program~fidasim fidasim program~fidasim->proc~neutron_f

Contents

Source Code


Source Code

subroutine neutron_f
    !+ Calculate neutron emission rate using a fast-ion distribution function F(E,p,r,z)
    integer :: ir, iz, ie, ip, iphi, nphi
    type(LocalProfiles) :: plasma
    type(LocalEMFields) :: fields
    real(Float64) :: eb,pitch,r,z
    real(Float64) :: erel, rate
    real(Float64), dimension(3) :: rg, ri
    real(Float64), dimension(3) :: vi
    real(Float64), dimension(3) :: uvw, uvw_vi
    real(Float64)  :: vnet_square, factor
    real(Float64)  :: maxcnt, inv_maxcnt, cnt

    allocate(neutron%weight(fbm%nenergy,fbm%npitch,fbm%nr,fbm%nz))
    neutron%weight = 0.d0

    nphi = 20
    maxcnt=fbm%nr*fbm%nz
    inv_maxcnt = 100.d0/maxcnt
    cnt=0.0
    !$OMP PARALLEL DO schedule(guided) private(fields,vi,ri,rg,pitch,eb,&
    !$OMP& ir,iz,ie,ip,iphi,plasma,factor,uvw,uvw_vi,vnet_square,rate,erel)
    z_loop: do iz = 1, fbm%nz
        r_loop: do ir=1, fbm%nr
            cnt = cnt+1
            !! Calculate position
            uvw(1) = fbm%r(ir)
            uvw(2) = 0.d0
            uvw(3) = fbm%z(iz)
            call uvw_to_xyz(uvw, rg)

            !! Get fields
            call get_fields(fields,pos=rg)
            if(.not.fields%in_plasma) cycle r_loop

            factor = 2*pi*fbm%r(ir)*fbm%dE*fbm%dp*fbm%dr*fbm%dz/nphi
            !! Loop over energy/pitch/phi
            pitch_loop: do ip = 1, fbm%npitch
                pitch = fbm%pitch(ip)
                energy_loop: do ie =1, fbm%nenergy
                    eb = fbm%energy(ie)
                    gyro_loop: do iphi=1, nphi
                        call gyro_correction(fields,eb,pitch,ri,vi)

                        !! Get plasma parameters at particle position
                        call get_plasma(plasma,pos=ri)
                        if(.not.plasma%in_plasma) cycle gyro_loop

                        !! Calculate effective beam energy
                        vnet_square=dot_product(vi-plasma%vrot,vi-plasma%vrot)  ![cm/s]
                        erel = v2_to_E_per_amu*inputs%ab*vnet_square ![kev]

                        !! Get neutron production rate
                        call get_neutron_rate(plasma, erel, rate)
                        neutron%weight(ie,ip,ir,iz) = neutron%weight(ie,ip,ir,iz) &
                                                    + rate*factor
                        rate = rate*fbm%f(ie,ip,ir,iz)*factor

                        !! Store neutrons
                        call store_neutrons(rate)
                    enddo gyro_loop
                enddo energy_loop
            enddo pitch_loop
            if (inputs%verbose.ge.2)then
              WRITE(*,'(f7.2,"% completed",a,$)') cnt*inv_maxcnt,char(13)
            endif
        enddo r_loop
    enddo z_loop
    !$OMP END PARALLEL DO

    if(inputs%verbose.ge.1) then
        write(*,'(T4,A,ES14.5," [neutrons/s]")') 'Rate:   ',sum(neutron%rate)
        write(*,'(30X,a)') ''
    endif

    call write_neutrons()

end subroutine neutron_f