nbi_spec Subroutine

public subroutine nbi_spec()

Calculates approximate neutral beam emission (full, half, third) from user supplied neutrals file

Arguments

None

Calls

proc~~nbi_spec~~CallsGraph proc~nbi_spec nbi_spec proc~ind2sub ind2sub proc~nbi_spec->proc~ind2sub proc~store_photons store_photons proc~nbi_spec->proc~store_photons proc~mc_nbi_cell mc_nbi_cell proc~nbi_spec->proc~mc_nbi_cell interface~randu randu proc~nbi_spec->interface~randu proc~in_plasma in_plasma proc~nbi_spec->proc~in_plasma interface~parallel_sum parallel_sum proc~nbi_spec->interface~parallel_sum proc~get_passive_grid_indices get_passive_grid_indices proc~store_photons->proc~get_passive_grid_indices proc~spectrum spectrum proc~store_photons->proc~spectrum proc~get_fields get_fields proc~store_photons->proc~get_fields proc~get_indices get_indices proc~store_photons->proc~get_indices proc~uvw_to_xyz uvw_to_xyz proc~store_photons->proc~uvw_to_xyz proc~mc_nbi_cell->interface~randu proc~cyl_to_uvw cyl_to_uvw proc~in_plasma->proc~cyl_to_uvw proc~xyz_to_uvw xyz_to_uvw proc~in_plasma->proc~xyz_to_uvw interface~interpol_coeff interpol_coeff proc~in_plasma->interface~interpol_coeff 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~parallel_sum_i3 parallel_sum_i3 interface~parallel_sum->proc~parallel_sum_i3 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~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_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_d0 parallel_sum_d0 interface~parallel_sum->proc~parallel_sum_d0 proc~xyz_to_cyl xyz_to_cyl proc~get_passive_grid_indices->proc~xyz_to_cyl proc~parallel_sum_d1->proc~parallel_sum_d1 mpi_allreduce mpi_allreduce proc~parallel_sum_d1->mpi_allreduce proc~parallel_sum_i2->proc~parallel_sum_i1 proc~parallel_sum_i2->mpi_allreduce proc~parallel_sum_i3->proc~parallel_sum_i1 proc~parallel_sum_i3->mpi_allreduce proc~parallel_sum_d3->proc~parallel_sum_d1 proc~parallel_sum_d3->mpi_allreduce proc~get_fields->proc~in_plasma proc~get_fields->proc~xyz_to_uvw proc~get_fields->proc~uvw_to_xyz proc~calc_perp_vectors calc_perp_vectors proc~get_fields->proc~calc_perp_vectors proc~get_position get_position proc~get_fields->proc~get_position proc~parallel_sum_d2->proc~parallel_sum_d1 proc~parallel_sum_d2->mpi_allreduce proc~parallel_sum_d4->proc~parallel_sum_d1 proc~parallel_sum_d4->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~parallel_sum_d0->mpi_allreduce proc~xyz_to_cyl->proc~xyz_to_uvw proc~uvw_to_cyl uvw_to_cyl proc~xyz_to_cyl->proc~uvw_to_cyl proc~cyl_to_xyz cyl_to_xyz proc~get_position->proc~cyl_to_xyz proc~cyl_to_xyz->proc~cyl_to_uvw proc~cyl_to_xyz->proc~uvw_to_xyz

Called by

proc~~nbi_spec~~CalledByGraph proc~nbi_spec nbi_spec program~fidasim fidasim program~fidasim->proc~nbi_spec

Contents

Source Code


Source Code

