neutron_f Subroutine

  1. fidasim.f90
  2. libfida
  3. neutron_f

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~time_string time_string proc~neutron_f->proc~time_string proc~my_rank my_rank proc~neutron_f->proc~my_rank proc~store_neutrons store_neutrons proc~neutron_f->proc~store_neutrons proc~get_plasma get_plasma proc~neutron_f->proc~get_plasma proc~write_neutrons write_neutrons proc~neutron_f->proc~write_neutrons proc~get_fields get_fields proc~neutron_f->proc~get_fields proc~get_dd_rate get_dd_rate proc~neutron_f->proc~get_dd_rate interface~parallel_sum parallel_sum proc~neutron_f->interface~parallel_sum proc~gyro_correction gyro_correction proc~neutron_f->proc~gyro_correction proc~xyz_to_uvw xyz_to_uvw proc~get_plasma->proc~xyz_to_uvw 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 h5close_f h5close_f proc~write_neutrons->h5close_f h5open_f h5open_f proc~write_neutrons->h5open_f h5fcreate_f h5fcreate_f proc~write_neutrons->h5fcreate_f h5ltmake_dataset_double_f h5ltmake_dataset_double_f proc~write_neutrons->h5ltmake_dataset_double_f interface~h5ltmake_compressed_dataset_double_f h5ltmake_compressed_dataset_double_f proc~write_neutrons->interface~h5ltmake_compressed_dataset_double_f h5ltmake_dataset_int_f h5ltmake_dataset_int_f proc~write_neutrons->h5ltmake_dataset_int_f h5fclose_f h5fclose_f proc~write_neutrons->h5fclose_f h5ltset_attribute_string_f h5ltset_attribute_string_f proc~write_neutrons->h5ltset_attribute_string_f proc~calc_perp_vectors calc_perp_vectors proc~get_fields->proc~calc_perp_vectors proc~get_fields->proc~xyz_to_uvw proc~get_fields->proc~uvw_to_xyz proc~get_fields->proc~in_plasma proc~get_fields->proc~get_position interface~interpol_coeff interpol_coeff proc~get_dd_rate->interface~interpol_coeff proc~parallel_sum_i3 parallel_sum_i3 interface~parallel_sum->proc~parallel_sum_i3 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_d5 parallel_sum_d5 interface~parallel_sum->proc~parallel_sum_d5 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_i2 parallel_sum_i2 interface~parallel_sum->proc~parallel_sum_i2 proc~parallel_sum_d4 parallel_sum_d4 interface~parallel_sum->proc~parallel_sum_d4 proc~parallel_sum_d3 parallel_sum_d3 interface~parallel_sum->proc~parallel_sum_d3 proc~parallel_sum_d2 parallel_sum_d2 interface~parallel_sum->proc~parallel_sum_d2 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 interface~randu randu proc~gyro_correction->interface~randu proc~parallel_sum_i3->proc~parallel_sum_i1 mpi_allreduce mpi_allreduce proc~parallel_sum_i3->mpi_allreduce proc~parallel_sum_d0->mpi_allreduce proc~parallel_sum_d1->proc~parallel_sum_d1 proc~parallel_sum_d1->mpi_allreduce proc~parallel_sum_d5->proc~parallel_sum_d1 proc~parallel_sum_d5->mpi_allreduce proc~parallel_sum_i1->proc~parallel_sum_i1 proc~parallel_sum_i1->mpi_allreduce proc~parallel_sum_i0->mpi_allreduce 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_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_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_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_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_7 h5ltmake_compressed_dataset_double_f_7 interface~h5ltmake_compressed_dataset_double_f->proc~h5ltmake_compressed_dataset_double_f_7 proc~parallel_sum_i2->proc~parallel_sum_i1 proc~parallel_sum_i2->mpi_allreduce proc~cross_product cross_product proc~gyro_step->proc~cross_product proc~in_plasma->proc~xyz_to_uvw proc~in_plasma->interface~interpol_coeff proc~cyl_to_uvw cyl_to_uvw proc~in_plasma->proc~cyl_to_uvw proc~parallel_sum_d4->proc~parallel_sum_d1 proc~parallel_sum_d4->mpi_allreduce proc~parallel_sum_d3->proc~parallel_sum_d1 proc~parallel_sum_d3->mpi_allreduce proc~cyl_to_xyz cyl_to_xyz proc~get_position->proc~cyl_to_xyz proc~parallel_sum_d2->proc~parallel_sum_d1 proc~parallel_sum_d2->mpi_allreduce proc~h5ltmake_compressed_dataset_double_f_4->h5ltmake_dataset_double_f proc~chunk_size chunk_size proc~h5ltmake_compressed_dataset_double_f_4->proc~chunk_size h5sclose_f h5sclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5sclose_f h5pset_chunk_f h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_chunk_f h5pcreate_f h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pcreate_f h5pset_shuffle_f h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_shuffle_f h5pset_deflate_f h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_4->h5pset_deflate_f h5dcreate_f h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_4->h5dcreate_f h5dwrite_f h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_4->h5dwrite_f h5dclose_f h5dclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5dclose_f h5pclose_f h5pclose_f proc~h5ltmake_compressed_dataset_double_f_4->h5pclose_f h5screate_simple_f h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_4->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_2->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_2->proc~chunk_size 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->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_2->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_2->h5dcreate_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->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_2->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_1->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_1->proc~chunk_size 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->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_1->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_1->h5dcreate_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->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_1->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_5->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_5->proc~chunk_size 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->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_5->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_5->h5dcreate_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->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_5->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_3->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_3->proc~chunk_size 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->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_3->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_3->h5dcreate_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->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_3->h5screate_simple_f proc~cyl_to_xyz->proc~uvw_to_xyz proc~cyl_to_xyz->proc~cyl_to_uvw proc~h5ltmake_compressed_dataset_double_f_6->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_6->proc~chunk_size proc~h5ltmake_compressed_dataset_double_f_6->h5sclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_chunk_f proc~h5ltmake_compressed_dataset_double_f_6->h5pcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_6->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_6->h5dcreate_f proc~h5ltmake_compressed_dataset_double_f_6->h5dwrite_f proc~h5ltmake_compressed_dataset_double_f_6->h5dclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_6->h5screate_simple_f proc~h5ltmake_compressed_dataset_double_f_7->h5ltmake_dataset_double_f proc~h5ltmake_compressed_dataset_double_f_7->proc~chunk_size 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->h5pset_shuffle_f proc~h5ltmake_compressed_dataset_double_f_7->h5pset_deflate_f proc~h5ltmake_compressed_dataset_double_f_7->h5dcreate_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->h5pclose_f proc~h5ltmake_compressed_dataset_double_f_7->h5screate_simple_f

