C-----------------------------------------------------------------
C    NS-BWR Equation
C    This program calculates properties of R22-R114 mixture
C    at buble and dew point at P=2.5[bar] and Z=0.5[kg/kg].
C-----------------------------------------------------------------
C
C---------------------------------------------
C    Definition of array PR1,PR2,CP1 and CP2.
C---------------------------------------------
      DIMENSION PR1(5),PR2(5),CP1(9),CP2(9)
C
C----------------------------------------
C    T=[C],P=[bar],quantity=[kg]
C    standard state=(saturated liquid at 273.15[K])
C-----------------------------------------
      call KPAMES(1,1)
      call STNKAS(1,1)
C
C----------------------------------------------------------------------------
C    Putting fundamental constants into PR1 and PR2 using DATA STATEMENT.
C    DATA(PR1(I),I=1,5)/MW1[KG/KMOL],TC1[K],dummy, VC1[M^3/KMOL],OMEGA1[-]/
C    DATA(PR2(I),I=1,5)/MW2[KG/KMOL],TC2[K],dummy, VC2[M^3/KMOL],OMEGA2[-]/
C----------------------------------------------------------------------------
      DATA(PR1(I),I=1,5)/ 86.469, 369.3,  -1.0, 0.168555, 0.205/
      DATA(PR2(I),I=1,5)/170.922, 418.78, -1.0, 0.296739, 0.234/
C
C----------------------------------------------------------------------------
C    Putting isobaric heat capacity at ideal gas state into CP1 and CP2
C    using DATA STATEMENT.
C    DATA(CP1(I),I=1,8)/EQUATION NUMBER,TMIN[K],TMAX[K], A,B,C,D,E/
C    DATA(CP2(I),I=1,8)/EQUATION NUMBER,TMIN[K],TMAX[K], A,B,C,D,E/
C----------------------------------------------------------------------------
      DATA(CP1(I),I=1,9)
     &  /3.0, 250.0, 350.0,
     &   -917.3,28.755,-2.030E-02,2.5351E-03,-2.185E-06, 7.1E-10/
C        [kJ/kg]
      DATA(CP2(I),I=1,9)
     &  /2.0, 250.0, 400.0,
     &   -0.12916, 4.7009E-03, -8.612E-06, 6.05E-09, 0.0, 0.0 /
C        [kJ/kg]
      do 10 I=4,9
C        [kJ/kg] -> [J/kmol]
         CP1(I) = CP1(I) * PR1(1) * 1000
         CP2(I) = CP2(I) * PR2(1) * 1000
 10   continue
C
C------------------------------------------
C    Setting interaction parametar AKIJ.
C------------------------------------------
      AKIJ=0.976
C
C--------------------------
C    Execution START2.
C--------------------------
      call START2(ECODE, PR1, PR2, CP1, CP2, AKIJ)
C
C--------------------------------------------
C   Calculation of properties at buble point and
C   dew point using SUBPB and output results.
C--------------------------------------------
      P=2.5
      Z=0.5
      CALL SUBTB(J,TB,P,Z,VB,HB,SB)
      CALL SUBTD(J,TD,P,Z,VD,HD,SD)
      WRITE(*,*) ' P[bar]     =',P
      WRITE(*,*) ' Z[kg/kg]   =',Z
      WRITE(*,1) TB, TD
      WRITE(*,2) VB, VD
      WRITE(*,3) HB/1000, HD/1000
      WRITE(*,4) SB/1000, SD/1000
 1    format('  TB[C]      =', F7.1,  '    TD =', F7.1)
 2    format('  VB[m^3/kg] =', E10.2, ' VB =', E10.2)
 3    format('  HB[kJ/kg]  =', F7.1,  '    HD =', F7.1)
 4    format('  SB[kJ/kg]  =', F7.3,  '    SD =', F7.3)
      STOP
      END
C     
Cกกกก result
C     
C  P[bar]     =  2.5
C  Z[kg/kg]   =  0.5
C  TB[C]      =  -14.0    TD =    3.5
C  VB[m^3/kg] =  0.72E-03 VB =  0.75E-01
C  HB[kJ/kg]  =  196.7    HD =  372.0
C  SB[kJ/kg]  =  1.025    SD =  1.683