c cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c cc fitgas routine ccc c cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c subroutine fitgas( a , z , tem , rhom , qgas ) real*8 mue implicit real*8 (a-h,o-z) common /wa/ n c mue = a/z rho = rhom*mue t8 = tem/1.d8 c a0 = 23.5d0 a1 = 6.83d4 a2 = 7.81d8 a3 = 230.d0 a4 = 6.70d5 a5 = 7.66d9 b1 = 1.47d0 b2 = 0.0329d0 b3 = 7.75d5*t8**1.5d0 + 247.d0*t8**3.85d0 b4 = 4.07d0 + 0.024d0*t8**1.4d0 b5 = 4.59d-5*t8**-0.11d0 fem= 56.d0/26.d0 c fem= 1.d0 c xi = fem*rhom/( 7.05d6*t8**1.5d0 + 5.12d4*t8**3.d0 ) fgas = 1.d0/(a0 + a1*t8**-2.d0 + a2*t8**-5.d0) & + ( 1.26d0*(1.d0 + xi**-1.d0) ) & / (1.d0 + b1*xi**-1.d0 + b2*xi**-2.d0) c ggas = 1.d0/(( 1.d0 + 1.d-9*fem*rhom) & *(a3 + a4*t8**-2.d0 + a5*t8**-5.d0)) & +1.d0/( b3/( fem*rhom ) + b4 + b5*fem*rhom**0.656d0) c c cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc c write(*,*) 'fgas & ggas = ',fgas,ggas c cccccccccccccccccccccccccccccccccccccccccccccccccccccccccccc dn = float( n ) snthtw = 0.2325d0 cv = 0.5d0 + 2.d0*snthtw ca = 0.5d0 cvp = 1.d0 - cv cap = 1.d0 - ca c qgas = 0.5738d0*z**2.d0/a*t8**6.d0*rho & *(0.5d0*(( cv **2.d0 + ca **2.d0 ) & + dn*( cvp**2.d0 + cap**2.d0 ))*fgas & - 0.5d0*(( cv **2.d0 - ca **2.d0 ) & + dn*( cvp**2.d0 - cap**2.d0 ))*ggas) c return end