get_flux_dist Subroutine

private subroutine get_flux_dist(g_in, fld, flx_dist)

Identical to get_flux except don't integrate over poloidal angle ! JRB

Arguments

Type IntentOptional Attributes Name
complex, intent(in), dimension (-ntgrid:,:,g_lo%llim_proc:) :: g_in

Input weighted distribution

complex, intent(in), dimension (-ntgrid:,:,:) :: fld
real, intent(inout), dimension (-ntgrid:,:,:,:) :: flx_dist

Contents

Source Code


Source Code

  subroutine get_flux_dist (g_in, fld, flx_dist)
    use theta_grid, only: ntgrid, grho, theta, bmag, gradpar
    use kt_grids, only: ntheta0, aky, naky
    use le_grids, only: integrate_moment
    use species, only: nspec
    use integration, only: trapezoidal_integration
    use gs2_time, only: woutunits
    use gs2_layouts, only: g_lo
    implicit none
    logical, parameter :: full_arr=moment_to_allprocs
    !> Input weighted distribution    
    complex, dimension (-ntgrid:,:,g_lo%llim_proc:), intent (in) :: g_in
    complex, dimension (-ntgrid:,:,:), intent (in) :: fld
    real, dimension (-ntgrid:,:,:,:), intent (in out) :: flx_dist
    complex, dimension (:,:,:,:), allocatable :: total
    real :: wgt
    integer :: ik, it, is

    allocate (total(-ntgrid:ntgrid,ntheta0,naky,nspec))
    ! EGH added final parameter 'all'
    ! to the call below for new parallel output
    ! This is temporary until distributed fields
    ! are implemented 1/2014
    call integrate_moment (g_in, total, moment_to_allprocs, full_arr)

    !Note that jacob = 1/(bmag*gradpar*drhodpsi). This is similar to the factor
    !we have below without the drhodpsi.
    wgt = 1.0/trapezoidal_integration(theta, grho / (bmag * gradpar))

    ! EGH for new parallel I/O everyone
    ! calculates fluxes
    ! Note: This doesn't include the field line jacobian (jacob) and as such
    ! it's not clear that the above wgt is correct (as it uses 1/(bmag*gradpar)
    ! which is jacob * drhodpsi, rather than jacob as one might expect).
    !if (proc0) then
       do is = 1, nspec
          do ik = 1, naky
             do it = 1, ntheta0
                flx_dist(:,it,ik,is) = aimag(total(:,it,ik,is)*conjg(fld(:,it,ik))) &
                     *aky(ik) * wgt / woutunits(ik)
             end do
          end do
       end do
       flx_dist = flx_dist*0.5
    ! end if

    deallocate (total)

  end subroutine get_flux_dist