pflux_vs_theta_vs_vpa Subroutine

public subroutine pflux_vs_theta_vs_vpa(vflx)

Uses

Diagnose particle flux contribution to toroidal momentum flux in the lab frame in terms of vpar and theta.

Arguments

Type IntentOptional Attributes Name
real, intent(out), dimension (-ntgrid:,:,:) :: vflx

Contents

Source Code


Source Code

  subroutine pflux_vs_theta_vs_vpa (vflx)
#ifdef LOWFLOW   
    use dist_fn_arrays, only: gnew, aj0, g_work
    use gs2_layouts, only: g_lo
    use gs2_layouts, only: it_idx, ik_idx, is_idx
    use theta_grid, only: Rplot !JPL
    use fields_arrays, only: phinew
    use kt_grids, only: aky
    use le_grids, only: integrate_volume, nlambda, negrid
    use le_grids, only: get_flux_vs_theta_vs_vpa
    use species, only:  nspec
    use lowflow, only: mach_lab    
#endif
    use theta_grid, only: ntgrid
    implicit none

    real, dimension (-ntgrid:,:,:), intent (out) :: vflx
#ifdef LOWFLOW  
    integer :: all = 1
    integer :: iglo, isgn, ig, it, ik, is

    real, dimension (:,:,:), allocatable :: g0r
    real, dimension (:,:,:,:,:), allocatable :: gavg
    allocate (g0r(-ntgrid:ntgrid,2,g_lo%llim_proc:g_lo%ulim_alloc))
    allocate (gavg(-ntgrid:ntgrid,nlambda,negrid,2,nspec))

    do iglo = g_lo%llim_proc, g_lo%ulim_proc
       it = it_idx(g_lo,iglo)
       ik = ik_idx(g_lo,iglo)
       is = is_idx(g_lo,iglo)
       do isgn = 1, 2
          do ig = -ntgrid, ntgrid
             g_work(ig,isgn,iglo) = gnew(ig,isgn,iglo)*aj0(ig,iglo) 

             g0r(ig,isgn,iglo) = aimag(g_work(ig,isgn,iglo)*conjg(phinew(ig,it,ik)))*aky(ik)*Rplot(ig)**2*mach_lab(is)

          end do
       end do
    end do
    
    call integrate_volume (g0r, gavg, all)
    call get_flux_vs_theta_vs_vpa (gavg, vflx)

    deallocate (gavg)
    deallocate (g0r)
#else
    vflx=0.
#endif
  end subroutine pflux_vs_theta_vs_vpa