iapwsif97水物性計算fortran90程序

1、iapws-if97水物性計算fortran90程序基于iapws-if97的水物性計算fortran90語言程序代碼，程序編寫為動態(tài)鏈接庫(dll)形式，以便于調(diào)用。下面具體講述。如上圖所示，水與水蒸氣可分為5個區(qū)，但本程序只適用于4個區(qū)，即圖中的第1、2、3、4區(qū)，不包括第5區(qū)。壓力與溫度范圍為：程序的輸入變量為壓力p與溫度t，單位分別為pa與k。輸出變量(12個)如下：-比容v內(nèi)能u比焓h比熵s定壓比熱容cp定容比熱容cv聲速wm3/kgj/kgj/kgj/kg*kj/kg*kj/kg*km/s定壓熱膨脹率等溫壓縮率導(dǎo)熱系數(shù)動力粘度表面張力1/k1/paw/m*kpa*sn/m-輸出變量

2、是以一個數(shù)組的形式存儲，例如數(shù)組為prop(12)，數(shù)組與變量的對應(yīng)順序為：vuhscpcvw單位如前所述?，F(xiàn)給出具體的調(diào)用方法：給定壓力p(pa)與溫度t(k)，計算物性，程序自己判斷所屬區(qū)域：call h2o(p,t,prop(12)若想得到給定壓力p(pa)或溫度t(k)下的飽和態(tài)物性，則可如下調(diào)用：a. 給定壓力call psatw(p,t)或給定溫度call tsatw(p,t)b. 飽和液態(tài)call h2o(p,1.001*t,prop(12)或飽和汽態(tài)call h2o(p,0.999*t,prop(12)當(dāng)然還有其他方法，讀者詳細(xì)了解程序內(nèi)容后，請自行發(fā)現(xiàn)。這里給出獲得動態(tài)鏈接庫

3、的方法：1. 新建一個動態(tài)鏈接庫工程2. 加入附錄所示源程序文件3. 編譯，獲得最重要的下面兩個文件：.lib與.dll4. 上面兩個文件就可用于其他程序的物性調(diào)用參考文獻(xiàn)wolfgang wagner, hans-joachim kretzschmar. properties ofwater and steam based on the industrial formulation iapws-if97, second edition.附錄：源程序(19個文件)【1】. const.f90module constimplicit nonereal*8:rg=0.461526e3 !j/kg/

4、kreal*8:tc=647.096 !kreal*8:pc=22.064e6 !pareal*8:rhoc=322.0 !kg/m3end module const【2】. reg1.f90module reg1implicit noneinteger:i(34),j(34)real*8:n(34)data i /&0 ,&0 ,&0 ,&0 ,&0 ,&0 ,&0 ,&0 ,&1 ,&1 ,&1 ,&1 ,&1 ,&1 ,&2 ,&2 ,&2 ,&2 ,&2 ,&3 ,&3 ,&3 ,&4 ,&4 ,&4 ,&5 ,&8 ,&8 ,&21 ,&23 ,&29 ,&30 ,&31 ,&32

5、&/data j /&-2 ,&-1 ,&0 ,&1 ,&2 ,&3 ,&4 ,&5 ,&-9 ,&-7 ,&-1 ,&0 ,&1 ,&3 ,&-3 ,&0 ,&1 ,&3 ,&17 ,&-4 ,&0 ,&6 ,&-5 ,&-2 ,&10 ,&-8 ,&-11 ,&-6 ,&-29 ,&-31 ,&-38 ,&-39 ,&-40 ,&-41 &/data n /&0.14632971213167 ,&-0.84548187169114 ,&-3.756360367204 ,&3.3855169168385 ,&-0.95791963387872 ,&0.15772038513228 ,&-0.

6、016616417199501 ,&8.1214629983568e-04 ,&2.8319080123804e-04 ,&-6.0706301565874e-04 ,&-0.018990068218419 ,&-0.032529748770505 ,&-0.021841717175414 ,&-5.283835796993e-05 ,&-4.7184321073267e-04 ,&-3.0001780793026e-04 ,&4.7661393906987e-05 ,&-4.4141845330846e-06 ,&-7.2694996297594e-16 ,&-3.1679644845054

7、e-05 ,&-2.8270797985312e-06 ,&-8.5205128120103e-10 ,&-2.2425281908e-06 ,&-6.5171222895601e-07 ,&-1.4341729937924e-13 ,&-4.0516996860117e-07 ,&-1.2734301741641e-09 ,&-1.7424871230634e-10 ,&-6.8762131295531e-19 ,&1.4478307828521e-20 ,&2.6335781662795e-23 ,&-1.1947622640071e-23 ,&1.8228094581404e-24 ,&

8、-9.3537087292458e-26 &/end module reg1【3】. reg2.f90module reg2implicit noneinteger:jo(9),i(43),j(43)real*8:no(9),n(43)data jo /&0 ,&1 ,&-5 ,&-4 ,&-3 ,&-2 ,&-1 ,&2 ,&3 &/data no /&-9.6927686500217 ,&10.086655968018 ,&-0.005608791128302 ,&0.071452738081455 ,&-0.40710498223928 ,&1.4240819171444 ,&-4.38

9、3951131945 ,&-0.28408632460772 ,&0.021268463753307 &/data i /&1 ,&1 ,&1 ,&1 ,&1 ,&2 ,&2 ,&2 ,&2 ,&2 ,&3 ,&3 ,&3 ,&3 ,&3 ,&4 ,&4 ,&4 ,&5 ,&6 ,&6 ,&6 ,&7 ,&7 ,&7 ,&8 ,&8 ,&9 ,&10 ,&10 ,&10 ,&16 ,&16 ,&18 ,&20 ,&20 ,&20 ,&21 ,&22 ,&23 ,&24 ,&24 ,&24 &/data j /&0 ,&1 ,&2 ,&3 ,&6 ,&1 ,&2 ,&4 ,&7 ,&36 ,&0

10、 ,&1 ,&3 ,&6 ,&35 ,&1 ,&2 ,&3 ,&7 ,&3 ,&16 ,&35 ,&0 ,&11 ,&25 ,&8 ,&36 ,&13 ,&4 ,&10 ,&14 ,&29 ,&50 ,&57 ,&20 ,&35 ,&48 ,&21 ,&53 ,&39 ,&26 ,&40 ,&58 &/data n /&-1.7731742473213e-03 ,&-0.017834862292358 ,&-0.045996013696365 ,&-0.057581259083432 ,&-0.05032527872793 ,&-3.3032641670203e-05 ,&-1.8948987

11、516315e-04 ,&-3.9392777243355e-03 ,&-0.043797295650573 ,&-2.6674547914087e-05 ,&2.0481737692309e-08 ,&4.3870667284435e-07 ,&-3.227767723857e-05 ,&-1.5033924542148e-03 ,&-0.040668253562649 ,&-7.8847309559367e-10 ,&1.2790717852285e-08 ,&4.8225372718507e-07 ,&2.2922076337661e-06 ,&-1.6714766451061e-11

12、,&-2.1171472321355e-03 ,&-23.895741934104 ,&-5.905956432427e-18 ,&-1.2621808899101e-06 ,&-0.038946842435739 ,&1.1256211360459e-11 ,&-8.2311340897998 ,&1.9809712802088e-08 ,&1.0406965210174e-19 ,&-1.0234747095929e-13 ,&-1.0018179379511e-09 ,&-8.0882908646985e-11 ,&0.10693031879409 ,&-0.33662250574171

13、 ,&8.9185845355421e-25 ,&3.0629316876232e-13 ,&-4.2002467698208e-06 ,&-5.9056029685639e-26 ,&3.7826947613457e-06 ,&-1.2768608934681e-15 ,&7.3087610595061e-29 ,&5.5414715350778e-17 ,&-9.436970724121e-07 &/end module reg2【4】. reg3.f90module reg3implicit noneinteger:i(40),j(40)real*8:n(40)data i /&0 ,&

14、0 ,&0 ,&0 ,&0 ,&0 ,&0 ,&0 ,&1 ,&1 ,&1 ,&1 ,&2 ,&2 ,&2 ,&2 ,&2 ,&2 ,&3 ,&3 ,&3 ,&3 ,&3 ,&4 ,&4 ,&4 ,&4 ,&5 ,&5 ,&5 ,&6 ,&6 ,&6 ,&7 ,&8 ,&9 ,&9 ,&10 ,&10 ,&11 &/data j /&0 ,&0 ,&1 ,&2 ,&7 ,&10 ,&12 ,&23 ,&2 ,&6 ,&15 ,&17 ,&0 ,&2 ,&6 ,&7 ,&22 ,&26 ,&0 ,&2 ,&4 ,&16 ,&26 ,&0 ,&2 ,&4 ,&26 ,&1 ,&3 ,&26 ,&0

15、 ,&2 ,&26 ,&2 ,&26 ,&2 ,&26 ,&0 ,&1 ,&26 &/data n /&1.0658070028513 ,&-15.732845290239 ,&20.944396974307 ,&-7.6867707878716 ,&2.6185947787954 ,&-2.808078114862 ,&1.2053369696517 ,&-8.4566812812502e-03 ,&-1.2654315477714 ,&-1.1524407806681 ,&0.88521043984318 ,&-0.64207765181607 ,&0.38493460186671 ,&-

16、0.85214708824206 ,&4.8972281541877 ,&-3.0502617256965 ,&0.039420536879154 ,&0.12558408424308 ,&-0.2799932969871 ,&1.389979956946 ,&-2.018991502357 ,&-8.2147637173963e-03 ,&-0.47596035734923 ,&0.0439840744735 ,&-0.44476435428739 ,&0.90572070719733 ,&0.70522450087967 ,&0.10770512626332 ,&-0.3291362325

17、8954 ,&-0.50871062041158 ,&-0.022175400873096 ,&0.094260751665092 ,&0.16436278447961 ,&-0.013503372241348 ,&-0.014834345352472 ,&5.7922953628084e-04 ,&3.2308904703711e-03 ,&8.0964802996215e-05 ,&-1.6557679795037e-04 ,&-4.4923899061815e-05 &/end module reg3【5】. reg4.f90module reg4implicit nonereal*8:

18、n(10)data n /&1167.0521452767 ,&-724213.16703206 ,&-17.073846940092 ,&12020.82470247 ,&-3232555.0322333 ,&14.91510861353 ,&-4823.2657361591 ,&405113.40542057 ,&-0.23855557567849 ,&650&/end module reg4【6】. mod_thcon.f90module mod_thcon!implicit nonereal*8:n0(4),n1(5),n2(10)!data n0 /&0.1

19、02811d-1 ,&0.299621d-1 ,&0.156146d-1 ,&-0.422464d-2 &/data n1 /&-0.397070 ,&0.400302 ,&1.06 ,&-0.171587 ,&2.39219 &/data n2 /&0.701309d-1 ,&0.118520d-1 ,&0.642857 ,&0.169937d-2 ,&-1.02 ,&-4.11717 ,&-6.17937 ,&0.0822994 ,&10.0932 ,&0.308976d-2 &/!end module mod_thcon【7】. mod_visc.f90module mod_visc!i

20、mplicit noneinteger:i(21),j(21)real*8:no(4),n(21)!data no /&0.167752e-1 ,&0.220462e-1 ,&0.6366564e-2 ,&-0.241605e-2 &/data i /&0 ,&0 ,&0 ,&0 ,&1 ,&1 ,&1 ,&1 ,&1 ,&2 ,&2 ,&2 ,&2 ,&2 ,&3 ,&3 ,&4 ,&4 ,&5 ,&6 ,&6 &/data j /&0 ,&1 ,&2 ,&3 ,&0 ,&1 ,&2 ,&3 ,&5 ,&0 ,&1 ,&2 ,&3 ,&4 ,&0 ,&1 ,&0 ,&3 ,&4 ,&3 ,&

21、5 &/data n /&0.520094 ,&0.850895e-1 ,&-0.108374e1 ,&-0.289555 ,&0.222531 ,&0.999115 ,&0.188797e1 ,&0.126613e1 ,&0.120573 ,&-0.281378 ,&-0.906851 ,&-0.772479 ,&-0.489837 ,&-0.257040 ,&0.161913 ,&0.257399 ,&-0.325372e-1 ,&0.698452e-1 ,&0.872102e-2 ,&-0.435673e-2 ,&-0.593264e-3 &/!end module mod_visc【8

22、】. b23.f90module b23implicit nonereal*8:n(5)data n /&348.05185628969 ,&-1.1671859879975 ,&1.0192970039326e-03 ,&572.54459862746 ,&13.91883977887 &/end module b23【9】. pb23.f90subroutine pb23(p,t)!dec$ attributes dllexport : pb23!use constuse b23!implicit nonereal*8:p,treal*8:pi,sitareal*8:pstar,tstar

23、!data pstar /1.0e6/ ! padata tstar /1.0/ ! k!if (t .lt. 623.15 .or. t .gt. 863.15) thenstop error: out b23 range !endif!sita=t/tstarpi=n(1)+n(2)*sita+n(3)*sita*sitap=pi*pstar!end subroutine pb23【10】. tb23.f90subroutine tb23(p,t)!dec$ attributes dllexport : tb23!use constuse b23!implicit nonereal*8:p

24、,treal*8:pi,sitareal*8:pstar,tstar!data pstar /1.0e6/ ! padata tstar /1.0/ ! k!if (p .lt. 16.529e6 .or. p .gt. 100.0e6) thenstop error: out b23 range !endif!pi=p/pstarsita=n(4)+sqrt(pi-n(5)/n(3)t=sita*tstar!end subroutine tb23【11】. psatw.f90subroutine psatw(p,t)!dec$ attributes dllexport : psatw!use

25、 constuse reg4!implicit nonereal*8:p,treal*8:pstar,tstar,sitareal*8:a,b,c!data pstar /1.0e6/ ! padata tstar /1.0/ ! k!if (t .lt. 273.15) .or. (t .gt. 647.096) stop error: reg4 (saturate) psatw !sita=t/tstar+n(9)/(t/tstar-n(10)!a=sita*sita+n(1)*sita+n(2)b=n(3)*sita*sita+n(4)*sita+n(5)c=n(6)*sita*sita

26、+n(7)*sita+n(8)!p=(2.0*c/(-b+sqrt(b*b-4.0*a*c)*4.0p=p*pstar!end subroutine psatw【12】. tsatw.f90subroutine tsatw(p,t)!dec$ attributes dllexport : tsatw!use constuse reg4!implicit nonereal*8:p,treal*8:pstar,tstar,sitareal*8:d,e,f,g!data pstar /1.0e6/ ! padata tstar /1.0/ ! k!if (p .lt. 611.657) .or. (

27、p .gt. 22.064e6) stop error: reg4 (saturate) tsatw!sita=(p/pstar)*0.25!e=sita*sita+n(3)*sita+n(6)f=n(1)*sita*sita+n(4)*sita+n(7)g=n(2)*sita*sita+n(5)*sita+n(8)d=2.0*g/(-f-sqrt(f*f-4.0*e*g)!t=(n(10)+d-sqrt(n(10)+d)*2.0-4.0*(n(9)+n(10)*d)/2.0t=t*tstar!end subroutine tsatw【13】. propreg1.f90subroutine p

28、ropreg1(p,t,prop)!dec$ attributes dllexport : propreg1!p: pa;!t: k ;!use constuse reg1!implicit noneinteger:iireal*8:p,treal*8:prop(12)real*8:gamma,pi,tau,gamma_pi,gamma_pipi,gamma_tau,gamma_tautau,gamma_pitaureal*8:pstar,tstar!data pstar /16.529e6/ ! padata tstar /1386.0/ ! k!pi=p/pstartau=tstar/t!

29、gamma=0.0do ii=1,34gamma=gamma+n(ii)*(7.1-pi)*i(ii)*(tau-1.222)*j(ii)enddogamma_pi=0.0do ii=1,34gamma_pi=gamma_pi-n(ii)*i(ii)*(7.1-pi)*(i(ii)-1.0)*(tau-1.222)*j(ii)enddogamma_pipi=0.0do ii=1,34gamma_pipi=gamma_pipi+n(ii)*i(ii)*(i(ii)-1.0)*(7.1-pi)*(i(ii)-2.0)*(tau-1.222)*j(ii)enddogamma_tau=0.0do ii

30、=1,34gamma_tau=gamma_tau+n(ii)*(7.1-pi)*i(ii)*j(ii)*(tau-1.222)*(j(ii)-1.0)enddogamma_tautau=0.0do ii=1,34gamma_tautau=gamma_tautau+n(ii)*(7.1-pi)*i(ii)*j(ii)*(j(ii)-1.0)*(tau-1.222)*(j(ii)-2.0)enddogamma_pitau=0.0do ii=1,34gamma_pitau=gamma_pitau-n(ii)*i(ii)*(7.1-pi)*(i(ii)-1.0)*j(ii)*(tau-1.222)*(

31、j(ii)-1.0)enddo!prop(1)=pi*gamma_pi*rg*t/pprop(2)=(tau*gamma_tau-pi*gamma_pi)*rg*tprop(3)=tau*gamma_tau*rg*tprop(4)=(tau*gamma_tau-gamma)*rgprop(5)=-tau*2.0*gamma_tautau*rgprop(6)=(-tau*2.0*gamma_tautau+(gamma_pi-tau*gamma_pitau)*2.0/gamma_pipi)*rgprop(7)=abs(gamma_pi)*sqrt(rg*t/(gamma_pi-tau*gamma_

32、pitau)*2.0/(tau*tau*gamma_tautau)-gamma_pipi)prop(8)=(1.0-tau*gamma_pitau/gamma_pi)/tprop(9)=-pi*gamma_pipi/(p*gamma_pi)!end subroutine propreg1【14】. propreg2.f90subroutine propreg2(p,t,prop)!dec$ attributes dllexport : propreg2!p: pa;!t: k ;!prop: v u h s cp cv w(sound speed) beta kappa thcon visc

33、sft!use constuse reg2!implicit noneinteger:iireal*8:p,treal*8:prop(12)real*8:pi,tau,gamma0,gamma0_pi,gamma0_pipi,gamma0_tau,gamma0_tautau,gamma0_pitaureal*8:gamma,gamma_pi,gamma_pipi,gamma_tau,gamma_tautau,gamma_pitaureal*8:pstar,tstar!data pstar /1.0e6/ ! padata tstar /540.0/ ! k!pi=p/pstartau=tsta

34、r/t!gamma0=log(pi)do ii=1,9gamma0=gamma0+no(ii)*tau*jo(ii)enddogamma0_pi=1.0/pigamma0_pipi=-1.0/pi*2.0gamma0_tau=0.0do ii=1,9gamma0_tau=gamma0_tau+no(ii)*jo(ii)*tau*(jo(ii)-1.0)enddogamma0_tautau=0.0do ii=1,9gamma0_tautau=gamma0_tautau+no(ii)*jo(ii)*(jo(ii)-1.0)*tau*(jo(ii)-2.0)enddogamma0_pitau=0.0

35、!gamma=0.0do ii=1,43gamma=gamma+n(ii)*pi*i(ii)*(tau-0.5)*j(ii)enddogamma_pi=0.0do ii=1,43gamma_pi=gamma_pi+n(ii)*i(ii)*pi*(i(ii)-1.0)*(tau-0.5)*j(ii)enddogamma_pipi=0.0do ii=1,43gamma_pipi=gamma_pipi+n(ii)*i(ii)*(i(ii)-1.0)*pi*(i(ii)-2.0)*(tau-0.5)*j(ii)enddogamma_tau=0.0do ii=1,43gamma_tau=gamma_ta

36、u+n(ii)*pi*i(ii)*j(ii)*(tau-0.5)*(j(ii)-1.0)enddogamma_tautau=0.0do ii=1,43gamma_tautau=gamma_tautau+n(ii)*pi*i(ii)*j(ii)*(j(ii)-1.0)*(tau-0.5)*(j(ii)-2.0)enddogamma_pitau=0.0do ii=1,43gamma_pitau=gamma_pitau+n(ii)*i(ii)*pi*(i(ii)-1.0)*j(ii)*(tau-0.5)*(j(ii)-1.0)enddo!prop(1)=pi*(gamma0_pi+gamma_pi)

37、*rg*t/pprop(2)=(tau*(gamma0_tau+gamma_tau)-pi*(gamma0_pi+gamma_pi)*rg*tprop(3)=tau*(gamma0_tau+gamma_tau)*rg*tprop(4)=(tau*(gamma0_tau+gamma_tau)-(gamma0+gamma)*rgprop(5)=-tau*2.0*(gamma0_tautau+gamma_tautau)*rgprop(6)=prop(5)-rg*(1.0+pi*gamma_pi-tau*pi*gamma_pitau)*2.0/(1.0-pi*pi*gamma_pipi)prop(7)

38、=(1.0+2.0*pi*gamma_pi+(pi*gamma_pi)*2.0)/(1.0-pi*pi*gamma_pipi)+(1.0+pi*gamma_pi-tau*pi*gamma_pitau)*2.0/(tau*tau*(gamma0_tautau+gamma_tautau)prop(7)=sqrt(prop(7)*rg*t)prop(8)=(1.0+pi*gamma_pi-tau*pi*gamma_pitau)/(t*(1.0+pi*gamma_pi)prop(9)=(1.0-pi*pi*gamma_pipi)/(p*(1.0+pi*gamma_pi)!end subroutine

39、propreg2【15】. propreg3.f90subroutine propreg3(p,t,prop)!dec$ attributes dllexport : propreg3!use constuse reg3!implicit noneinteger:iireal*8:p,t,rho,p1,rho1,rho2real*8:prop(12)real*8:delta,tau,fi,fi_delta,fi_deltadelta,fi_tau,fi_tautau,fi_deltataureal*8:rhostar,tstar!rhostar=rhoctstar=tctau=tstar/t!

40、rho1=1.0e-3rho2=2.98*rhocrho=rhocdodelta=rho/rhostarfi_delta=n(1)/deltado ii=2,40fi_delta=fi_delta+n(ii)*i(ii)*delta*(i(ii)-1.0)*tau*j(ii)enddop1=rho*rg*t*delta*fi_deltaif (abs(p1-p) .le. 1.0) thenexitendifif (p1 .gt. p) thenrho2=rhoelserho1=rhoendifrho=0.5*(rho1+rho2)enddo!delta=rho/rhostarfi=n(1)*

41、log(delta)do ii=2,40fi=fi+n(ii)*delta*i(ii)*tau*j(ii)enddofi_delta=n(1)/deltado ii=2,40fi_delta=fi_delta+n(ii)*i(ii)*delta*(i(ii)-1.0)*tau*j(ii)enddofi_deltadelta=-n(1)/delta*2.0do ii=2,40fi_deltadelta=fi_deltadelta+n(ii)*i(ii)*(i(ii)-1.0)*delta*(i(ii)-2.0)*tau*j(ii)enddofi_tau=0.0do ii=2,40fi_tau=f

42、i_tau+n(ii)*delta*i(ii)*j(ii)*tau*(j(ii)-1.0)enddofi_tautau=0.0do ii=2,40fi_tautau=fi_tautau+n(ii)*delta*i(ii)*j(ii)*(j(ii)-1.0)*tau*(j(ii)-2.0)enddofi_deltatau=0.0do ii=2,40fi_deltatau=fi_deltatau+n(ii)*i(ii)*delta*(i(ii)-1.0)*j(ii)*tau*(j(ii)-1.0)enddo!prop(1)=1.0/rhoprop(2)=tau*fi_tau*rg*tprop(3)

43、=tau*fi_tau-fiprop(4)=tau*fi_tau+delta*fi_deltaprop(6)=-tau*tau*fi_tautauprop(5)=prop(5)+(delta*fi_delta-delta*tau*fi_deltatau)*2.0/(2.0*delta*fi_delta+delta*delta*fi_deltadelta)prop(7)=2.0*delta*fi_delta+delta*delta*fi_deltadelta-(delta*fi_delta-delta*tau*fi_deltatau)*2.0/(tau*tau*fi_tautau)prop(7)

44、=sqrt(prop(7)*rg*t)prop(8)=(fi_delta-tau*fi_deltatau)/(2.0*fi_delta+delta*fi_deltadelta)prop(9)=1.0/(2.0*delta*fi_delta+delta*delta*fi_deltadelta)!end subroutine propreg3【16】. thcon.f90subroutine thcon(rho,t,lamna)!dec$ attributes dllexport : thcon!calculate thermal conductivity: w/(m*k)!use mod_thc

45、on!implicit noneinteger:iireal*8:rho,treal*8:lamna,delta,sita,dsitareal*8:rhostar,tstar,lamnastarreal*8:a0,a1,a2,sumreal*8:aa,bb!lamnastar=1.0 ! w/(m*k)rhostar=317.7tstar=647.26!delta=rho/rhostarsita=t/tstardsita=abs(sita-1.0)+n2(10)!sum=0.0do ii=1,4sum=sum+n0(ii)*sita*(ii-1)enddoa0=sita*0.5*sum!a1=

46、n1(1)+n1(2)*delta+n1(3)*exp(n1(4)*(delta+n1(5)*2.0)!aa=2.0+n2(8)*dsita*(-0.6)if (dsita .ge. 1.0) thenbb=1.0/dsitaelsebb=n2(9)*dsita*(-0.6)endifa2=(n2(1)*sita*(-10.0)+n2(2)*delta*1.8*exp(n2(3)*(1.0-delta*2.8)a2=a2+n2(4)*aa*delta*bb*exp(bb/(bb+1.0)*(1.0-delta*(1.0+bb)a2=a2+n2(5)*exp(n2(6)*sita*1.5+n2(7)*delta*(-5.0)!lamna=a0+a1+a2lamna=lamna*lamnastar!end subroutine thcon【17】. visc.f90subroutine visc(rho,t,eta)!dec$ attributes dllexport : visc!calculate viscosity: pa*s!use constuse mod_visc!implicit noneinteger:iireal*8:rho,treal*8:eta,delta,sitareal*8:etastar,rhostar,tstarreal*8:fi0,fi1,s