beamfus Subroutine

  subroutine beamfus(beamfus0,betbm0,bp,bt,cnbeam,dene,deni,dlamie, &
       ealphadt,enbeam,fdeut,ftrit,ftritbm,sigvdt,ten,tin,vol,zeffai, &
       betanb,dnbeam2,palpnb)

    !! Routine to calculate beam slowing down properties
    !! author: P J Knight, CCFE, Culham Science Centre
    !! beamfus0: input real : multiplier for beam-background fusion calculation
    !! betbm0 : input real :  leading coefficient for neutral beam beta fraction
    !! bp     : input real :  poloidal field (T)
    !! bt     : input real :  toroidal field on axis (T)
    !! cnbeam : input real :  neutral beam current (A)
    !! dene   : input real :  electron density (/m3)
    !! deni   : input real :  fuel ion density (/m3)
    !! dlamie : input real :  ion-electron coulomb logarithm
    !! ealphadt : input real :  alpha particle birth energy (D-T) (keV)
    !! enbeam : input real :  neutral beam energy (keV)
    !! fdeut  : input real :  deuterium fraction of main plasma
    !! ftrit  : input real :  tritium fraction of main plasma
    !! ftritbm: input real :  tritium fraction of neutral beam
    !! sigvdt : input real :  profile averaged <sigma v> for D-T (m3/s)
    !! ten    : input real :  density weighted average electron temperature (keV)
    !! tin    : input real :  density weighted average ion temperature (keV)
    !! vol    : input real :  plasma volume (m3)
    !! zeffai : input real :  mass weighted plasma effective charge
    !! betanb : output real : neutral beam beta component
    !! dnbeam2: output real : hot beam ion density (/m3)
    !! palpnb : output real : alpha power from hot neutral beam ions (MW)
    !! This routine calculates the beam slowing down properties.
    !!     !
    ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    implicit none

    !  Arguments

    real(dp), intent(in) :: beamfus0, betbm0, bp, bt, cnbeam, &
         dene, deni, dlamie, ealphadt, enbeam, fdeut, ftrit, ftritbm, &
         sigvdt, ten, tin, vol, zeffai
    real(dp), intent(out) :: betanb, dnbeam2, palpnb

    !  Local variables

    real(dp) :: denid,denit,ecritd,ecritt,ehotnb,palpdb, &
         palptb,tausl

    ! !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

    !  Velocity slowing down time

    tausl = 1.99D19 * (2.0D0*(1.0D0-ftritbm) + 3.0D0*ftritbm) * &
         ten**1.5D0 / dene / dlamie

    !  Critical energy for electron/ion slowing down of the beam ion
    !  (deuterium and tritium neutral beams, respectively) (keV)

    ecritd = 14.8D0 * ten * 2.0D0 * zeffai**0.6666D0 * &
         (dlamie+4.0D0)/dlamie
    ecritt = ecritd * 1.5D0

    !  Deuterium and tritium ion densities

    denid = deni * fdeut
    denit = deni * ftrit

    !  Perform beam calculations

    call beamcalc(denid,denit,ealphadt,enbeam,ecritd,ecritt,tausl, &
         ftritbm,cnbeam,tin,vol,sigvdt,palpdb,palptb,dnbeam2,ehotnb)

    !  Neutral beam alpha power

    palpnb = beamfus0 * (palpdb + palptb)

    !  Neutral beam beta

    betanb = betbm0 * 4.03D-22 * 0.66666D0 * dnbeam2 * ehotnb / &
         (bt**2 + bp**2)

  end subroutine beamfus