c-----------------------------------------------------------------------
c # Update temperature in 1d case
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@amroc.net
c
subroutine TUpdate1(q,mx,lb,ub,lbr,ubr,shaper,meqn)
c
implicit double precision (a-h,o-z)
c
include "ck.i"
c
integer meqn, mx
double precision q(meqn,mx)
integer lb(1), ub(1), lbr(1), ubr(1), shaper(1), getindx
integer i, ic, k, stride
c
if (meqn.le.4) return
c
stride = (ub(1) - lb(1))/(mx-1)
c
do 10 ic=lbr(1), ubr(1), stride
i = getindx(ic, lb(1), stride)
c
rho = 0.d0
rhoW = 0.d0
do k = 1, Nsp
rho = rho + q(k,i)
rhoW = rhoW + q(k,i)/Wk(k)
enddo
rhoek = 0.5d0*q(Nsp+1,i)**2/rho
rhoe = q(Nsp+2,i) - rhoek
call SolveTrhok(q(Nsp+3,i), rhoe, rhoW,
& q(1,i), Nsp, ifail)
if (ifail.ne.0) then
p = rhoW*RU*q(Nsp+3,i)
q(Nsp+2,i) = rhoek +
& avgtabip(T,q(1,i),hms,Nsp) - p
endif
c
10 continue
c
return
end
c