1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
|
SUBROUTINE ALG04(H,S,VW,R1,R2,X1,X2,VM,EPS,SCLFAC,G,EJ,HMIN,VMIN,
1PSMID,NSTRMS,LOG2,LNCT,IFAIL)
C
DIMENSION H(1),S(1),VW(1),R1(1),R2(1),X1(1),X2(1),VM(1)
DIMENSION VZUP(21),VZDN(21),PSDN(21),HUP(21),VWUP(21),SUP(21),RMID
1(20),DELR(20),PSUP(21),HSUP(21),VWFUN(21),VZFUN(21),SDN(21),HSDN(2
21),VWDN(21),HDN(21),XX1(21),XX2(21),R(21)
C
DO 80 J=1,NSTRMS
80 R(J)=(R1(J)+R2(J))*0.5
Q1=R(NSTRMS)-R(1)
Q2=VM(1)
DO 90 J=2,NSTRMS
IF(R(J)-R(J-1).LT.Q1)Q1=R(J)-R(J-1)
IF(VM(J).LT.Q2)Q2=VM(J)
90 CONTINUE
DELZ=Q2*Q1**2/(EPS*SCLFAC)*0.25
Q1=(X2(1)+X2(NSTRMS)-X1(1)-X1(NSTRMS))*0.5
ISTEP=Q1/DELZ+1.0
DELZ=Q1/FLOAT(ISTEP)
VM2=VMIN**2
ITUB=NSTRMS-1
IMID=NSTRMS/2+1
DO 110 J=1,NSTRMS
PSUP(J)=PSMID
HUP(J)=H(J)
VWUP(J)=VW(J)
110 SUP(J)=S(J)
DO 120 J=1,ITUB
RMID(J)=(R(J)+R(J+1))*0.5
120 DELR(J)=R(J+1)-R(J)
IFAIL=0
KSTEP=1
130 CALL ALG29(VWUP,R,XX2,NSTRMS)
DO 150 J=1,NSTRMS
150 VWFUN(J)=EPS/R(J)*(XX2(J)-VWUP(J)/R(J))*SCLFAC
IF(KSTEP.GT.1)GO TO 280
JSTEP=1
J1=IMID
190 J2=J1+JSTEP
JJ=J1
IF(JSTEP.EQ.-1)JJ=J2
Q1=((VWUP(J1)+VWUP(J2))*0.5)**2/RMID(JJ)
Q1=DELR(JJ)*Q1*FLOAT(JSTEP)
X3=(SUP(J1)+SUP(J2))*0.5
K=1
200 Q2=ALG2(X3,(PSUP(J1)+PSUP(J2))*0.5)
IF(Q2.GE.HMIN)GO TO 210
IFAIL=1
GO TO 600
210 Q2=ALG5(Q2,X3)/G
X4=PSUP(J2)
PSUP(J2)=PSUP(J1)+Q1*Q2
IF(ABS(X4/PSUP(J2)-1.0).LE.1.0E-5)GO TO 220
K=K+1
IF(K.LE.10)GO TO 200
IFAIL=2
GO TO 600
220 IF(J2.EQ.1)GO TO 240
IF(J2.EQ.NSTRMS)GO TO 230
J1=J2
GO TO 190
230 JSTEP=-1
J1=IMID
GO TO 190
240 DO 260 J=1,NSTRMS
HSUP(J)=ALG2(SUP(J),PSUP(J))
IF(HSUP(J).GE.HMIN)GO TO 250
IFAIL=3
GO TO 600
250 Q1=2.0*G*EJ*(HUP(J)-HSUP(J))-VWUP(J)**2
IF(Q1.GE.VM2)GO TO 260
IFAIL=4
GO TO 600
260 VZUP(J)=SQRT(Q1)
FLOW=0.0
DO 270 J=1,ITUB
270 FLOW=FLOW+(R(J+1)**2-R(J)**2)*(VZUP(J)+VZUP(J+1))*ALG5((HSUP(J)+HS
1UP(J+1))*0.5,(SUP(J)+SUP(J+1))*0.5)
280 CALL ALG29(HSUP,R,XX2,NSTRMS)
DO 300 J=1,NSTRMS
300 HDN(J)=HUP(J)+DELZ/VZUP(J)*EPS/R(J)*XX2(J)*SCLFAC
DO 310 J=1,NSTRMS
310 VWDN(J)=VWUP(J)+DELZ/VZUP(J)*VWFUN(J)
CALL ALG29(VZUP,R,VZFUN,NSTRMS)
DO 330 J=1,NSTRMS
VZFUN(J)=DELZ*EPS*SCLFAC*VZFUN(J)/R(J)
SDN(J)=SUP(J)
330 PSDN(J)=PSUP(J)
KK=1
340 J1=IMID
JSTEP=1
350 J2=J1+JSTEP
JJ=J1
IF(JSTEP.EQ.-1)JJ=J2
Q1=((VWDN(J1)+VWDN(J2))*0.5)**2/RMID(JJ)
Q1=DELR(JJ)*Q1*FLOAT(JSTEP)
X3=(SDN(J1)+SDN(J2))*0.5
K=1
360 Q2=ALG2(X3,(PSDN(J1)+PSDN(J2))*0.5)
IF(Q2.GE.HMIN)GO TO 370
IFAIL=5
GO TO 600
370 Q2=ALG5(Q2,X3)/G
X4=PSDN(J2)
PSDN(J2)=PSDN(J1)+Q1*Q2
IF(ABS(X4/PSDN(J2)-1.0).LE.1.0E-5)GO TO 380
K=K+1
IF(K.LE.10)GO TO 360
IFAIL=6
GO TO 600
380 IF(J2.EQ.1)GO TO 400
IF(J2.EQ.NSTRMS)GO TO 390
J1=J2
GO TO 350
390 J1=IMID
JSTEP=-1
GO TO 350
400 DO 410 J=1,NSTRMS
VZDN(J)=VZUP(J)+(VZFUN(J)-(PSDN(J)-PSUP(J))/ALG5(HSUP(J),SUP(J))*G
1)/VZUP(J)
HSDN(J)=HDN(J)-(VZDN(J)**2+VWDN(J)**2)/(2.0*G*EJ)
IF(HSDN(J).GE.HMIN)GO TO 410
IFAIL=7
GO TO 600
410 SDN(J)=ALG3(PSDN(J),HSDN(J))
XX1(1)=0.0
DO 420 J=1,ITUB
420 XX1(J+1)=XX1(J)+(R(J+1)**2-R(J)**2)*(VZDN(J+1)+VZDN(J))*ALG5((HSDN
1(J)+HSDN(J+1))*0.5,(SDN(J)+SDN(J+1))*0.5)
Q1=XX1(NSTRMS)
IF(ABS(Q1/FLOW-1.0).LE.1.0E-5.AND.KK.GT.1)GO TO 450
IF(KK.LE.15)GO TO 430
IFAIL=8
GO TO 600
430 Q2=ALG9(HSDN(IMID),SDN(IMID),VZDN(IMID)**2)
Q1=(Q1-FLOW)*PSDN(IMID)*Q2/(FLOW*(1.0-Q2))
DO 440 J=1,NSTRMS
440 PSDN(J)=PSDN(J)+Q1
KK=KK+1
GO TO 340
450 IF(KSTEP.EQ.ISTEP)GO TO 510
DO 500 J=1,NSTRMS
PSUP(J)=PSDN(J)
HSUP(J)=HSDN(J)
VZUP(J)=VZDN(J)
VWUP(J)=VWDN(J)
HUP(J)=HDN(J)
500 SUP(J)=SDN(J)
KSTEP=KSTEP+1
GO TO 130
510 DO 520 J=1,NSTRMS
H(J)=HDN(J)
S(J)=SDN(J)
520 VW(J)=VWDN(J)
RETURN
600 CALL ALG03(LNCT,1)
WRITE(LOG2,610)IFAIL
610 FORMAT(5X,30HMIXING CALCULATION FAILURE NO.,I2)
RETURN
END
|
