10000 REM HGK+Pitzer et al. for NaCl(aq) 10100 DEFDBL A-H, M-Z 10150 DIM HGKG(40),II(40),JJ(40),BP(10),BQ(10) 10200 DIM ATZ(4),ADZ(4),AAT(4),AAD(4) 10250 DIM BV(10),A(8),C(18) 10300 DIM QR(11),QT(10),QZR(9),QZT(9) 10350 DIM FFD(2),FFP(5),NNT$(2),NND$(2),NNP$(5) 10400 DIM DU(10),ZPIT(60),ZPITZ(120) 10450 GOSUB *BLOCKDATA 10500 GOSUB *UNIT 10850 GOSUB *PARAMETERS 10900 LPRINT 10950 INPUT"Pressure? If end, input 0";X 11000 IF X=0 THEN GOTO 13750 11050 INPUT"Temperature";TT 11100 T=TT 11150 GOSUB *TTTT 11200 T=TTT 11250 IF T=<338.15# THEN PA$="L" ELSE PA$="O" 11300 IF PA$="L" THEN GOTO 11400 11350 FOR I=1 TO 53 : ZPIT(I)=ZPITZ(I) : NEXT I : GOTO 11450 11400 FOR I=1 TO 53 : ZPIT(I)=ZPITZ(I+53) : NEXT I 11450 IF PA$="L" THEN LPRINT "Parameters from low-temperature fits" : GOTO 11550 11500 IF PA$="O" THEN LPRINT "Parameters from overall fits" 11550 RT=GASCON*T 11600 GOSUB *BBT 11650 INPUT "Molality";MOL 11660 IF MOL=<0 THEN GOTO 11650 11700 PRES=X 11750 PINPUT=PRES/FP 11800 DGSS=PINPUT/(T*.4#) 11810 IF T>=647.126# THEN GOTO 12250 11850 DLL=0 11900 DVV=0 12000 DL=DLL:DV=DVV 12050 GOSUB *PCORRTPDLDV 12100 IF ABS((PINPUT-P)/P)=<5D-005 THEN PPP=PINPUT : GOTO 12900 12150 IF PINPUT>P THEN DGSS=DL:GOTO 12250 12200 IF PINPUT