subroutine nbi_spec
    !+ Calculates approximate neutral beam emission (full, half, third)
    !+ from user supplied neutrals file
    integer :: ic, i, j, k, it
    real(Float64), dimension(3) :: ri, vnbi, random3, rc
    integer,dimension(3) :: ind
    !! Determination of the CX probability
    real(Float64) :: nbif_photons, nbih_photons, nbit_photons
    real(Float64) :: f_wght, h_wght, t_wght
    real(Float64) :: f_tot, h_tot, t_tot
    real(Float64), dimension(n_stark,inputs%nlambda,spec_chords%nchan) :: full, half, third
    real(Float64), dimension(:,:,:,:), allocatable :: fullstokes, halfstokes, thirdstokes
    logical :: inp
    integer :: n = 10000

    !$OMP PARALLEL DO schedule(dynamic,1) private(i,j,k,ic,ind, &
    !$OMP& nbif_photons, nbih_photons, nbit_photons, rc, ri,inp, vnbi,&
    !$OMP& random3,f_tot,h_tot,t_tot,f_wght,h_wght,t_wght,&
    !$OMP& full,half,third, &
    !$OMP& fullstokes, halfstokes, thirdstokes)
    loop_over_cells: do ic = istart, spec_chords%ncell, istep
        call ind2sub(beam_grid%dims,spec_chords%cell(ic),ind)
        i = ind(1) ; j = ind(2) ; k = ind(3)

        nbif_photons = neut%full%dens(3,i,j,k)*tables%einstein(2,3)
        nbih_photons = neut%half%dens(3,i,j,k)*tables%einstein(2,3)
        nbit_photons = neut%third%dens(3,i,j,k)*tables%einstein(2,3)
        if((nbif_photons + nbih_photons + nbit_photons).le.0.0) then
            cycle loop_over_cells
        endif

        rc = [beam_grid%xc(i), beam_grid%yc(j), beam_grid%zc(k)]

        !Find a point in cell that is also in the plasma
        ri = rc
        call in_plasma(ri, inp)
        do while (.not.inp)
            call randu(random3)
            ri = rc + beam_grid%dr*(random3 - 0.5)
            call in_plasma(ri,inp)
        enddo
        allocate(fullstokes(n_stark,4,inputs%nlambda,spec_chords%nchan))
        allocate(halfstokes(n_stark,4,inputs%nlambda,spec_chords%nchan))
        allocate(thirdstokes(n_stark,4,inputs%nlambda,spec_chords%nchan))
        f_tot = 0.0 ; h_tot = 0.0 ; t_tot = 0.0
        full  = 0.0 ; half  = 0.0 ; third = 0.0
        fullstokes  = 0.0 ; halfstokes  = 0.0 ; thirdstokes = 0.0
        do it=1, n
            !! Full Spectra
            call mc_nbi_cell(ind, nbif_type, vnbi, f_wght)
            f_tot = f_tot + f_wght
            call store_photons(ri, vnbi, beam_lambda0, f_wght*nbif_photons, full, fullstokes)

            !! Half Spectra
            call mc_nbi_cell(ind, nbih_type, vnbi, h_wght)
            h_tot = h_tot + h_wght
            call store_photons(ri, vnbi, beam_lambda0, h_wght*nbih_photons, half, halfstokes)

            !! Third Spectra
            call mc_nbi_cell(ind, nbit_type, vnbi, t_wght)
            t_tot = t_tot + t_wght
            call store_photons(ri, vnbi, beam_lambda0, t_wght*nbit_photons, third, thirdstokes)
        enddo
        !$OMP CRITICAL(nbi_spec_1)
        spec%full = spec%full + full/f_tot
        spec%half = spec%half + half/h_tot
        spec%third = spec%third + third/t_tot
        spec%fullstokes = spec%fullstokes + fullstokes/f_tot
        spec%halfstokes = spec%halfstokes + halfstokes/h_tot
        spec%thirdstokes = spec%thirdstokes + thirdstokes/t_tot
        !$OMP END CRITICAL(nbi_spec_1)
        deallocate(fullstokes)
        deallocate(halfstokes)
        deallocate(thirdstokes)
    enddo loop_over_cells
    !$OMP END PARALLEL DO

#ifdef _MPI
    !! Combine Spectra
    call parallel_sum(spec%full)
    call parallel_sum(spec%half)
    call parallel_sum(spec%third)
#endif
end subroutine nbi_spec