Called by

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

Contents

Source Code

neutron_f

Source Code

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

    if(.not.any(thermal_mass.eq.H2_amu)) then
        write(*,'(T2,a)') 'NEUTRON_F: Thermal Deuterium is not present in plasma'
        return
    endif
    if(any(thermal_mass.eq.H3_amu)) then
        write(*,'(T2,a)') 'NEUTRON_F: D-T neutron production is not implemented'
    endif
    if(beam_mass.ne.H2_amu) then
        write(*,'(T2,a)') 'NEUTRON_F: Fast-ion species is not Deuterium'
        return
    endif

    if(inputs%calc_neutron.ge.2) then
        allocate(neutron%weight(fbm%nenergy,fbm%npitch,fbm%nr,fbm%nz,fbm%nphi))
        neutron%weight = 0.d0
        allocate(neutron%emis(fbm%nr,fbm%nz,fbm%nphi))
        neutron%emis = 0.d0
    endif

    ngamma = 20
    rate = 0
    !$OMP PARALLEL DO schedule(guided) private(fields,vi,ri,pitch,eb,&
    !$OMP& ir,iphi,iz,ie,ip,igamma,plasma,factor,uvw,uvw_vi,vnet_square,rate,erel,s,c)
    z_loop: do iz = istart, fbm%nz, istep
        r_loop: do ir=1, fbm%nr
            phi_loop: do iphi = 1, fbm%nphi
                !! Calculate position
                if(fbm%nphi.eq.1) then
                    s = 0.d0
                    c = 1.d0
                else
                    s = sin(fbm%phi(iphi))
                    c = cos(fbm%phi(iphi))
                endif

                uvw(1) = fbm%r(ir)*c
                uvw(2) = fbm%r(ir)*s
                uvw(3) = fbm%z(iz)

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

                factor = fbm%r(ir)*fbm%dE*fbm%dp*fbm%dr*fbm%dz*fbm%dphi/ngamma
                !! Loop over energy/pitch/gamma
                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 igamma=1, ngamma
                            call gyro_correction(fields,eb,pitch,fbm%A,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*fbm%A*vnet_square ![kev]

                            !! Get neutron production rate
                            call get_dd_rate(plasma, erel, rate, branch=2)
                            if(inputs%calc_neutron.ge.2) then
                                neutron%weight(ie,ip,ir,iz,iphi) = neutron%weight(ie,ip,ir,iz,iphi) &
                                                                 + rate * factor
                                !$OMP CRITICAL(neutron_emis)
                                neutron%emis(ir,iz,iphi) = neutron%emis(ir,iz,iphi) &
                                                                 + rate * fbm%f(ie,ip,ir,iz,iphi) &
                                                                 * factor
                                !$OMP END CRITICAL(neutron_emis)
                            endif

                            rate = rate*fbm%f(ie,ip,ir,iz,iphi)*factor
                            !! Store neutrons
                            call store_neutrons(rate)
                        enddo gyro_loop
                    enddo energy_loop
                enddo pitch_loop
            enddo phi_loop
        enddo r_loop
    enddo z_loop
    !$OMP END PARALLEL DO

#ifdef _MPI
    call parallel_sum(neutron%rate)
    if(inputs%calc_neutron.ge.2) then
        call parallel_sum(neutron%weight)
        call parallel_sum(neutron%emis)
    endif
#endif

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

#ifdef _MPI
    if(my_rank().eq.0) call write_neutrons()
#else
    call write_neutrons()
#endif


end subroutine neutron_f