39
CALL OMTAU(RO,T)
IF(Z.EQ.0D0) GO TO 134
Z=1.D0+AO
134 RETURN
END
С Подпрограмма, реализующая итерационный процесс определения
С плотности из уравнения состояния (метод Ньютона)
SUBROUTINE FUN(X)
IMPLICIT REAL*8(A-H,О-Z)
COMMON/P/P/RM/RM/T/T/AI/AO,A1
ITER=1
1 CONTINUE
NPRIZ=0
IF(ITER.NE.l) NPRIZ=1
CALL COMPL(X,T,NPRIZ)
Z=1.D0+AO
FX=1.D6*(P-(1.D-3*RM*T*Z*X))
F=1.D3*RM*T*(1.D0+A1)
DR=FX/F
X=X+DR
IF(ITER.GT.10) GO TO 4
ITER=ITER+1
IF(DABS(DR/X).GT.1.D-6) GO TO 1
4 CALL COMPL(X,T,NPRIZ)
RETURN
END
SUBROUTINE OMTAU(RO,T)
IMPLICIT REAL*8(A-H,O-Z)
COMMON/PARCM/TCM,VCM/Z/Z
Z=-1D0
TR=T/TCM
ROR=RO*VCM
IF(TR.LT.1.05D0) Z=0D0
IF(ROR.LT.0.D0.OR.ROR.GT.3.D0) Z=0D0
RETURN
END
SUBROUTINE COMPL(RO,T,NPRIZ)
IMPLICIT REAL*8(A-H,O-Z)
DIMENSION B(10,8),BK(10)
COMMON/PARCM/TCM,VCM/B/B/AI/AO,A1
IF(NPRIZ.NE.0) GO TO 7
TR=T/TCM
DO 1 I=1,10
BK(I)=0
DO 1 J=1,8
1 BK(I)=BK(I)+B(I,J)/TR**(J-1)
7 ROR=RO*VCM
AO=0.D0
A1=0.D0
DO 33 I=1,10
D=BK(I)*ROR**I
AO=AO+D