Барилович В.А. Основы термогазодинамики двухфазных потоков и их численное решение
Подождите немного. Документ загружается.
310
PROGRAM SCCA
RO(A)=999.9-0.1475*(A-273.)-2.695E-3*(A-273.)**2
OPEN(2,FILE='SCCNDA.DAT')
OPEN(8,FILE='REZSCA.DAT')
READ(2,*) GSM,Z0,D,PI,P1,DB1
F=PI*D**2/4.
TV1=5.632*P1**0.25+273.
ROV1=1.213E-5*((TV1-273.)/5.677)**3.75
ROL1=RO(TV1)
ROSM1=ROL1/2.
WSM1=GSM/ROSM1/F
X1=ROV1/(ROL1-ROV1)
SB=6.*GSM*X1/PI/(DB1**3)/ROV1
ROSG=ROL1*(1.-PI*(2.**0.5)/6.)
P=P1
J=0
WRITE(8,*)' Z P SB ROSM WSM DB
& FV ROSG TV XV FL ROV ALFA BETA'
DO 2 I=1,201
Z=(I-1)*1.E-4
DB=DB1*(1.-Z/Z0)
TV=5.632*P**0.25+273.
ROV=1.213E-5*((TV1-273.)/5.677)**3.75
XV=ROV*PI*DB**3*SB/GSM/6.
ROSM=ROL1*ROV/((1.-XV)*ROV+XV*ROL1)
ALFA=XV*ROSM/ROV
WSM=WSM1*ROSM1/ROSM
P=P1+ROSM1*WSM1**2-ROSM*WSM**2
FV=PI*DB**3*SB/WSM/6.
FL=F-FV
BETA=FV/F
WRITE(8,200)
Z,P,SB,ROSM,WSM,DB,FV,ROSG,TV,XV,FL,ROV,ALFA,BETA
J=J+1
IF(J.EQ.10) THEN
WRITE(8,*)' Z P SB ROSM WSM DB
& FV ROSG TV XV FL ROV ALFA BETA'
J=0
END IF
IF(Z.EQ.19.9E-3) GOTO 3
2 CONTINUE
200 FORMAT(F6.4,F10.0,E10.2,2F8.3,2E10.2,F8.3,F6.0,2E10.2,F8.3,
311
&2F6.3)
3 STOP
END
Файл исходных данных
0.415 20.E-3 4.5E-3 3.1416 36.E3 8.E-5
GSM Z0 D PI P1 DB1
PROGRAM SCC2
RO(A)=999.9-0.1475*(A-273.)-2.695E-3*(A-273.)**2
YIL(A)=(12.475+3.935*(A-273.)+1.3125E-3*(A-273.)**2)*1.E3
SIG(A)=(75.54-0.168*(A-273.))*1.E-3
SQV(A)=(2447.7-1.283*(A-273.)-6.26E-3*(A-273.)**2)*1.E3
CPL(A)=4275.5-1.2209*(A-273.)+0.59823E-2*(A-273.)**2
VI(A)=(2461.1+2.636*(A-273.)-4.88E-3*(A-273.)**2)*1.E3
OPEN(2,FILE='SCCOND2.DAT')
OPEN(8,FILE='REZSCND2.DAT')
READ(2,*) GSM,Z0,D,PI,P1,DB1,TL1,AMV
F=PI*D**2/4.
TV1=5.6322*P1**0.25+273.
ROV1=1.213E-5*((TV1-273.)/5.6322)**3.75
ROL1=RO(TL1)
X1=ROV1/(ROL1-ROV1)
GVB1=GSM*X1
GL1=GSM-GVB1
PB1=P1+4.*SIG(TL1)/DB1
TVB1=5.6322*PB1**0.25+273.
ROVB1=1.213E-5*((TVB1-273.)/5.6322)**3.75
WSM1=2.*GSM/ROL1/F
ISTAR1=GVB1*(VI(TVB1)+WSM1**2/2.)+GL1*(YIL(TL1)
& +WSM1**2/2.)
SB=6.*GSM*X1/PI/(DB1**3)/ROVB1
ROSG=ROL1*(1.-PI*(2.**0.5)/6.)
TL=TL1
J=0
WRITE(8,*)' Z WSM P ROSM PB DB
& FV FL XV ROV ROSG ROVB TV TVB TL
& R0 DLT'
DO 2 I=1,201
J=J+1
Z=(I-1)*1.E-4
312
DB=DB1*(1.-Z/Z0)
WSM=GSM*(1.+(1.-Z/Z0)**3)/ROL1/F
ROSM=GSM/F/WSM
FV=PI*DB**3*SB/WSM/6.
FL=F-FV
P=P1+GSM*(WSM1-WSM)/F
TV=5.6322*P**0.25+273.
ROV=1.213E-5*((TV-273.)/5.6322)**3.75
PB=P+4.*SIG(TL)/DB
TVB=5.6322*PB**0.25+273.
ROVB=1.213E-5*((TVB-273.)/5.6322)**3.75
GVB=PI*DB**3*ROVB*SB/6.
GL=GSM-GVB
XV=GVB/GSM
ALFA=FV/F
TL=((ISTAR1-(VI(TVB)+WSM**2/2.)*GVB)/GL-
WSM**2/2.)/CPL(TL)+273.
PSTAR=P+ROSM*WSM**2/2.
R0=PB*AMV/ROVB/TVB
DLT=TL-TL1
IF(J.EQ.10)THEN
WRITE(8,200)Z,WSM,P,ROSM,PB,DB,FV,FL,XV,ROV,ROSG,ROVB,
TV,TVB,TL,
& R0,DLT
J=0
END IF
IF(Z.EQ.19.9E-3) GOTO 3
2 CONTINUE
200 FORMAT(F6.3,F8.2,F10.0,F8.3,F10.0,3E10.3,F6.3,3F8.3,
&3F6.0,F8.1,F8.3)
3 STOP
END
Файл исходных данных
0.415 20.E-3 4.5E-3 3.1416 36.E3 8.E-5
346.36 18.
GSM Z0 D PI P1 DB1
TL1 AMV
PROGRAM SCC2PI
DIMENSION B(10),C(10)
REAL LAMLS,IMP1,IMP,MAX
313
VIS(A)=(613.-4.227*(A-273.)+9.22E-3*(A-273.)**2)*1E-6
RO(A)=999.9-0.1475*(A-273.)-2.695E-3*(A-273.)**2
YIL(A)=(12.475+3.935*(A-273.)+1.3125E-3*(A-273.)**2)*1.E3
SIG(A)=(75.54-0.168*(A-273.))*1.E-3
SQV(A)=(2447.7-1.283*(A-273.)-6.26E-3*(A-273.)**2)*1.E3
CPL(A)=4275.5-1.2209*(A-273.)+0.59823E-2*(A-273.)**2
VI(A)=(2461.1+2.636*(A-273.)-4.88E-3*(A-273.)**2)*1.E3
PS(A)=1.E6*EXP(82.86568+0.01028003*A-7821.541/A
& -11.48776*ALOG(A))
OPEN (2,FILE='SC2PI.DAT')
OPEN (8,FILE='REZSC2PI.DAT')
READ(2,*) N,DX,GSM,TL1,P1,DB1,
& DC1,DCM,DCSR,ZCC,ZCM,ZCD,UGC,AMV,PI,B1,B2,
& LAMLS,KK
N=2
Z=ZCC+ZCM+ZCD
FC1=PI*DC1**2/4.
ROL1=RO(TL1)
SIGL1=SIG(TL1)
PB1=P1+4.*SIGL1/DB1
TVB1=5.6322*PB1**0.25+273.
ROVB1=1.213E-5*((TVB1-273.)/5.6322)**3.75
XV1=ROVB1/(ROL1-ROVB1)
GVB1=XV1*GSM
ROSM1=ROVB1*ROL1/(ROVB1+GVB1*(ROL1-ROVB1)/GSM)
WSM1=GSM/ROSM1/FC1
IMP1=GSM*WSM1+FC1*P1
SB=6.*GSM*XV1/PI/(DB1**3)/ROVB1
GL1=GSM-GVB1
ISTAR1=GVB1*(VI(TVB1)+WSM1**2/2.)+GL1*(YIL(TL1)+
& WSM1**2/2.)
P=P1
ROL=ROL1
TL=TL1
DB=DB1
WSM=WSM1
ROSM=ROSM1
GVB=GVB1
XV=XV1
FC=FC1
IMP=IMP1
GCV=UGC/AMV
314
IS=0
B(1)=DB
B(2)=IMP
J=0
K=0
X=0
WRITE(8,*)' X P WSM GVB TL ROSM D
&TPL DB TVB XV AV ALFA ASM MAX RJDL'
120 X=X+DX
K=K+1
DO 110 LL=1,4
SIGL=SIG(TL)
VISL=VIS(TL)
ROL=RO(TL)
PB=P+4.*SIGL/DB
CAP=B2/(PB**0.56)
ROVB=1.213E-5*PB**0.9375
TVB=5.6322*PB**0.25+273.
GVB=SB*PI*DB**3*ROVB/6.
XV=GVB/GSM
ROSM=ROVB*ROL/(ROVB+GVB*(ROL-ROVB)/GSM)
IF(X.LT.ZCC) DC=DC1-(DC1-DCM)*X/ZCC
IF(X.GT.ZCC.AND.X.LE.(ZCC+ZCM)) DC=DCM
IF(X.GT.(ZCC+ZCM)) DC=DCM+(DCSR-DCM)*(X-
(ZCC+ZCM))/ZCD
FC=PI*DC**2/4.
IF(X.LT.ZCC)PRFC=-PI*(DC1-(DC1-DCM)*X/ZCC)*(DC1-
DCM)/2./ZCC
IF(X.GE.ZCC.AND.X.LE.(ZCC+ZCM))PRFC=0.
IF(X.GT.(ZCC+ZCM))PRFC=PI*(DCM+(DCSR-DCM)*(X-
(ZCC+ZCM)
DTPL=TL-TVB
SQ=SQV(TVB)
CPLIT=CPL(TL)
WSM=GSM/ROSM/FC
REC=WSM*DC*ROSM/VISL
CF=0.0791/(REC**0.25)
TW=CF*ROSM*WSM**2/2.
GL=GSM-GVB
TL=((ISTAR1-(VI(TVB)+WSM**2/2.)*GVB)/GL-
WSM**2/2.)/CPLIT+273.
CAP=B2/(P**0.56)
315
RJDL=CAP*(PB/(TVB**0.5)-
PS(TL)/(TL**0.5))/((2.*PI*GCV)**0.5)
Y=(TL-273.)/100.
AV=(1.25*PB/ROVB)**0.5+13.36363636-0.03636363636*(TVB-
273.)
AL=1429.00787+595.633911*Y-749.205444*Y**2+394.216553*
& Y**3-109.585632*Y**4+11.354808*Y**5
ALFA=(ROL-ROVB*ROL/(ROVB+XV*(ROL-ROVB)))/(ROL-
ROVB)
ASM=1./(((ALFA*ROVB+(1.-
ALFA)*ROL)*(ALFA/(ROVB*AV**2)+
& (1.-ALFA)/(ROL*AL**2)))**0.5)
DALFDP=ALFA*(1.-ALFA)*(1./(ROL*AL**2)-
1./(ROVB*AV**2))
ASM1=1./(((ROVB-ROL)*DALFDP+ALFA/(AV**2)+(1.-ALFA)/
& (AL**2))**0.5)
ASM2=AV*(ROVB/(ALFA*(ROL-(ROL-ROVB)*ALFA)+(1.-
ALFA)*ROVB*
& AV**2/(AL**2)))**0.5
ALFMIN=(1-ROVB*AV**2/(ROL*AL**2))/(2.*(1.-ROVB/ROL))
MAX=WSM/ASM
OPR=WSM**2/(ASM**2)-1.
C(1)=8.*LAMLS*CPLIT*RO(TVB)*DTPL*ABS(DTPL)*B1/
& (SQ**2*ROVB**2*DB*WSM)
C(2)=-TW*PI*DC
CALL RKM(A,DX,N,IS,B,C)
DB=B(1)
110 IMP=B(2)
P=(IMP-GSM**2/FC/ROSM)/FC
WSM=GSM/FC/ROSM
PSTAR=P+ROSM*WSM**2/2.
IF(K.EQ.KK)THEN
WRITE(8,200) X,P,WSM,GVB,TL,ROSM,DTPL,DB,TVB,XV,AV,
& ALFA,ASM,MAX,RJDL
J=J+1
K=0
IF(J.EQ.10)THEN
J=0
WRITE(8,*)' X P WSM GVB TL ROSM D
&TPL DB TVB XV AV ALFA ASM MAX RJDL'
END IF
END IF
316
IF(GVB.EQ.0) GOTO 16
IF(X.LT.Z) GOTO 120
IF(X.GT.Z) GOTO 16
200 FORMAT(F6.4,F10.0,F6.1,E10.3,F6.0,F8.3,E12.3,E10.3,
& F6.0,E10.3,F6.0,F6.3,F8.1,2F7.3)
16 STOP
END
Файл исходных данных
2 0.1E-3 0.415 346.36 36.E3 8.E-5
4.5E-3 4.5E-3 1.E-2 0.0 1.E-2 0.05
8314. 18. 3.1416 0.00195 1. 0.68 2
N DX GSM TL0 P0 DB0
DC0 DCM DCSR ZCC ZCM ZCD
UGC AMV PI B1 B2 LAMLS KK
PROGRAM TAU_BUBBLE
C Программа позволяет определить время схлопывания
C пузырька под действием заданной разности давлений в пу-
C зырьке и воде (Рп меньше Рж).
C
ROL(A)=1008.7-0.288*(A-273.)-2.16E-3*(A-273.)**2
P(T0)=1.E6*EXP(82.86568+0.01028003*T0-7821.541/T0
& -11.48776*ALOG(T0))
F(C)=C**1.5*((1.-C)/(1.-C**3))**0.5
OPEN(2,FILE='TAU.DAT')
OPEN(4,FILE='PRN')
OPEN(6,FILE='CON')
OPEN(8,FILE='TAUREZ.DAT')
READ(2,*) RB0,DLP,T0
ROV=1.1203E-5*PB0**0.9435
PL=P(T0)
PB0=PL-DLP
C0=RB0*(3.*ROL(T0)/DLP/2.)**0.5
RB=2E-6
RBOT=RB/RB0
PRINT*, ' DLP=',DLP,' PL=',PL,' C0=',C0,' RBOT=',RBOT
A=RBOT
B=0.999
S1=F(A)
N=100
D=(B-A)/N
S=(F(A)+F(B))/2
317
X=A
DO 2 I=1,99
X=X+D
2 S=S1+X**1.5*((1.-X)/(1.-X**3))**0.5
S=S*D
TAU=S*C0
WRITE(8,*) ' TAU=',TAU,' C0=',C0,' PL=',PL,' S=',S
STOP
END
Файл TAU.DAT
2.E-3 0.5E5 373.
RB0 DLP T0
Файл TAUREZ.DAT
TAU= 1.9325557E-06 C0= 3.3911061E-04 PL= 100779.5
S= 5.6988946E-06
PROGRAM CAV_1
C Программа позволяет выполнить расчет
C роста пузырька под действием инерционных, поверхностных
C и вязкостных сил при заданной разности давлений в пу-
C зырьке и воде (Рп больше Рж).
Идентификаторы программы
N- число дифференциальных уравнений системы
DTAU- шаг интегрирования по времени
τ
Δ
,с
DLP- перепад давления между пузырьком и жидкостью
RB- радиус пузырька
SIG- коэффициент поверхностного натяжения жидкости, Н/м
TAU- время, с
TL- температура жидкости, К
UB- скорость перемещения поверхности пузырька, м/с
VISL- коэффициент абсолютной вязкости жидкости, Па⋅с
PI- число “пи”
AMV- молекулярная масса пара,
μ
г
, кг/кмоль
UGC-универсальная газовая постоянная, 8314 Дж/(кмоль К)
KK- интервал для выдачи результатов расчета (через каждые КК
шагов)
DIMENSION B(10),C(10)
REAL IVTV,ISM, LAMV,LAM,IV,LBFI
YILTL(A)=(12.475+3.935*(A-273.)+1.3125E-3*(A-273.)**2)*1.E3
SIG(A)=(75.54-0.168*(A-273.))*1.E-3
318
IV(A)=(2461.1+2.636*(A-273.)-4.88E-3*(A-273.)**2)*1.E3
CPL(A)=4331.-3.09*(A-273.)+1.98E-2*(A-273.)**2
VISL(A)=(613.-4.227*(A-273.)+9.22E-3*(A-273.)**2)*1.E-6
VISP(A)=(8.35+3.42E-2*(A-273.)+2.E-5*(A-273.)**2)*1.E-6
SQV(A)=(2447.7-1.283*(A-273.)-6.26E-3*(A-273.)**2)*1.E3
PR(A)=1.28-4.4E-3*(A-273.)+2.4E-5*(A-273.)**2
LAM(A)=(1.801+2.69E-3*(A-273.)+3.02E-5*(A-273.)**2)*1.E-2
LAML(A)=0.57225+1.92E-3*(A-273.)-8.1E-6*(A-273.)**2
CVL(A)=4244.1-1.583*(A-273.)+1.314e-2*(A-273.)**2
CVV(A)=135950-799.5*(A-273.)+1.502*(A-273.)**2
ROP(A)=6.7962-0.129299*(A-273.)+6.7314e-4*(A-273.)**2
PS(A)=1.E6*EXP(82.86568+0.01028003*A-7821.541/
& A-11.48776*ALOG(A))
RO(A)=999.9-0.1475*(A-273.)-2.695E-3*(A-273.)**2
OPEN (2,FILE='CAV.DAT')
OPEN (4,FILE='GRAF.DAT')
OPEN (8,FILE='REZCAV.DAT')
READ(2,*)
N,TAU0,DTAU,DLP,UB0,RB0,TL,UGC,AMV,PI,TAUM,KK
GCA=UGC/AMV
N=2
IS=0
UB=UB0
RB=RB0
TAU=TAU0
B(1)=URB
B(2)=RB
J=0
K=0
120 TAU=TAU+DTAU
K=K+1
DO 110 LL=1,4
C(1)=DLP/RO(TL)/RB-3.*C(2)**2/RB/2.-2.*SIG(TL)/(RO(TL)*
& RB**2)-4.*VISL(TL)*C(2)/(RO(TL)*RB**2)
C(2)=URB
110 CALL RKM(TAU,DTAU,N,IS,B,C)
URB=B(1)
RB=B(2)
IF(K.EQ.KK) THEN
PRINT *, " TAU=",TAU," URB=",URB," RB=",RB
K=0
J=J+1
319
IF(J.EQ.10) THEN
J=0
ENDIF
END IF
WRITE(4,*) TAU,URB,RB
WRITE(8,*) ' TAU=',TAU,' URB=',URB,' RB=',RB ,
& ' TAUM=',TAUM,' C(1)=',C(1),' C(2)=',C(2)
IF(RB.LT.1E-5) GOTO 130
IF(TAU.LT.TAUM) GOTO 120
130 STOP
END
SUBROUTINE RKM(A,DX,N,IS,B,C)
DIMENSION F(10),G(40),C(10),B(10)
IS=IS+1
GOTO(10,30,60,80),IS
10 E=A
DO 20 I=1,N
F(I)=B(I)
G(4*I-3)=C(I)*DX
20 B(I)=F(I)+G(4*I-3)/2.
GOTO 50
30 DO 40 I=1,N
G(4*I-2)=C(I)*DX
40 B(I)=F(I)+G(4*I-2)/2.
50 A=E+DX/2.
GOTO 100
60 DO 70 I=1,N
G(4*I-1)=C(I)*DX
70 B(I)=F(I)+G(4*I-1)
A=E+DX
GOTO 100
80 DO 90 I=1,N
G(4*I)=C(I)*DX
B(I)=G(4*I-3)+2.*(G(4*I-2)+G(4*I-1))
90 B(I)=(B(I)+G(4*I))/6.+F(I)
IS=0
100 RETURN
END
Файл CAV.DAT
2 0.0 0.1E-6 0.4E5 0.0 1.E-5 293. 8314. 18. 3.1416 0.01 10
N TAU0 DTAU DLP URB0 RB0 TL UGC AMV PI TAUM KK