subroutine brempde(dens,t,gamma,qbremde)
*
* this subroutine returns the neutrino pair bremsstrahlung rxn rates
* for partially degenerate electrons for the liquid metal regime.
* From Munkata, Kohyama, & Itoh 1987, ApJ, 316, 708.
* Valid in the range:  Gamma 1 to 178. (DO NOT use below Gamma = 1)
*                      Temp  10**8 to 10**11
*                      Dens  1     to 10**12
* Z and AT are only places where composition dependence shows up.
*
      implicit double precision (a-z)

      common/comp/amhyhe,amheca,x(4)
      dimension a(6),b(5)

      t8 = 10.**(t - 8.0)
      d = 10.**(dens)
      at = 1. /(x(2)*4. + x(3)*12. + x(4)*16.)
      z = (x(2)*2. + x(3)*6. + x(4)*8.)
      en = 2.
      s2theta = 0.217
      cv = 0.5 + 2.*s2theta
      ca = 0.5
      cvp = 1. - cv
      cap = 1. - ca
*
* load coefficient arrays. NOTE: a(1) = a_0 etc.
*
      a(1) = 23.5 
      a(2) = 6.83e+04
      a(3) = 7.81e+04
      a(4) = 230.
      a(5) = 6.70e+05
      a(6) = 7.66e+09
      b(1) = 1.47
      b(2) = 0.0329
      b(3) = (7.75e+05*t8**1.5) + (247.*t8**3.85)
      b(4) = 4.07 + (0.0240*t8**1.40)
      b(5) = 4.59e-05*t8**-0.110
*
* crank it out.
*
      xi = d / ((7.05e+06*t8**1.5) + (5.12e+04*t8**3))
      x2 = xi*xi
 
      f1 = 1./(a(1)+(a(2)/t8**2)+(a(3)/t8**5))
      f2 = 1.26*(1.+(1./xi))/(1.+(b(1)/xi)+(b(2)/x2))
      fsum = f1 + f2
      g1 = (1. + 1.e-9*d)*(a(4)+(a(5)/t8**2))+(a(6)/t8**5)
      g2 = (b(3)/d) + b(4) + b(5)*d**0.656
      gsum = (1./g1) + (1./g2)
 
      first = 0.5*(cv*cv + ca*ca + en*cvp*cvp + en*cap*cap)*fsum
      sec = 0.5*(cv*cv - ca*ca + en*cvp*cvp - en*cap*cap)*gsum
      qbremde = 0.5738*(z*z/at)*(t8**6)*(first - sec)

      return
      end 

************************************************************************