gyro_range Subroutine

  1. fidasim.f90
  2. libfida
  3. gyro_range

public subroutine gyro_range(b, gs, gyrange, nrange)

Calculates the range(s) of gyro-angles that would land within a bounded plane

Arguments

TypeIntentOptionalAttributesName
type(BoundedPlane), intent(in) :: b

Bounded Plane
type(GyroSurface), intent(in) :: gs

Gyro-surface
real(kind=Float64), intent(out), dimension(2,4):: gyrange

(theta, dtheta) values
integer, intent(out) :: nrange

Number of ranges. 1 <= nrange <= 4

Calls

proc~~gyro_range~~CallsGraph proc~gyro_range gyro_range proc~in_boundary in_boundary proc~gyro_range->proc~in_boundary proc~line_plane_intersect line_plane_intersect proc~gyro_range->proc~line_plane_intersect proc~line_gyro_surface_intersect line_gyro_surface_intersect proc~gyro_range->proc~line_gyro_surface_intersect proc~boundary_edge boundary_edge proc~gyro_range->proc~boundary_edge proc~gyro_surface_coordinates gyro_surface_coordinates proc~gyro_range->proc~gyro_surface_coordinates proc~in_gyro_surface in_gyro_surface proc~gyro_range->proc~in_gyro_surface

Called by

proc~~gyro_range~~CalledByGraph proc~gyro_range gyro_range proc~npa_gyro_range npa_gyro_range proc~npa_gyro_range->proc~gyro_range proc~npa_f npa_f proc~npa_f->proc~npa_gyro_range proc~pnpa_f pnpa_f proc~pnpa_f->proc~npa_gyro_range proc~npa_mc npa_mc proc~npa_mc->proc~npa_gyro_range proc~pnpa_mc pnpa_mc proc~pnpa_mc->proc~npa_gyro_range program~fidasim fidasim program~fidasim->proc~npa_f program~fidasim->proc~pnpa_f program~fidasim->proc~npa_mc program~fidasim->proc~pnpa_mc

Contents

Source Code

gyro_range

Source Code

subroutine gyro_range(b, gs, gyrange, nrange)
    !+ Calculates the range(s) of gyro-angles that would land within a bounded plane
    type(BoundedPlane), intent(in)             :: b
        !+ Bounded Plane
    type(GyroSurface), intent(in)              :: gs
        !+ Gyro-surface
    real(Float64), dimension(2,4), intent(out) :: gyrange
        !+ (theta, dtheta) values
    integer, intent(out) :: nrange
        !+ Number of ranges. `1 <= nrange <= 4`

    integer :: nb, i, j, ninter
    logical :: in_gs
    logical, dimension(8) :: cross = .False.
    real(Float64) :: t_p, th1, th2, dth
    real(Float64), dimension(2) :: u_cur, t_i
    real(Float64), dimension(3) :: rc, p_pre, p_cur, v0, ri
    real(Float64), dimension(2,8) :: u
    real(Float64), dimension(3,50) :: bedge

    nrange = 0
    gyrange = 0.d0

    call line_plane_intersect(gs%center, gs%basis(:,3), b%origin, b%basis(:,3), rc, t_p)
    if(t_p.eq.0.0) return

    call boundary_edge(b, bedge, nb)
    p_pre = bedge(:,1)
    in_gs = in_gyro_surface(gs, p_pre)

    ninter = 0
    u = 0.d0
    boundary_loop: do i=1,nb
        p_cur = bedge(:,modulo(i,nb)+1)
        v0 = p_cur - p_pre
        call line_gyro_surface_intersect(p_pre, v0, gs, t_i)
        do j=1,2
            if((t_i(j).gt.0.0).and.(t_i(j).lt.1.0)) then
                ri = p_pre + t_i(j)*v0
                call gyro_surface_coordinates(gs, ri, u_cur)
                if(u_cur(2).gt.0.0) then
                    in_gs = .not.in_gs
                    ninter = ninter + 1
                    cross(ninter) = in_gs
                    u(:,ninter) = u_cur
                endif
            endif
        enddo
        p_pre = p_cur
    enddo boundary_loop

    if(ninter.eq.0) then
        if(in_boundary(b, rc)) then
            nrange = 1
            gyrange(:,1) = [0.d0,2*pi]
        endif
        return
    endif

    do i=1, ninter
        if(cross(i)) then
            th1 = u(1,i)
            j = modulo(i,ninter) + 1
            th2 = u(1,j)
            dth = th2-th1
            nrange = nrange + 1
            if(dth.gt.0.0) then
                gyrange(:,nrange) = [th1, dth]
            else
                gyrange(:,nrange) = [th2, -dth]
            endif
        endif
        !! OpenMP with multiple threads is duplicating gyro-ranges for some markers
        !! causing double counting and I don't know why.
        !! It should be very unlikely for multiple gyro-ranges to occur so for
        !! now I'm including this cludge to force only one gyro-range when using
        !! OpenMP.
#ifdef _OMP
        if(nrange.eq.1) exit
#endif
    enddo

end subroutine gyro_range