c
c =========================================================
subroutine flx1(maxmx,meqn,mbc,mx,q,maux,aux,f)
c =========================================================
c
c # Copyright (C) 2002 Ralf Deiterding
c # Brandenburgische Universitaet Cottbus
c
c # Copyright (C) 2003-2007 California Institute of Technology
c # Ralf Deiterding, ralf@cacr.caltech.edu
c
implicit double precision (a-h,o-z)
include "ck.i"
c
dimension q(1-mbc:maxmx+mbc, meqn)
dimension aux(1-mbc:maxmx+mbc, maux)
dimension f(1-mbc:maxmx+mbc, meqn)
dimension rk(LeNsp)
c
mu = Nsp+1
mE = Nsp+2
mT = Nsp+3
c
do i=1-mbc,mx+mbc
rho = 0.d0
rhoW = 0.d0
do k = 1, Nsp
rk(k) = q(i,k)
rho = rho + q(i,k)
rhoW = rhoW + q(i,k)/Wk(k)
enddo
u = q(i,mu)/rho
rhoe = q(i,mE)-0.5d0*rho*u**2
call SolveTrhok(q(i,mT),rhoe,rhoW,rk,Nsp,ifail)
if (ifail.ne.0) write(6,100) i,q(i,mT)
p = rhoW*RU*q(i,mT)
do k = 1, Nsp
f(i,k) = u*q(i,k)
enddo
f(i,mu) = u*q(i,mu)+p
f(i,mE) = u*(q(i,mE)+p)
f(i,mT) = 0.d0
enddo
c
100 format('flx1eurhok q(',i2,'): T out of range using: ',f16.8)
return
end
c