c ===================================================== subroutine combl() c ===================================================== c c Create and initialize application specific common-blocks. c c Copyright (C) 2003-2007 California Institute of Technology c Ralf Deiterding, ralf@amroc.net c implicit double precision (a-h,o-z) c include "ck.i" include "cuser.i" c dimension Xk(leNsp), Yk(leNsp) c call cheminit() c c SI-units c rhol = 0.072d0 ul = 0.d0 vl = 0.d0 wl = 0.d0 pl = 7123.d0 c rhor = 0.18075d0 ur = 0.d0 vr = 0.d0 wr = 0.d0 pr = 35594.d0 c c Convert kg/m**3 into g/cm**3 rhol = rhol * 1.d-3 rhor = rhor * 1.d-3 c Convert Pa = J/m**2 into dynes/cm**2 pl = pl * 1.d1 pr = pr * 1.d1 c Convert m/sec into cm/sec ul = ul * 1.d2 ur = ur * 1.d2 c Convert m/sec into cm/sec vl = vl * 1.d2 vr = vr * 1.d2 u = u * 1.d2 c do k=1,Nsp Xk(k) = 0.d0 enddo Xk(1) = 1.d0 Xk(6) = 2.d0 Xk(9) = 7.d0 c sum = 0.d0 do k = 1, Nsp sum = sum + Xk(k) enddo do k=1,Nsp Xk(k) = Xk(k)/sum enddo c call ckxty (Xk,iwork,rwork,Yk) call ckmmwy (Yk,iwork,rwork,W) c c # data in left state: Tl = (pl*W)/(rhol*RU) rhoul = rhol * ul rhovl = rhol * vl rhowl = rhol * wl el = rhol * (0.5d0*(ul**2+vl**2+wl**2) + & avgtabip(Tl,Yk,hms,Nsp)) - pl c do k=1, Nsp qin(k) = rhol*Yk(k) enddo qin(Nsp+1) = rhoul qin(Nsp+2) = rhovl qin(Nsp+3) = rhowl qin(Nsp+4) = el qin(Nsp+5) = Tl c do k=1,Nsp Xk(k) = 0.d0 enddo Xk(2) = 1.d0 Xk(3) = 2.d0 Xk(9) = 7.d0 c sum = 0.d0 do k = 1, Nsp sum = sum + Xk(k) enddo do k=1,Nsp Xk(k) = Xk(k)/sum enddo c call ckxty (Xk,iwork,rwork,Yk) call ckmmwy (Yk,iwork,rwork,W) c c # data in right state: Tr = (pr*W)/(rhor*RU) rhour = rhor * ur rhovr = rhor * vr rhowr = rhor * wr er = rhor * (0.5d0*(ur**2+vr**2+wr**2) + & avgtabip(Tr,Yk,hms,Nsp)) - pr c write(6,*) 'left state' write(6,*) rhol, ul, vl, wl, pl, Tl write(6,*) 'right state' write(6,*) rhor, ur, vr, wr, pr, Tr c x0 = 0.d0 y0 = 0.d0 z0 = 0.d0 ra = 6.0d0 c do k=1, Nsp qout(k) = rhor*Yk(k) enddo qout(Nsp+1) = rhour qout(Nsp+2) = rhovr qout(Nsp+3) = rhowr qout(Nsp+4) = er qout(Nsp+5) = Tr c pi = 4.d0*atan(1.d0) idisc = 2 c c Geometry for reclecting boundary conditions c angpl = pi/4.d0 rdi = 0.5d0 c c Outer radius - extrapolation c rd = 11.5d0 c return end