33
DO 5 I=1,NC-1
DO 5 J=I+1,NC
UM=UM+2.*Y(I)*Y(J)*(UIJ(I,J)**5-1D0)*(EI(I)*EI(J))**2.5
GM=GM+2.*Y(I)*Y(J)*(GIJ(I,J)-1D0)*(GI(I)+GI(J))
5 KM=KM+2.*Y(I)*Y(J)*(KIJ(I,J)**5-1D0)*(KI(I)*KI(J))**2.5
KM=KM**.6
UM=UM**.2
DO 7 N=1,13
B1(N)=0D0
DO 9 I=1,NC
9 В1(N)=B1(N)+Y(I)*Y(I)(GIJM(I,I)+ 1D0-GN(N))**GN(N)*
*(QI(I)*QI(I)+1D0-QN(N))**QN(N)*(FI(I)+1D0-FN(N))*FN(N)*
*EIJM(I,I)"UN(N)*KI(I)*KI(I)*KI(I)
DO 11 I=1,NC-1
DO 11 J=I+1,NC
11 В1(N)=B1(N)+2.*Y(I)*Y(J)(GIJМ(I,J)+1D0-GN(N))**GN(N)*
*(QI(I)*QI(J)+1D0-QN(N))**QN(N)*((FI(I)*FI(J))**.5+
1D0-FN(N))**FN(N)*EIJM(I,J)**UN(N)*(KI(I)*KI(J)**1.5
7 CONTINUE
DO 13 N=8,53
13 C1(N)=AN(N)*(GM+1D0-GN(N))**GN(N)*(QM**2+1D0-QN(N))**
*QN(N)*(FM+1D0-FN(N))**FN(N)*UM**UN(N)
RETURN
END
SUBROUTINE PHASE1
IMPLICIT REAL*8(A-H,O-Z)
COMMON/Z/Z/RM/RM/T/T/P/P/AI1/AO,A1/AN/AN(53)
*/COEF1/B1(13),C1(53)/COEF2/B,C(53)/UN/UN(53)
CALL PCONT1(P,T)
IF(Z.EQ.0D0) GO TO 134
B=0D0
DO 1 N=1,13
1 B=B+AN(N)/T**UN(N)*B1(N)
DO 3 N=8,53
3 C(N)=C1(N)/T**UN(N)
PR=P/5.
RO=9D3*P/(RM*T*(1.1*PR+0.7D0))
CALL FUN1(RO)
Z=1D0+AO
134 RETURN
END
С Подпрограмма, реализующая итерационный процесс определения
С плотности из уравнения состояния (метод Ньютона)
SUBROUTINE FUN1(X)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/P/P/RM/RM/T/T/AI1/AO,A1
ITER=1
1 CONTINUE
CALL COMPL1(X)
Z=1.D0+AO
FX= 1 .D6*(P-(1.D-3*RM*T*Z*X))
F= 1 .D3*RM*(1.D0+A1)