c c ========================================================== subroutine flg3eurhok(q,mx,my,mz,lb,ub,qf,mxo,myo,mzo,lbo,ubo, & lbr,ubr,shaper,meqn,mcomp,time) 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 integer meqn, mx, my, mz, mxo, myo, mzo dimension q(meqn,mx,my,mz), qf(mxo,myo,mzo) c integer lb(3), ub(3), lbo(3), ubo(3), lbr(3), ubr(3), shaper(3), & mresult, stride, imin(3), imax(3), i, j, k, d, getindx c stride = (ub(1) - lb(1))/(mx-1) do 5 d = 1, 3 imin(d) = max(lb(d), lbr(d)) imax(d) = min(ub(d), ubr(d)) if (mod(imin(d)-lb(d),stride) .ne. 0) then imin(d) = imin(d) + stride - mod(imin(d)-lb(d),stride) endif imin(d) = getindx(imin(d), lb(d), stride) if (mod(imax(d)-lb(d),stride) .ne. 0) then imax(d) = imax(d) - mod(imax(d)-lb(d),stride) endif imax(d) = getindx(imax(d), lb(d), stride) 5 continue do 10 k = imin(3), imax(3) do 10 j = imin(2), imax(2) do 10 i = imin(1), imax(1) c Compute rho if (mcomp.eq.1) then qf(i,j,k) = 0.d0 do m = 1, Nsp qf(i,j,k) = qf(i,j,k) + q(m,i,j,k) enddo c # Convert g/cm**3 into kg/m**3 if (ckunits) qf(i,j,k) = qf(i,j,k)*1.d3 c Compute p else if (mcomp.eq.2) then rhoW = 0.d0 do m = 1, Nsp rhoW = rhoW + q(m,i,j,k)/Wk(m) enddo T = q(Nsp+5,i,j,k) qf(i,j,k) = rhoW*RU*T c # Convert dynes/cm**2 into Pa = J/m**2 if (ckunits) qf(i,j,k) = qf(i,j,k)*1.d-1 c Compute gamma else if (mcomp.eq.3) then rhoW = 0.d0 do m = 1, Nsp rhoW = rhoW + q(m,i,j,k)/Wk(m) enddo T = q(Nsp+5,i,j,k) rhoCp = avgtabip(T,q(1,i,j,k),cpk,Nsp) qf(i,j,k) = RU/(rhoCp/rhoW-RU)+1.d0 c Compute Yk else if (mcomp.ge.5.and.mcomp.lt.5+Nsp) then rho = 0.d0 do m = 1, Nsp rho = rho + q(m,i,j,k) enddo qf(i,j,k) = q(mcomp-4,i,j,k)/rho else qf(i,j,k) = 0.d0 endif c 10 continue c return end