646.3# THEN DOUT=D : GOTO 13100 13000 GOSUB *DFINDDOUTPDTDPD 13050 D=DOUT 13100 GOSUB *THERMDT 13150 DD=DOUT/FD 13200 H=HD*RT 13250 S=SD*GASCON 13300 PRES=PPP*FP 13350 G=GD*RT 13400 LPRINT"Liquid phase at the vapor-saturated pressure" 13450 LPRINT USING"T(& &)=+####.#### P(& &)= +#.#####^^^^^ D(water)=+#.#####^^^^^(& &)"; NT$,TT,NP$,PRES,DD,ND$ 13500 GOSUB *SECDERIVP 13550 GOSUB *DEBYEHUCKEL 13600 GOSUB *NACLAQ 13650 INPUT"Will you continue the calculation? Input Y(or y) or N(or n)";CAL$ 13700 IF CAL$="Y" OR CAL$="y" THEN LPRINT : LPRINT : PRINT : GOTO 10900 13750 END 13800 *DFINDDOUTPDTDPD 13850 DD=D 14000 LL=0 14050 LL=LL+1 14100 IF DD=<1D-008 THEN DD=1D-008 14150 IF DD>1.9# THEN DD=1.9# 14200 D=DD 14250 GOSUB *QQTD 14300 Q0=Q 14350 GOSUB *BASEDT 14400 PP=RT*DD*BASEF+Q0 14450 DPD=RT*(Z+Y*DZB)+Q5 : DQ=DPD 14500 IF DPD>0 THEN GOTO 14700 14550 IF D>=.2967# THEN DD=DD*1.02 14600 IF D<.2967# THEN DD=DD*.98# 14650 IF LL=<10 GOTO 14050 14700 DPDX=DPD*1.1# 14750 IF DPDX<.1# THEN DPDX=.1# 14800 DP=ABS(1#-PP/PPP) 14850 IF DP<1D-009 THEN GOTO 15250 14900 IF D>.3# AND DP<1D-008 THEN GOTO 15250 14950 IF D>.7# AND DP<1D-007 THEN GOTO 15250 15000 XP=(PPP-PP)/DPDX 15050 IF ABS(XP)>.1# THEN XP=XP*.1#/ABS(XP) 15100 DD=DD+XP 15150 IF DD=<0 THEN DD=1D-008 15200 IF LL=<30 THEN GOTO 14050 15250 DOUT=DD 15300 RETURN 15350 *CORRTPDLDVDELG 15400 IF T>646.3# THEN GOTO 16300 15450 DLIQ=DLL 15500 IF DLL=<0 THEN DLIQ=1.11#-.0004*T 15550 DLL=DLIQ:D=DLIQ 15600 GOSUB *DFINDDOUTPDTDPD 15650 D=DOUT:DL=DOUT 15700 GOSUB *THERMDT 15750 GL=GD 15800 DVAP=DVV 15850 IF DVV=<0 THEN DVAP=PPP/RT 15900 D=DVAP:DVV=DVAP 15950 GOSUB *DFINDDOUTPDTDPD 16000 IF DOUT<5D-007 THEN DOUT=5D-007 16050 D=DOUT:DV=DOUT 16100 GOSUB *THERMDT 16150 GV=GD 16200 DELG=GL-GV 16250 RETURN 16300 PPP=0 16350 IF T>647.126# THEN RETURN 16400 DELG=0 16450 TAUC=.657128#*(1#-T/647.126#)^.325# 16500 DL=.322#+TAUC 16550 DV=.322#-TAUC 16600 D=DV 16650 GOSUB *BASEDT 16700 GOSUB *QQTD 16750 ZB=BASEF 16800 PPP=RT*DV*ZB+Q 16850 RETURN 16900 *BBT 16950 BV(1)=1# 17000 FOR I=2 TO 10 17050 BV(I)=BV(I-1)*TZ/T 17100 NEXT I 17150 B1=BP(1)+BP(2)*LOG(1#/BV(2)) 17200 B2=BQ(1) 17250 B1T=BP(2)*BV(2)/TZ 17300 B2T=0 17350 B1TT=0 17400 B2TT=0 17450 FOR I=3 TO 10 17500 B1=B1+BP(I)*BV(I-1) 17550 B2=B2+BQ(I)*BV(I-1) 17600 B1T=B1T-CDBL((I-2))*BP(I)*BV(I-1)/T 17650 B2T=B2T-CDBL((I-2))*BQ(I)*BV(I-1)/T 17700 B1TT=B1TT+BP(I)*CDBL((I-2))*CDBL((I-2))*BV(I-1)/(T*T) 17750 B2TT=B2TT+BQ(I)*CDBL((I-2))*CDBL((I-2))*BV(I-1)/(T*T) 17800 NEXT I 17850 B1TT=B1TT-B1T/T 17900 B2TT=B2TT-B2T/T 17950 RETURN 18000 *BASEDT 18050 Y=.25#*B1*D 18100 XX=1#-Y 18150 Z0=(1#+ALPHA*Y+BETA*Y*Y)/(XX*XX*XX) 18200 Z=Z0+4#*Y*(B2/B1-GAMMA) 18250 DZ0=(ALPHA+2#*BETA*Y)/(XX*XX*XX)+3#*(1#+ALPHA*Y+BETA*Y*Y)/(XX*XX*XX*XX) 18300 DZB=DZ0+4#*(B2/B1-GAMMA) 18350 AB=(-1#)*LOG(XX)-(BETA-1#)/XX+28.16666667#/(XX*XX)+4#*Y*(B2/B1-GAMMA)+15.166666667#+LOG(D*RT/.101325) 18450 BASEF=Z 18500 BB2TT=T*T*B2TT 18550 UB=(-1#)*T*B1T*(Z-1#-D*B2)/B1-D*T*B2T 18650 CVB=2#*UB+(Z0-1#)*((T*B1T/B1)*(T*B1T/B1)-T*T*B1TT/B1)-D*(BB2TT-GAMMA*B1TT*T*T)-(T*B1T/B1)*(T*B1T/B1)*Y*DZ0 18700 DPDTB=BASEF/T+BASEF*D/Z*(DZB*B1T/4#+B2T-B2/B1*B1T) 18750 SB=UB-AB 18800 RETURN 18850 *QQTD 18900 QR(1)=0 18950 Q5=0 19000 Q=0 19050 AR=0 19100 DADT=0 19150 CVR=0 19200 DPDTR=0 19250 E=EXP((-1#)*AA*D) 19300 Q10=D*D*E 19350 Q20=1#-E 19400 QR(2)=Q10 19450 QV=TZ/T 19500 QT(1)=T/TZ 19550 FOR I=2 TO 10 19600 QR(I+1)=QR(I)*Q20 19650 QT(I)=QT(I-1)*QV 19700 NEXT I 19750 FOR I=1 TO INC 19800 K=II(I)+1 19850 L=JJ(I) 19900 QK=CDBL(K) : QL=CDBL(L) 20000 QZR(K-1)=QR(K+1):QZT(L)=QT(L+1):QZR(K)=QR(K+2):QZT(L+1)=QT(L+2) 20050 QP=HGKG(I)*AA*QZR(K-1)*QZT(L) 20100 Q=Q+QP 20150 Q5=Q5+AA*(2#/D-AA*(1#-E*(QK-1#)/Q20))*QP 20200 AR=AR+HGKG(I)*QZR(K)*QZT(L)/(Q10*QK*RT) 20250 DFDT=Q20^QK*(1#-QL)*QZT(L+1)/(TZ*QK) 20300 D2F=QL*DFDT 20350 DPT=DFDT*Q10*AA*QK/Q20 20400 DADT=DADT+HGKG(I)*DFDT 20450 DPDTR=DPDTR+HGKG(I)*DPT 20500 CVR=CVR+HGKG(I)*D2F/GASCON 20550 NEXT I 20600 QP=0 20650 Q2A=0 20700 FOR J=37 TO 40 20750 IF HGKG(J)=0 THEN GOTO 22000 20800 K=II(J) 20850 KM=JJ(J) 20900 QK=CDBL(K) : QKM=CDBL(KM) 20950 DDZ=ADZ(J-36) 21000 DEL=D/DDZ-1# 21050 IF ABS(DEL)<1D-010 THEN DEL=1D-010 21150 EX1=(-1#)*AAD(J-36)*DEL^QK 21200 DEX=EXP(EX1)*DEL^QKM 21250 ATT=AAT(J-36) 21300 TX=ATZ(J-36) 21350 TAU=T/TX-1# 21400 EX2=(-1#)*ATT*TAU*TAU 21450 TEX=EXP(EX2) 21500 Q10=DEX*TEX 21550 QM=QKM/DEL-QK*AAD(J-36)*DEL^(QK-1#) 21600 FCT=QM*D*D*Q10/DDZ 21650 Q5T=FCT*(2#/D+QM/DDZ)-(D/DDZ)*(D/DDZ)*Q10*(QKM/(DEL*DEL)+QK*(QK-1#)*AAD(J-36)*DEL^(QK-2#)) 21700 Q5=Q5+Q5T*HGKG(J) 21750 QP=QP+HGKG(J)*FCT 21800 DADT=DADT-2#*HGKG(J)*ATT*TAU*Q10/TX 21850 DPDTR=DPDTR-2#*HGKG(J)*ATT*TAU*FCT/TX 21900 Q2A=Q2A+T*HGKG(J)*(4#*ATT*EX2+2#*ATT)*Q10/(TX*TX) 21950 AR=AR+Q10*HGKG(J)/RT 22000 NEXT J 22050 SR=(-1#)*DADT/GASCON 22100 UR=AR+SR 22150 CVR=CVR+Q2A/GASCON 22200 Q=Q+QP 22250 RETURN 22300 *THERMDT 22350 GOSUB *IDEALT 22400 GOSUB *BASEDT 22450 GOSUB *QQTD 22500 QPQ=Q:QDPQ=Q5 22550 Z=BASEF+QPQ/(RT*D) 22600 DPDD=RT*(BASEF+Y*DZB)+QDPQ 22650 AD=AB+AR+AI-UREF/T+SREF 22700 GD=AD+Z 22750 UD=UB+UR+UI-UREF/T 22800 DPDT=RT*D*DPDTB+DPDTR 22850 CVDX=CVB+CVR+CVIX 22900 CPD=CVDX+T*DPDT*DPDT/(D*D*DPDD*GASCON) 22950 HD=UD+Z 23000 SD=SB+SR+SI-SREF 23050 RETURN 23100 *SECDERIVP 23250 D2PDDDD1=0 : D2PDD2=0 : D2PDTDD=0 : D2PDTDT1=0 : D2PDT2=0 : ZX=0 23300 D2PRESIDDDDD1=0 : D2PRESIDDDDT=0 : D2PRESIDDT2=0 : D2PRESIDDD2=0 23310 D2DDT2A=0 : D2DDT2=0 : DDDT=0 23400 D2PDDDD1=3#+ALPHA+3#*Y+4#*ALPHA*Y+3#*BETA*Y+ALPHA*Y*Y+3#*BETA*Y*Y 23450 D2PDDDD1=D2PDDDD1*B1/(2#*XX*XX*XX*XX*XX)+2#*B1*(B2/B1-GAMMA) 23500 D2PDD2=D2PDDDD1*RT 23550 D2PDTDD=Z0+8#*Y*(B2/B1-GAMMA)+(Y+B1T*D*T/2#)*((ALPHA+2#*BETA*Y)/(XX*XX*XX)+3#*Z0/XX) 23600 ZX=6#+3#*ALPHA+BETA+3#*ALPHA*Y+4#*BETA*Y+BETA*Y*Y 23650 D2PDTDD=D2PDTDD+B1T*D*T*Y*ZX/(2#*XX*XX*XX*XX*XX) 23700 D2PDTDD=D2PDTDD+2#*B2T*D*T-2#*B1T*D*T*GAMMA 23750 D2PDTDD=D2PDTDD*GASCON 23800 D2PDTDT1=(ALPHA+2#*BETA*Y)/(XX*XX*XX)+3#*Z0/XX 23850 D2PDTDT1=D2PDTDT1*(2#*B1T*D+B1TT*D*T)/4# 23900 D2PDTDT1=D2PDTDT1+B1T*B1T*D*D*T*ZX/(8#*XX*XX*XX*XX*XX)+2#*B2T*D 23950 D2PDTDT1=D2PDTDT1-2#*B1T*D*GAMMA+B2TT*D*T-B1TT*D*T*GAMMA 24000 D2PDT2=D2PDTDT1*GASCON*D 24050 FOR I=1 TO INC 24100 K=II(I)+1 24150 L=JJ(I) 24200 QK=CDBL(K) : QL=CDBL(L) 24250 D2PRESIDDDDD1=2#/(D*D)-4#/D+4#*(QK-1#)*E/(D*Q20)+1#-3#*(QK-1#)*E/Q20+(QK-1#)*(QK-2#)*E*E/(Q20*Q20) 24300 D2PRESIDDD2=D2PRESIDDD2+D2PRESIDDDDD1*HGKG(I)*QT(L+1)*QR(K+1) 24350 D2PRESIDDDDT=D2PRESIDDDDT-(QL-1#)*HGKG(I)*QT(L+1)*QR(K+1)*(2#/D-1#+(QK-1#)*E/Q20)/T 24400 D2PRESIDDT2=D2PRESIDDT2+(QL-1#)*QL*HGKG(I)*QR(K+1)*QT(L+1)/(T*T) 24450 NEXT I 24500 D2PRESIDDD2A=0 : D2PRESIDDDDTC=0 24550 FOR J=37 TO 40 24560 IF HGKG(J)=0 THEN GOTO 25900 24600 K=II(J) 24650 KM=JJ(J) 24700 QK=CDBL(K) : QKM=CDBL(KM) 24750 DDZ=ADZ(J-36) 24800 DEL=D/DDZ-1# 24850 IF ABS(DEL)<1D-010 THEN DEL=1D-010 24900 EX1=(-1#)*AAD(J-36)*DEL^QK 24950 DEX=EXP(EX1)*DEL^QKM 25000 ATT=AAT(J-36) 25050 TX=ATZ(J-36) 25100 TAU=T/TX-1# 25150 EX2=(-1#)*ATT*TAU*TAU 25200 TEX=EXP(EX2) 25250 Q10=DEX*TEX 25300 QM=QKM/DEL-QK*AAD(J-36)*DEL^(QK-1#) 25350 D2PRESIDDD2A=QM*(2#/(D*D)+4#*QKM/(D*DDZ*DEL)+4#*QK*EX1/(D*DDZ*DEL)) 25400 D2PRESIDDD2B=QM*(QKM*(QKM-1#)+2#*QK*QKM*EX1+QK*(QK-1#)*EX1+QK*QK*EX1*EX1) 25450 D2PRESIDDD2A=D2PRESIDDD2A/DDZ+D2PRESIDDD2B/(DDZ*DDZ*DDZ*DEL*DEL) 25500 D2PRESIDDD2C=(QKM-QK*(QK-1#)*EX1)*(4#/D+2*QKM/(DDZ*DEL)+2#*QK*EX1/(DDZ*DEL))/(DDZ*DDZ*DEL*DEL) 25550 D2PRESIDDD2A=D2PRESIDDD2A-D2PRESIDDD2C 25600 D2PRESIDDD2A=D2PRESIDDD2A+(2#*QKM+QK*(QK-1#)*(QK-2#)*EX1)/(DDZ*DDZ*DDZ*DEL*DEL*DEL) 25650 D2PRESIDDD2=D2PRESIDDD2+HGKG(J)*Q10*D2PRESIDDD2A*D*D 25700 D2PRESIDDDDTC=2#*QM+D*QKM*QM/(DDZ*DEL)+D*QK*EX1*QM/(DDZ*DEL) 25750 D2PRESIDDDDTC=D2PRESIDDDDTC+D*(QK*(QK-1#)*EX1/(DEL*DEL)-QKM/(DEL*DEL))/DDZ 25800 D2PRESIDDDDT=D2PRESIDDDDT-2#*D*HGKG(J)*ATT*TAU*Q10*D2PRESIDDDDTC/(TX*DDZ) 25850 D2PRESIDDT2=D2PRESIDDT2-2#*D*D*HGKG(J)*ATT*(1#+2#*EX2)*Q10*QM/(DDZ*TX*TX) 25900 NEXT J 25950 D2PDD2=D2PDD2+D2PRESIDDD2 26000 D2PDTDD=D2PDTDD+D2PRESIDDDDT 26050 D2PDT2=D2PDT2+D2PRESIDDT2 26100 D2DDT2A=DPDD*DPDD*D2PDT2-2#*DPDT*DPDD*D2PDTDD+DPDT*DPDT*D2PDD2 26150 D2DDT2=(-1#)*D2DDT2A/(DPDD*DPDD*DPDD) 26200 DDDT=(-1#)*DPDT/DPDD 26250 RETURN 26300 *PST 26350 IF T>314# THEN GOTO 26550 26400 PL=6.3573118#-8858.843#/T+607.56335*T^(-.6#) 26450 PS=.1#*EXP(PL) 26500 RETURN 26550 TR=T/647.25# 26600 W=ABS(1#-TR) 26650 BPST=0 26700 FOR I=1 TO 8 26750 ZPST=CDBL(I) 26800 BPST=BPST+A(I)*W^((ZPST+1#)/2#) 26850 NEXT I 26900 QPST=BPST/TR 26950 PS=22.093*EXP(QPST) 27000 RETURN 27050 *IDEALT 27100 TIDEAL=T/100 27150 TL=LOG(TIDEAL) 27200 GI=(-1#)*(C(1)/TIDEAL+C(2))*TL 27250 HI=(C(2)+C(1)*(1#-TL)/TIDEAL) 27300 CPI=C(2)-C(1)/TIDEAL 27350 FOR I=3 TO 18 27400 GI=GI-C(I)*TIDEAL^CDBL(I-6) 27450 HI=HI+C(I)*CDBL((I-6))*TIDEAL^CDBL(I-6) 27500 CPI=CPI+C(I)*CDBL((I-6))*CDBL((I-5))*TIDEAL^CDBL(I-6) 27550 NEXT I 27600 AI=GI-1# 27650 UI=HI-1# 27700 CVIX=CPI-1# 27750 SI=UI-AI 27800 RETURN 27850 *PCORRTPDLDV 27900 GOSUB *PST 27950 PPP=PS 28000 GOSUB *CORRTPDLDVDELG 28050 DP=0 28100 DP=DELG*RT/(1#/DV-1#/DL) 28150 PPP=PPP+DP 28200 IF ABS(DELG)<1D-005 THEN GOTO 28300 28250 DLL=DL:DVV=DV:GOTO 28000 28300 P=PPP 28350 RETURN 28400 *UNIT 28450 PRINT"*******************" 28500 PRINT"* Enter units *" 28550 PRINT"*******************" 28600 PRINT A1$ 28650 PRINT"Choose from 1=deg K, 2=deg C" 28700 INPUT IT 28750 IF IT<1 OR IT>2 THEN GOTO 28650 28800 NT$=NNT$(IT) 28850 PRINT A2$ 28900 PRINT"Choose from 1=kg/m3, 2=g/cm3" 28950 INPUT ID 29000 IF ID>2 OR ID<1 THEN GOTO 28900 29050 ND$=NND$(ID) 29100 FD=FFD(ID) 29150 PRINT A3$ 29200 PRINT"Choose from 1=MPa, 2=bar" 29250 INPUT IP 29300 IF IP>2 OR IP<1 THEN GOTO 29200 29350 NP$=NNP$(IP) 29400 FP=FFP(IP) 29750 RETURN 29800 *TTTT 29850 ON IT GOTO 29900, 30050 29900 TTT=T 30000 GOTO 30150 30050 TTT=T+273.15# 30150 RETURN 30200 *BLOCKDATA 30250 FOR I=1 TO 4:READ ATZ(I):NEXT I 30300 DATA 640#,640#,641.6#,270# 30350 FOR I=1 TO 4:READ ADZ(I):NEXT I 30400 DATA 0.319#,0.319#,0.319#,1.55# 30450 FOR I=1 TO 4:READ AAT(I):NEXT I 30500 DATA 2.0D+004,2.0D+004,4.0D+004,25.0# 30550 FOR I=1 TO 4:READ AAD(I):NEXT I 30600 DATA 34.0#,40.0#,30.0#,1.05D+003 30650 GASCON=.461518# : TZ=647.073 : AA=1# : INC=36 30710 UREF=-4328.454977# : SREF=7.618072# 30750 ALPHA=11#:BETA=44.333333333333#:GAMMA=3.5# 30800 FOR I=1 TO 10:READ BP(I):NEXT I 30850 DATA 0.7478629#,-0.3540782#,0.0#,0.0#,0.007159876#,0.0#,-0.003528426#,0.0#,0.0#,0.0# 30900 FOR I=1 TO 10:READ BQ(I):NEXT I 30950 DATA 1.1278334#,0.0#,-0.5944001#,-5.010996#,0.0#,0.63684256#,0.0#,0.0#,0.0#,0.0# 31000 FOR I=1 TO 40:READ HGKG(I):NEXT I 31050 DATA -5.3062968529023D+002,2.2744901424408D+003,7.8779333020687D+002 31100 DATA -6.9830527374994D+001,1.7863832875422D+004,-3.9514731563338D+004 31150 DATA 3.3803884280753D+004,-1.3855050202703D+004,-2.5637436613260D+005 31200 DATA 4.8212575981415D+005,-3.4183016969660D+005,1.2223156417448D+005 31250 DATA 1.1797433655832D+006,-2.1734810110373D+006,1.0829952168620D+006 31300 DATA -2.5441998064049D+005,-3.1377774947767D+006,5.2911910757704D+006 31350 DATA -1.3802577177877D+006,-2.5109914369001D+005,4.6561826115608D+006 31400 DATA -7.2752773275387D+006,4.1774246148294D+005,1.4016358244614D+006 31450 DATA -3.1555231392127D+006,4.7929666384584D+006,4.0912664781209D+005 31500 DATA -1.3626369388386D+006,6.9625220862664D+005,-1.0834900096447D+006 31550 DATA -2.2722827401688D+005,3.8365486000660D+005,6.8833257944332D+003 31600 DATA 2.1757245522644D+004,-2.6627944829770D+003,-7.0730418082074D+004 31650 DATA -0.225#,-1.68#,0.055#,-93.0# 31660 FOR I=37 TO 40 : HGKG(I)=0 : NEXT I 31700 FOR I=1 TO 40:READ II(I):NEXT I 31750 DATA 0,0,0,0,1,1,1,1,2,2,2,2,3,3,3,3,4,4,4,4,5,5,5,5,6,6,6,6,8,8,8,8,2,2,0,4,2,2,2,4 31800 FOR I=1 TO 40:READ JJ(I):NEXT I 31850 DATA 2,3,5,7,2,3,5,7,2,3,5,7,2,3,5,7,2,3,5,7,2,3,5,7,2,3,5,7,2,3,5,7,1,4,4,4,0,2,0,0 31900 FOR I=1 TO 8:READ A(I):NEXT I 31950 DATA -7.8889166#,2.5514255#,-6.716169#,33.239495# 32000 DATA -105.38479#,174.35319#,-148.39348#,48.631602# 32050 FOR I=1 TO 18:READ C(I):NEXT I 32100 DATA 1.9730271018D+001,2.09662681977D+001,-4.83429455355D-001,6.05743189245D+000 32150 DATA 2.256023885D+001,-9.87532442D+000,-4.3135538513D+000,4.58155781D-001 32200 DATA -4.7754901883D-002,4.1238460633D-003,-2.7929052852D-004 32250 DATA 1.4481695261D-005,-5.6473658748D-007,1.6200446D-008,-3.303822796D-010 32300 DATA 4.51916067368D-012,-3.70734122708D-014,1.37546068238D-016 32350 FOR I=1 TO 2:READ FFD(I):NEXT I 32400 DATA 1.0D-003,1.0# 32450 FOR I=1 TO 2:READ FFP(I):NEXT I 32500 DATA 1.0#,10.0# 32650 FOR I=1 TO 2:READ NNT$(I):NEXT I 32700 DATA "K","deg C" 32750 FOR I=1 TO 2:READ NND$(I):NEXT I 32800 DATA "kg/m3","g/cm3" 32850 FOR I=1 TO 2:READ NNP$(I):NEXT I 32900 DATA "MPa","bar" 33050 A1$="TEMPERATURE":A2$="DENSITY":A3$="PRESSURE" 33100 RETURN 33150 *PARAMETERS 33200 FOR I=1 TO 106 : READ ZPITZ(I) : NEXT I 33250 REM Overall fit (Pitzer, 1987) 33300 DATA -71637.203#,2.2209012#,-7.7991396D-005,-4.8099272D-009,624.68125# 33350 DATA 6.0159787D-004,3.4069074D-007,2.1962044D-011,-110.74702#,0.039494473# 33400 DATA -6.5313475D-007,-6.4781894D-010,-1.5842012D-005,3.2452006D-009,516.99706# 33450 DATA -5.9960301D+006,-656.81518#,24.86912950#,5.381275267D-005,-5.588746990D-008 33500 DATA 6.589326333D-012,-4.4640952#,0.01110991383#,-2.657339906D-007,1.746006963D-010 33550 DATA 1.046261900D-014,-5.307012889D-006,8.634023325D-010,-4.178596200D-013,-1.579365943# 33600 DATA 2.202282079D-003,-1.310550324D-007,-6.381368333D-011,9.706578079#,-2.686039622D-002 33650 DATA 1.534474401D-005,-3.215398267D-009,119.31966#,-0.48309327#,1.4068095D-003 33700 DATA -4.2345814#,-6.1084589#,0.40217793#,2.2902837D-005,-0.075354649# 33750 DATA 1.531767295D-004,-9.0550901D-008,-1.538600820D-008,8.6926600D-011,0.353104136# 33800 DATA -4.3314252D-004,-0.09187145529#,5.1904777D-004 33850 REM Low-temperature fit 33900 DATA -71659.531#,2.3483335#,-8.3668484D-005,2.4018168D-009,624.88208# 33950 DATA -5.3697119D-004,3.5126966D-007,0#,-110.74702#,0.038900801# 34000 DATA 2.6973456D-006,-6.2746876D-010,-1.5267612D-005,0#,516.99706# 34050 DATA -5.9960301D+006,-656.81518#,24.879183#,-2.1552731D-005,5.0166855D-008 34100 DATA 0#,-4.4640952#,0.011087099#,-6.4479761D-008,-2.3234032D-010 34150 DATA 0#,-5.2194871D-006,2.4445210D-010,2.8527066D-013,-1.5696231# 34200 DATA 2.2337864D-003,-6.3933891D-007,4.5270573D-011,5.4151933#,0# 34250 DATA 0#,0#,119.31966#,-0.48309327#,1.4068095D-003 34300 DATA -4.2345814#,-6.1084589#,0.40743803#,-6.8152430D-006,-0.075354649# 34350 DATA 1.2609014D-004,6.2480692D-008,1.8994373D-008,-1.0731284D-010,0.32136572# 34400 DATA -2.5382945D-004,0#,0# 34450 RVGAS=83.144# : RGAS=8.3144# : MW=18.01534 34500 FOR I=1 TO 9: READ DU(I) : NEXT I 34550 DATA 3.4279D+002,-5.0866D-003,9.4690D-007,-2.0525#,3.1159D+003 34600 DATA -1.8289D+002,-8.0325D+003,4.2142D+006,2.1417# 34650 EE=4.803242D-010 : BC=1.380662D-016 34700 NACL=58.4428# : MREF=5.550825 : YPPB=10 34750 RETURN 34800 *DEBYEHUCKEL 34850 PRES=(PRES/FP)*(FFP(2)/FFP(1)) 34950 DPDD=DPDD*(FFP(2)/FFP(1)) 35000 DPDT=DPDT*(FFP(2)/FFP(1)) 35150 EPS=DU(1)*EXP(DU(2)*T+DU(3)*T*T) 35200 E=1#+(PRES-1000)/(DU(7)+DU(8)/T+DU(9)*T+1000) 35250 E=LOG(E) 35300 EPS=EPS+(DU(4)+(DU(5)/(DU(6)+T)))*E 35350 APHI=SQR(2#*3.14159265#*6.022045D+023*D/1000)/3# 35400 APHI=APHI*EE*EE*EE/(BC*SQR(BC)*T*SQR(T)*EPS*SQR(EPS)) 35450 DET=DU(7)+DU(8)/T+DU(9)*T+PRES 35500 DEPSDP=DU(4)+DU(5)/(DU(6)+T) 35550 DEPSDP=DEPSDP/DET 35600 DRHODPDD=1#/(D*DPDD) 35650 AV=2#*RVGAS*T*APHI*(3#*DEPSDP/EPS-DRHODPDD) 35700 ALPH=DPDT/(D*DPDD) 35750 DE=DU(7)+DU(8)/T+DU(9)*T 35800 DEPS=DU(1)*(DU(2)+2#*DU(3)*T)*EXP(DU(2)*T+DU(3)*T*T) 35850 DEPS=DEPS-DU(5)*LOG(1#+(PRES-1000)/(DE+1000))/((DU(6)+T)*(DU(6)+T)) 35900 DEPS=DEPS+(DU(4)+DU(5)/(DU(6)+T))*(1000-PRES)*(DU(9)-DU(8)/(T*T))/((DE+PRES)*(DE+1000)) 35950 DEPS=DEPS/EPS 36000 AH=1#+T*DEPS+T*ALPH/3# 36050 AH=AH*(-6#)*APHI*RGAS*T 36100 D2EPS=DU(1)*(DU(2)+2#*DU(3)*T)*(DU(2)+2#*DU(3)*T)*EXP(DU(2)*T+DU(3)*T*T) 36150 D2EPS=D2EPS+2#*DU(1)*DU(3)*EXP(DU(2)*T+DU(3)*T*T) 36200 D2EPS=D2EPS+(2#*DU(5)/((DU(6)+T)*(DU(6)+T)*(DU(6)+T)))*LOG(1#+(PRES-1000)/(DE+1000)) 36250 D2EPS=D2EPS-(2#*DU(5)/((DU(6)+T)*(DU(6)+T)))*(DU(9)-DU(8)/(T*T))*(1#/(DE+PRES)-1#/(DE+1000)) 36300 D2EPS=D2EPS+(DU(4)+DU(5)/(DU(6)+T))*(2#*DU(8)/(T*T*T))*(1#/(DE+PRES)-1#/(DE+1000)) 36350 D2EPSX=(DU(9)-DU(8)/(T*T))*(DU(9)-DU(8)/(T*T))*(1#/((DE+PRES)*(DE+PRES))-1#/((DE+1000)*(DE+1000))) 36400 D2EPS=D2EPS-(DU(4)+DU(5)/(DU(6)+T))*D2EPSX 36450 D2EPS=D2EPS/EPS 36500 DWDDDT=DDDT/D 36550 DWD2DDT2=D2DDT2/D 36600 AJ=2#*DWD2DDT2-DWDDDT*DWDDDT-2#*DWDDDT/T-6#*D2EPS+15#*DEPS*DEPS+6#*DEPS/T 36650 AJ=AJ-6#*DWDDDT*DEPS+3#/(T*T) 36700 AJ=AJ*APHI*RGAS*T*T 36750 LPRINT USING"APHI= +#.####";APHI 36800 LPRINT USING"AH/RT= +##.###";AH/(RGAS*T) 36850 LPRINT USING"AJ/R= +###.##";AJ/RGAS 36900 LPRINT USING"AV= +##.###";AV 36950 RETURN 37000 *NACLAQ 37050 MIREF=MREF 37051 S=S*MW : H=H*MW : G=G*MW 37052 S=S+SREF*RGAS : H=H+UREF*RGAS : G=G+UREF*RGAS-T*SREF*RGAS 37310 CPW=CPD*RGAS 37400 TL=T-227# 37450 TH=680-T 37500 HFREF=LOG(1#+1.2#*SQR(MIREF))/(2#*1.2#) 37550 BM0V=ZPIT(19)+ZPIT(20)*2#*PRES+ZPIT(21)*3#*PRES*PRES 37600 BM0V=BM0V+(ZPIT(24)+ZPIT(25)*2#*PRES+ZPIT(26)*3#*PRES*PRES)*T 37650 BM0V=BM0V+(ZPIT(28)+ZPIT(29)*2#*PRES)*T*T 37700 BM0V=BM0V+(ZPIT(31)+ZPIT(32)*2#*PRES+ZPIT(33)*3#*PRES*PRES)/TL 37750 BM0V=BM0V+(ZPIT(35)+ZPIT(36)*2#*PRES+ZPIT(37)*3#*PRES*PRES)/TH 37800 CMV=ZPIT(44) 37850 CMV=CMV+ZPIT(47)*T+ZPIT(49)*T*T 37900 CMV=CMV+ZPIT(51)/TL+ZPIT(53)/TH 37950 CMV=CMV/2# 38000 VNACL=RVGAS*(ZPIT(2)+2#*ZPIT(3)*PRES+3#*ZPIT(4)*PRES*PRES) 38050 VNACL=VNACL+RVGAS*T*(ZPIT(6)+2#*ZPIT(7)*PRES+3#*ZPIT(8)*PRES*PRES) 38100 VNACL=VNACL+RVGAS*T*T*(ZPIT(11)+2#*ZPIT(12)*PRES)+RVGAS*T*T*T*ZPIT(14) 38150 VNACL=VNACL-YPPB*MW/D-2#*AV*HFREF-2#*RVGAS*T*MREF*BM0V-2#*RVGAS*T*MREF*MREF*CMV 38200 LPRINT 38250 LPRINT USING"V(water)=+##.### VNaCl=+###.##";MW/D,VNACL 38300 LPRINT 38350 BM0=ZPIT(17)/T+ZPIT(18)+ZPIT(19)*PRES+ZPIT(20)*PRES*PRES+ZPIT(21)*PRES*PRES*PRES+ZPIT(22)*LOG(T) 38400 BM0=BM0+(ZPIT(23)+ZPIT(24)*PRES+ZPIT(25)*PRES*PRES+ZPIT(26)*PRES*PRES*PRES)*T 38450 BM0=BM0+(ZPIT(27)+ZPIT(28)*PRES+ZPIT(29)*PRES*PRES)*T*T 38500 BM0=BM0+(ZPIT(30)+ZPIT(31)*PRES+ZPIT(32)*PRES*PRES+ZPIT(33)*PRES*PRES*PRES)/TL 38550 BM0=BM0+(ZPIT(34)+ZPIT(35)*PRES+ZPIT(36)*PRES*PRES+ZPIT(37)*PRES*PRES*PRES)/TH 38600 BM1=ZPIT(38)/T+ZPIT(39)+ZPIT(40)*T+ZPIT(41)/TL 38650 CM=ZPIT(42)/T+ZPIT(43)+ZPIT(44)*PRES+ZPIT(45)*LOG(T) 38700 CM=CM+(ZPIT(46)+ZPIT(47)*PRES)*T+(ZPIT(48)+ZPIT(49)*PRES)*T*T 38750 CM=CM+(ZPIT(50)+ZPIT(51)*PRES)/TL+(ZPIT(52)+ZPIT(53)*PRES)/TH 38800 BM0L=(-1#)*ZPIT(17)/(T*T)+ZPIT(22)/T 38850 BM0L=BM0L+(ZPIT(23)+ZPIT(24)*PRES+ZPIT(25)*PRES*PRES+ZPIT(26)*PRES*PRES*PRES) 38900 BM0L=BM0L+(ZPIT(27)+ZPIT(28)*PRES+ZPIT(29)*PRES*PRES)*2#*T 38950 BM0L=BM0L-(ZPIT(30)+ZPIT(31)*PRES+ZPIT(32)*PRES*PRES+ZPIT(33)*PRES*PRES*PRES)/(TL*TL) 39000 BM0L=BM0L+(ZPIT(34)+ZPIT(35)*PRES+ZPIT(36)*PRES*PRES+ZPIT(37)*PRES*PRES*PRES)/(TH*TH) 39050 BM1L=(-1#)*ZPIT(38)/(T*T)+ZPIT(40)-ZPIT(41)/(TL*TL) 39100 CML=(-1#)*ZPIT(42)/(T*T)+ZPIT(45)/T 39150 CML=CML+(ZPIT(46)+ZPIT(47)*PRES)+(ZPIT(48)+ZPIT(49)*PRES)*2#*T 39200 CML=CML-(ZPIT(50)+ZPIT(51)*PRES)/(TL*TL)+(ZPIT(52)+ZPIT(53)*PRES)/(TH*TH) 39250 CML=CML/2# 39300 BM0J=ZPIT(22)/(T*T)+2#*(ZPIT(23)+ZPIT(24)*PRES+ZPIT(25)*PRES*PRES+ZPIT(26)*PRES*PRES*PRES)/T 39350 BM0J=BM0J+6#*(ZPIT(27)+ZPIT(28)*PRES+ZPIT(29)*PRES*PRES) 39400 BM0J=BM0J+454#*(ZPIT(30)+ZPIT(31)*PRES+ZPIT(32)*PRES*PRES+ZPIT(33)*PRES*PRES*PRES)/(TL*TL*TL*T) 39450 BM0J=BM0J+1360*(ZPIT(34)+ZPIT(35)*PRES+ZPIT(36)*PRES*PRES+ZPIT(37)*PRES*PRES*PRES)/(T*TH*TH*TH) 39500 BM1J=2#*ZPIT(40)/T+454#*ZPIT(41)/(T*TL*TL*TL) 39550 CMJ=ZPIT(45)/(2#*T*T)+(ZPIT(46)+ZPIT(47)*PRES)/T+3#*(ZPIT(48)+ZPIT(49)*PRES) 39600 CMJ=CMJ+227#*(ZPIT(50)+ZPIT(51)*PRES)/(T*TL*TL*TL) 39650 CMJ=CMJ+680*(ZPIT(52)+ZPIT(53)*PRES)/(T*TH*TH*TH) 39700 EAIREF=EXP((-2#)*SQR(MIREF)) 39750 GEXREF=(-8#)*APHI*HFREF+2#*MREF*BM0+(4#*BM1/4#)*(1#-(1#+2#*SQR(MIREF))*EAIREF)+MREF*MREF*CM 39800 GNACL=(-1#)*YPPB*G/(RGAS*T)-GEXREF+(ZPIT(1)+ZPIT(2)*PRES+ZPIT(3)*PRES*PRES+ZPIT(4)*PRES*PRES*PRES)/T 39850 GNACL=GNACL+ZPIT(5)+ZPIT(6)*PRES+ZPIT(7)*PRES*PRES 39900 GNACL=GNACL+ZPIT(8)*PRES*PRES*PRES+ZPIT(9)*LOG(T) 39950 GNACL=GNACL+(ZPIT(10)+ZPIT(11)*PRES+ZPIT(12)*PRES*PRES)*T+(ZPIT(13)+ZPIT(14)*PRES)*T*T 40000 GNACL=GNACL+ZPIT(15)/(T*TL)+ZPIT(16)/(T*TH*TH*TH) 40050 GNACL=GNACL*RGAS*T 40100 BLREF=BM0L+2#*BM1L*(1#-(1#+2#*SQR(MIREF))*EAIREF)/(4#*MIREF) 40150 PHILREF=2#*AH*HFREF-2#*RGAS*T*T*(MREF*BLREF+MREF*MREF*CML) 40200 HNACL=(-1#)*YPPB*H/(RGAS*T)-PHILREF/(RGAS*T) 40250 HNACL=HNACL+(ZPIT(1)+ZPIT(2)*PRES+ZPIT(3)*PRES*PRES+ZPIT(4)*PRES*PRES*PRES)/T-ZPIT(9) 40300 HNACL=HNACL-(ZPIT(10)+ZPIT(11)*PRES+ZPIT(12)*PRES*PRES)*T 40350 HNACL=HNACL-2#*(ZPIT(13)+ZPIT(14)*PRES)*T*T 40400 HNACL=HNACL+ZPIT(15)*(2#*T-227#)/(T*TL*TL) 40450 HNACL=HNACL+ZPIT(16)*(680-4#*T)/(T*TH*TH*TH*TH) 40500 HNACL=HNACL*RGAS*T 40550 SNACL=(HNACL-GNACL)/T 40600 BJREF=BM0J+2#*BM1J*(1#-(1#+2#*SQR(MIREF))*EAIREF)/(4#*MIREF) 40650 DPHILREFDT=2#*HFREF*AJ-2#*MREF*RGAS*T*T*BJREF-2#*MREF*MREF*RGAS*T*T*CMJ 40700 CPNACL=(-1#)*YPPB*CPW-DPHILREFDT-RGAS*ZPIT(9)-2#*RGAS*T*(ZPIT(10)+ZPIT(11)*PRES+ZPIT(12)*PRES*PRES) 40750 CPNACL=CPNACL-6#*RGAS*T*T*(ZPIT(13)+ZPIT(14)*PRES)-2#*RGAS*T*ZPIT(15)/(TL*TL*TL) 40800 CPNACL=CPNACL-12#*RGAS*T*ZPIT(16)/(TH*TH*TH*TH*TH) 40850 LPRINT USING"G/RT= +##.#### GNaCl/RT= +##.###";G/(RGAS*T),GNACL/(RGAS*T) 40900 LPRINT 40950 LPRINT USING"H/RT= +##.#### HNaCl/RT= +##.###";H/(RGAS*T),HNACL/(RGAS*T) 41000 LPRINT 41050 LPRINT USING"S/R= +##.#### SNaCl/R= +##.###";S/RGAS,SNACL/RGAS 41100 LPRINT 41150 LPRINT USING"Cp/R= +##.### CpNaCl/R= +###.##";CPW/RGAS,CPNACL/RGAS 41200 LPRINT 41250 MI=MOL 41300 HF=LOG(1#+1.2#*SQR(MI))/(2#*1.2#) 41350 VPHI=VNACL+2#*AV*HF+2#*RVGAS*T*(MOL*BM0V+MOL*MOL*CMV) 41400 VSOLN=VPHI*MOL+(1000/MW)*(MW/D) 41450 DSOLN=(1000+MOL*NACL)/VSOLN 41500 OSC=0 41550 GM=0 41600 OSC=1#-APHI*SQR(MI)/(1#+1.2#*SQR(MI)) 41650 EAI=EXP((-2#)*SQR(MI)) 41700 OSC=OSC+MOL*(BM0+BM1*EAI)+MOL*MOL*CM 41750 GM=(-1#)*APHI*(SQR(MI)/(1#+1.2#*SQR(MI))+4#*HF) 41800 GM1=1#-(1#+2#*SQR(MI)-4#*MI/2#)*EAI 41850 GM=GM+MOL*(2#*BM0+(2#*BM1/(4#*MI))*GM1)+(3#*MOL*MOL/2#)*CM 41900 GM=EXP(GM) 41950 BL=BM0L+2#*BM1L*(1#-(1#+2#*SQR(MI))*EAI)/(4#*MI) 42000 PHIL=2#*AH*HF-2#*RGAS*T*T*(MOL*BL+MOL*MOL*CML) 42050 BJ=BM0J+2#*BM1J*(1#-(1#+2#*SQR(MI))*EAI)/(4#*MI) 42100 CPX=2#*HF*AJ-2#*MOL*RGAS*T*T*BJ-2#*MOL*MOL*RGAS*T*T*CMJ 42150 BMX=BM0+2#*BM1*(1#-(1#+2#*SQR(MI))*EAI)/(4#*MI) 42200 GEX=(-8#)*APHI*HF*MI+2#*MOL*MOL*BMX+2#*MOL*MOL*MOL*(CM/2#) 42250 GEX=GEX*RGAS*T 42300 SX=(PHIL-GEX/MOL)/T 42350 SSPEC=S*(1000/MW)+MOL*(SNACL+SX)+2#*MOL*RGAS*(1#-LOG(MOL)) 42400 SSPEC=SSPEC/(1000+MOL*NACL) 42450 HSPEC=H*(1000/MW)+MOL*(HNACL+PHIL) 42500 HSPEC=HSPEC/(1000+MOL*NACL) 42550 CPSPEC=CPW*(1000/MW)+MOL*(CPNACL+CPX) 42600 CPSPEC=CPSPEC/(1000+MOL*NACL) 42650 LPRINT USING"m=#.##### Density(g/cm3)= +#.#####";MOL,DSOLN 42700 LPRINT USING" Osmotic coeff= +#.###";OSC 42750 LPRINT USING" Activity coeff= +#.###";GM 42800 LPRINT USING" phiL/RT= +##.###";PHIL/(RGAS*T) 42850 LPRINT USING" Ex entr/R= +##.###";SX/RGAS 42900 LPRINT USING" phiCp/R= +###.##";(CPNACL+CPX)/RGAS 42950 LPRINT USING" Hspecific(J/g)= +#.####^^^^^";HSPEC 43000 LPRINT USING" Sspecific(J/g K)= +#.###";SSPEC 43050 LPRINT USING" Cpspecific(J/g K)= +#.###";CPSPEC 43100 LPRINT 43150 RETURN