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) include "cuser.i" c parameter( lin = 5, lmechout = 11 ) character *16 cwork dimension cwork(2) c cwork(1)= 'Product' cwork(2)= 'Reactant' c open(unit=lin,status='old',form='formatted',file='init.dat') read (lin, *) gamma, qr read (lin, *) vr, pr read (lin, *) sloc, Nden, NCJ read (lin, *) Wk(1), Wk(2), RU, PA read (lin, *) Nbnd close (lin) c open(unit=lmechout, status='unknown', form='formatted', & file='chem.dat') write (lmechout,400) RU write (lmechout,401) PA write (lmechout,402) (k,cwork(k),k=1,2) write (lmechout,403) (k,Wk(k),k=1,2) close (lmechout) c 400 format('RU ',e16.8) 401 format('PA ',e16.8) 402 format('Sp(',i2.2,') ',a16) 403 format('W(',i2.2,') ',e16.8) c gamma1 = gamma-1.d0 c Yr = 1.d0 ur = 0.d0 if (Nden.eq.1) then rhor = vr else Tr = vr Wr = 1.d0/((1.d0-Yr)/Wk(1) + Yr/Wk(2)) rhor = (pr*Wr)/(RU*Tr) end if V0 = 1.d0/rhor P0 = pr c if (NCJ.eq.1) then Dj = qr Dj = Dj/dsqrt(P0*V0) qn = 0.5d0*((Dj+gamma/Dj)**2-4.d0*gamma)/ & (gamma**2-1.d0) q0 = qn*P0*V0 else q0 = qr qn=q0/(P0*V0) Dj=dsqrt((gamma**2-1.d0)*qn/2.d0) & +dsqrt((gamma**2-1.d0)*qn/2.d0+gamma) end if Vj=gamma*(1.d0+Dj**2)/((gamma+1.d0)*Dj**2) Pj=(1.d0+Dj**2)/(gamma+1.d0) Uj=-Dj*Vj+Dj c write (6,*) 'Normalized CJ state:', & Dj,Pj,1.d0/Vj,Uj,qn c Dj = Dj*dsqrt(P0*V0) Pj = Pj*P0 Vj = Vj*V0 Uj = Uj*dsqrt(P0*V0) c write (6,*) 'Actual CJ state:', & Dj,Pj,1.d0/Vj,Uj,q0 c rhol = 1.d0/Vj ul = Uj pl = Pj Yl = 0.d0 c Yact = 1.d-8 c return end