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
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
|
10 ' **************************************************************************
11 ' * VLMU.BAS *
12 ' * PENG-ROBINSON EQUATION OF STATE PHASE EQUILIBRIUM PROGRAM FOR MIXTURES *
13 ' * FOR USE IN CHAPTER 8 OF THE SECOND EDITION OF *
14 ' * CHEMICAL AND ENGINEERING THERMODYNAMICS BY S. I. SANDLER *
15 ' * PUBLISHED BY J. WILEY & SONS, INC., N. Y. *
16 ' * PROGRAM LAST REVISED 12/6/87 *
17 ' * *
18 ' **************************************************************************
19 KEY OFF: CLS
20 PRINT TAB(10);"***************************************************************"
21 PRINT TAB(10);"* PROGRAM FOR CALCULATION OF PHASE EQUILIBRIUM IN MIXTURES *"
22 PRINT TAB(10);"* USING THE PENG-ROBINSON EQUATION OF STATE *"
23 PRINT TAB(10);"* BY STANLEY I. SANDLER, DEPT. OF CHEMICAL ENGINEERING *"
24 PRINT TAB(10);"* UNIVERSITY OF DELAWARE, NEWARK, DE 19716 U. S. A. *"
25 PRINT TAB(10);"* COPYRIGHT SEPTEMBER, 1987 BY STANLEY I. SANDLER *"
26 PRINT TAB(10);"* *"
27 PRINT TAB(10);"* DUPLICATION OR DISTRIBUTION IN ANY FORM WITHOUT THE *"
28 PRINT TAB(10);"* WRITTEN PERMISSION OF AUTHOR AND COPYRIGHT HOLDER IS *"
29 PRINT TAB(10);"* PROHIBITED AND IN VIOLATION OF FEDERAL COPYRIGHT LAWS. *"
30 PRINT TAB(10);"***************************************************************"
31 L42$="####.##": L14$="#.####": L14E$="#.####^^^^": L24$="##.####"
32 L1$="-------------------------------------------------------------------"
33 L2$= "\ \ #.#### #.####^^^^ #.####^^^^"
34 L3$= "\ \ #.#### #.#### #.####^^^^"
35 L5$= "#.#### #.####"
36 L6$= "\ \ #.#### #.#### #.#### #.####^^^^"
37 L4$= "Phase compressibilities "
38 PI#=3.1415926541#
40 DIM KK1$(10),TC#(10),PC#(10),WC#(10),TB#(10)
45 DIM DEL#(10,10),X1#(10),X2#(10),Z#(10),ZX#(10),AC#(10)
50 DIM Y#(10),Y1#(10),Y2#(10),YY#(10),YZ#(10),ZZZ#(10),XX#(10),A#(10)
55 DIM C#(3),SA#(10),S#(2),D1#(10),B#(10),AIJ#(10,10)
60 DIM XLF#(10),XK1#(10),RAT1#(10),PS#(10)
70 DIM FOX#(10),F#(10),F1#(10),F2#(10)
80 R#=8.314# 'R = J/Mol K
90 RG#=.00008314# 'R = bar m^3/mol K
100 PRINT: PRINT: INPUT"Enter input device number, 0-keyboard, 1-diskfile";INPF
110 IF INPF=1 GOTO 330
120 PRINT "Enter disk filename to save critical properties"
130 INPUT "Example, a:co2&nc4.dat";FILEOUT$
140 PRINT "This file can be used for data input on rerunning program"
150 OPEN FILEOUT$ FOR OUTPUT AS 2
160 INPUT "Enter the number of components";NC
170 FOR I=1 TO NC
180 INPUT "Enter species name (alphabetic, 7 letters)";KK1$(I)
190 INPUT "Enter Tc(K),Pc(bar),wc, and Tboil(K)";TC#(I),PC#(I),WC#(I),TB#(I)
200 NEXT I
210 PRINT #2, NC
220 FOR I=1 TO NC
230 PRINT #2,KK1$(I)
240 PRINT #2,TC#(I),PC#(I),WC#(I),TB#(I)
250 NEXT I
260 CLOSE #2
270 FOR I=1 TO NC
280 FOR J=1 TO NC
290 DEL#(I,J)=0#
300 NEXT J
310 NEXT I
320 GOTO 410
330 INPUT "Enter input diskfile, ex. a:co2&nc4.dat";FILEIN$
340 OPEN FILEIN$ FOR INPUT AS 1
350 INPUT #1,NC
360 FOR I=1 TO NC
370 INPUT #1,KK1$(I)
380 INPUT #1,TC#(I),PC#(I),WC#(I),TB#(I)
390 NEXT I
400 CLOSE #1
410 IF NC=1 GOTO 600
420 PRINT "Type 0 for all binary interaction parameters equal to zero"
430 INPUT " 1 otherwise "; INTFLG
440 IF INTFLG=0 GOTO 530
450 FOR I=1 TO NC
460 FOR J=I TO NC
470 IF I=J GOTO 510
480 PRINT "k(";I;",";J;") = "
490 INPUT DEL#(I,J)
500 DEL#(J,I)=DEL#(I,J)
510 NEXT J
520 NEXT I
530 IF INTFLG=0 GOTO 600
540 PRINT "matrix of binary interaction parameters"
550 FOR I=1 TO NC
560 FOR J=I TO NC
570 PRINT DEL#(I,J)
580 NEXT J
590 NEXT I
600 PRINT:PRINT "Type 1 for compressibility calculation"
610 PRINT " 2 for species fugacity calculation"
620 PRINT " 3 for bubble point temperature calculation"
630 PRINT " 4 for dew point temperature calculation"
640 PRINT " 5 for bubble point pressure calculation"
650 PRINT " 6 for dew point pressure calculation"
660 PRINT " 7 for isothermal flash calculation"
670 PRINT " 0 to leave"
680 INPUT ICODE
690 IF ICODE=0 THEN END
700 PRINT
710 ON ICODE GOSUB 1950,3300,3630,4960,6250,7600,2260
720 GOTO 600
730 '=======================================================================
740 'subroutine prfuga
750 C1#=SQR(2#)
760 C2#=1#+C1#
770 C3#=C1#-1#
780 YMAX#=0#
790 YMIN#=0#
800 IF NC>1 GOTO 840
810 AA#=A#(1)
820 BB#=B#(1)
830 GOTO 1000
840 FOR I=1 TO NC
850 SA#(I)=0#
860 NEXT I
870 AA#=0#
880 BB#=0#
890 FOR I=1 TO NC
900 BB#=BB#+ZX#(I)*B#(I)
910 FOR M=1 TO NC
920 IF I=M GOTO 960
930 AA#=AA#+ZX#(I)*ZX#(M)*AIJ#(I,M)
940 SA#(M)=SA#(M)+ZX#(I)*AIJ#(I,M)
950 GOTO 980
960 AA#=AA#+ZX#(I)*ZX#(I)*A#(I)
970 SA#(M)=SA#(M)+ZX#(M)*A#(M)
980 NEXT M
990 NEXT I
1000 CA#=AA#*P#/((RG#*T#)^2)
1010 CB#=BB#*P#/(RG#*T#)
1020 C#(1)=CB#-1#
1030 C#(2)=CA#-CB#*(3#*CB#+2#)
1040 C#(3)=CB#*(CB#*CB#+CB#-CA#)
1050 SSA# = (3#*C#(2)-C#(1)*C#(1))/3#
1060 SB# = (2#*(C#(1)^3)-9#*C#(1)*C#(2)+27#*C#(3))/27#
1070 SR# = ((SSA#/3#)^3)+((SB#/2#)^2)
1080 IF SR#<0# THEN GOTO 1310
1090 IF SR#>0# THEN GOTO 1150
1100 AZ#=(-SB#/2#)^(1#/3#)
1110 YZ#(1)=AZ#+AZ#
1120 YZ#(2)=-YZ#(1)/2#
1130 YZ#(3)=YZ#(2)
1140 GOTO 1400
1150 XA#=(SQR(SR#))-(SB#/2#)
1160 IF XA#<0# GOTO 1170 ELSE 1200
1170 XA#=-XA#
1180 AZ#=-1#*((XA#)^(1#/3#))
1190 GOTO 1210
1200 AZ#=(XA#)^(1#/3#)
1210 XB#=-(SB#/2#)-SQR(SR#)
1220 IF XB#>0 GOTO 1260
1230 XB#=-XB#
1240 AB#=-(XB#^(1#/3#))
1250 GOTO 1270
1260 AB#=(XB#)^(1#/3#)
1270 YZ#(1)=AZ#+AB#
1280 YZ#(2)=YZ#(1)
1290 YZ#(3)=YZ#(1)
1300 GOTO 1400
1310 TEMP#=((SSA#/3#)^3#)
1320 TEMP#=SQR(-TEMP#)
1330 CSPHI#=-SB#/(2#*TEMP#)
1340 TANPHI#=(SQR(1#-CSPHI#^2))/CSPHI#
1350 PHI#=ATN(TANPHI#)
1360 T11# = 2#*SQR(-SSA#/3#)
1370 YZ#(1) = T11#*COS(PHI#/3#)
1380 YZ#(2) = T11#*COS((PHI# + 2#*PI#)/3#)
1390 YZ#(3) = T11#*COS((PHI# + 4#*PI#)/3#)
1400 FOR I = 1 TO 3
1410 YZ#(I) = YZ#(I)-C#(1)/3#
1420 NEXT I
1430 IF YZ#(1)>YZ#(2) THEN 1460
1440 IF YZ#(2)>YZ#(3) THEN YMAX#=YZ#(2) ELSE YMAX#=YZ#(3)
1450 GOTO 1470
1460 IF YZ#(1)>YZ#(3) THEN YMAX#=YZ#(1) ELSE YMAX#=YZ#(3)
1470 IF YZ#(1)<YZ#(2) THEN 1500
1480 IF YZ#(2)<YZ#(3) THEN YMIN#=YZ#(2) ELSE YMIN#=YZ#(3)
1490 GOTO 1510
1500 IF YZ#(1)<YZ#(3) THEN YMIN#=YZ#(1) ELSE YMIN#=YZ#(3)
1510 Z0#=YMAX#
1520 Z1#=YMIN#
1530 V1# = Z1#*RG#*T#/P#
1540 IF V1#<BB# THEN Z1#=Z0#
1550 ' 0 = VAPOR AND 1 = LIQUID
1580 IF ICOMP=0 THEN RETURN
1590 IF IFASE = 0 THEN ZZ# = Z0# ELSE ZZ# = Z1#
1600 AG1#=(ZZ#+C2#*CB#)/(ZZ#-C3#*CB#)
1610 AG1#=LOG(AG1#)
1620 AG2#=CA#/(2#*CB#*C1#)
1630 IF NC>1 GOTO 1680
1640 FOX#(1)=ZZ#-1#-LOG(ZZ#-CB#)-AG1#*AG2#
1650 FOX#(1)=EXP(FOX#(1))
1660 F#(1)=P#*FOX#(1)
1670 GOTO 1740
1680 FOR M=1 TO NC
1690 AG3#=(2#*SA#(M)/AA#)-(B#(M)/BB#)
1700 FOX#(M)=(B#(M)*(ZZ#-1#)/BB#)-LOG(ZZ#-CB#)-AG1#*AG2#*AG3#
1710 FOX#(M)=EXP(FOX#(M))
1720 F#(M)=ZX#(M)*P#*FOX#(M)
1730 NEXT M
1740 RETURN
1750 '========================================================================
1760 'subroutine prcons
1770 FOR I=1 TO NC
1780 AC#(I)=.457235529#*((RG#*TC#(I))^2)/PC#(I)
1790 B#(I)=.077796074#*RG#*TC#(I)/PC#(I)
1800 XK#=.37464#+(1.54226#-.26992#*WC#(I))*WC#(I)
1810 TR#=T#/TC#(I)
1820 ALSQR#=1#+XK#*(1#-SQR(TR#))
1830 ALPHA#=ALSQR#*ALSQR#
1840 A#(I)=ALPHA#*AC#(I)
1850 NEXT I
1860 IF NC=1 GOTO 1930
1870 FOR I=1 TO NC-1
1880 FOR J=I+1 TO NC
1890 AIJ#(I,J)=(1#-DEL#(I,J))*SQR(A#(I)*A#(J))
1900 AIJ#(J,I)=AIJ#(I,J)
1910 NEXT J
1920 NEXT I
1930 RETURN
1940 '=======================================================================
1950 'subroutine comp
1960 SUM#=0#
1970 INPUT "Enter feed temperature (K), pressure (bar)";T#,P#
1980 IF NC=1 THEN Z#(1)=1#
1990 IF NC=1 GOTO 2080
2000 FOR I=1 TO NC
2010 PRINT "Enter feed mole fraction of ";KK1$(I);
2020 INPUT Z#(I)
2030 SUM#=SUM#+Z#(I)
2040 NEXT I
2050 FOR I=1 TO NC
2060 Z#(I)=Z#(I)/SUM#
2070 NEXT I
2080 GOSUB 1760 'call prcons
2090 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
2100 IFASE=0:ICOMP=0
2110 GOSUB 740 'call prfuga
2120 PRINT:PRINT L1$
2130 PRINT "Compressibility calculation using the Peng-Robinson EOS": PRINT
2140 PRINT "Temperature= ";:PRINT USING L42$;T#;:PRINT " K"
2150 PRINT "Pressure = ";:PRINT USING L42$;P#;:PRINT " bar": PRINT
2160 PRINT,:PRINT"Compressibility",:PRINT"Specific volume (m^3/mol)"
2170 PRINT "Vapor",:PRINT USING L14$;Z0#;:PRINT,,:PRINT USING L14E$;Z0#*RG#*T#/P#
2180 PRINT "Liquid",:PRINT USING L14$;Z1#;:PRINT,,:PRINT USING L14E$;Z1#*RG#*T#/P#
2190 PRINT L1$
2200 PRINT
2210 PRINT "Type 1 to run again with same components and kijs"
2220 PRINT " 0 to return to main menu";:INPUT IENDFG
2230 IF IENDFG<>0 GOTO 1960
2240 RETURN
2250 '=======================================================================
2260 'subroutine flash
2270 ICOMP = 1
2280 IF NC<=1 THEN PRINT"Need 2 or more components to do an isothermal flash calculation"
2290 IF NC<=1 GOTO 3280
2300 EPS2#=.0001#
2310 KVALUE=0
2320 SUM#=0#:ICONV=0
2330 INPUT "Enter feed temperature (K) ";T#
2340 INPUT "Enter feed pressure (bar) ";P#
2350 FOR I=1 TO NC
2360 PRINT "Enter feed mole fraction of ";KK1$(I);
2370 INPUT Z#(I)
2380 SUM#=SUM#+Z#(I)
2390 NEXT I
2400 FOR I=1 TO NC
2410 Z#(I)=Z#(I)/SUM#
2420 NEXT I
2430 GOSUB 1760 'call prcons
2440 IF KVALUE=1 GOTO 2570
2450 FOR I=1 TO NC
2460 IF T#>(1.1#*TC#(I)) GOTO 2520
2470 DT1#=(1#/T#)-(1#/TB#(I))
2480 DT2#=(1#/TC#(I))-(1#/TB#(I))
2490 DLNP#=LOG(PC#(I))
2500 PS#(I)=EXP(DLNP#*DT1#/DT2#)
2510 GOTO 2550
2520 TR#=T#/TC#(I)
2530 AXP#=7.224#-2.598#*LOG(TR#)-7.534#/TR#
2540 PS#(I)=PC#(I)*EXP(AXP#)
2550 XK1#(I)=PS#(I)/P#
2560 NEXT I
2570 XLF#(1)=1#/2#
2580 K=0
2590 NLAP=0
2600 FOR I=1 TO NC
2610 D1#(I)=1#-XK1#(I)
2620 NEXT I
2630 NLAP=NLAP+1
2640 SXA#=0#
2650 SXV#=0#
2660 S1A#=0#
2670 FOR I=1 TO NC
2680 DENO#=XLF#(1)*D1#(I)+XK1#(I)
2690 X1#(I)=Z#(I)/DENO#
2700 Y#(I)=X1#(I)*XK1#(I)
2710 SXA#=SXA#+X1#(I)
2720 SXV#=SXV#+Y#(I)
2730 DS1A#=-X1#(I)*D1#(I)*D1#(I)/DENO#
2740 S1A#=S1A#+DS1A#
2750 NEXT I
2760 IF NLAP>20 GOTO 2840
2770 XMX#=SXA#-SXV#
2780 IF ABS(XMX#)<EPS2# GOTO 2840
2790 IF S1A#=0# GOTO 3070
2800 DA#=XMX#/S1A#
2810 XLF#(1)=XLF#(1)-DA#
2820 XLF#(2)=1#-XLF#(1)
2830 GOTO 2630
2840 K=K+1: IF K>75 GOTO 3210
2850 FOR I=1 TO NC:ZX#(I)=X1#(I):NEXT I
2860 IFASE=1
2870 GOSUB 740 'call prfuga
2880 FOR I=1 TO NC:F1#(I)=F#(I):NEXT I
2890 ZZ1#=ZZ#
2900 FOR I=1 TO NC:ZX#(I)=Y#(I):NEXT I
2910 IFASE=0
2920 GOSUB 740 'call prfuga
2930 FOR I=1 TO NC:F2#(I)=F#(I):NEXT I
2940 ZZ2#=ZZ#
2950 IFLAG=0
2970 FOR I=1 TO NC
2980 RAT1#(I)=F1#(I)/F2#(I)
2990 ZZZ#=ABS(RAT1#(I)-1#)
3000 IF ZZZ#<.0001# GOTO 3030
3010 XK1#(I)=XK1#(I)*RAT1#(I)
3020 IFLAG=IFLAG+1
3030 NEXT I
3040 IF IFLAG=0 GOTO 3060
3050 GOTO 2590
3060 IF (XLF#(1)<0#) OR (XLF#(1)>1#) GOTO 3210
3070 PRINT: PRINT L1$
3080 PRINT "Isothermal flash calculation using the PENG-ROBINSON equation of state"
3090 PRINT "Temperature = ";:PRINT USING L42$;T#;:PRINT" K"
3100 PRINT "Pressure = ";:PRINT USING L42$;P#;:PRINT" bar":PRINT
3110 PRINT "Component Feed Liquid Vapor K"
3120 FOR I=1 TO NC
3130 PRINT USING L6$; KK1$(I), Z#(I), X1#(I), Y#(I), XK1#(I)
3140 NEXT I
3150 PRINT
3160 PRINT "Split",:PRINT,:PRINT USING L14$;XLF#(1);:PRINT,:PRINT USING L14$;XLF#(2)
3170 PRINT:PRINT L4$;:PRINT,:PRINT USING L14$;ZZ1#;:PRINT,:PRINT USING L14$;ZZ2#
3180 PRINT:PRINT" Number of iterations = " ;K
3190 PRINT L1$: PRINT: ICONV=1
3200 GOTO 3220
3210 PRINT "No flash calculation possibile at these conditions"
3220 PRINT "Type 1 to run again with same components and kijs"
3230 PRINT " 0 to return to main menu";:INPUT IENDFG
3240 IF IENDFG=0 GOTO 3280
3245 IF ICONV=0 GOTO 2320
3250 PRINT "Type 1 to use last K values as initial guess"
3260 INPUT " 0 to reinitialize K values ";KVALUE
3270 GOTO 2320
3280 RETURN
3290 '========================================================================
3300 'subroutine spfug
3310 ICOMP=1
3320 SUM#=0#
3330 INPUT "Enter feed temperature (K), pressure (bar)";T#,P#
3340 IF NC=1 THEN Z#(1)=1#
3350 IF NC=1 GOTO 3440
3360 FOR I=1 TO NC
3370 PRINT "Enter feed mole fraction of ";KK1$(I);
3380 INPUT Z#(I)
3390 SUM#=SUM#+Z#(I)
3400 NEXT I
3410 FOR I=1 TO NC
3420 Z#(I)=Z#(I)/SUM#
3430 NEXT I
3440 INPUT "Enter 0 for vapor fugacity, 1 for liquid";IFASE
3450 GOSUB 1760 'call prcons
3460 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
3470 GOSUB 740 'call prfuga
3480 PRINT:PRINT L1$
3490 PRINT "Fugacity calculation using the Peng-Robinson equation of state"
3500 PRINT "Temperature = ";:PRINT USING L42$;T#;:PRINT " K"
3510 PRINT "Pressure = ";:PRINT USING L42$;P#;:PRINT" bar":PRINT
3520 PRINT "Component Feed Fugacity (bar) (f/xP) "
3530 FOR I=1 TO NC
3540 PRINT USING L2$;KK1$(I), ZX#(I), F#(I), FOX#(I)
3550 NEXT I
3560 PRINT:PRINT L4$;:PRINT USING L14$;ZZ#
3570 PRINT L1$: PRINT
3580 PRINT "Type 1 to run again with same components and kijs"
3590 PRINT " 0 to return to main menu";:INPUT IENDFG
3600 IF IENDFG<>0 GOTO 3320
3610 RETURN
3620 '=======================================================================
3630 'subroutine bubptt
3640 K=0: ICOMP=1: KVALUE=0
3650 SUM#=0#: ICONV=0
3660 INPUT "Enter feed pressure (bar)";P#
3670 PRINT "Enter estimated bubble point temperature (K)"
3680 INPUT " or 0 for default estimate";T#
3690 IF NC=1 THEN Z#(1)=1#
3700 IF NC=1 GOTO 3790
3710 FOR I=1 TO NC
3720 PRINT "Enter feed mole fraction of ";KK1$(I);
3730 INPUT Z#(I)
3740 SUM#=SUM#+Z#(I)
3750 NEXT I
3760 FOR I=1 TO NC
3770 Z#(I)=Z#(I)/SUM#
3780 NEXT I
3790 IF T#<>0 GOTO 3870
3800 TBG#=0
3810 TCG#=0
3820 FOR I=1 TO NC
3830 TBG#=TBG#+Z#(I)*TB#(I)
3840 TCG#=TCG#+Z#(I)*TC#(I)
3850 NEXT I
3860 T#=.5#*(TBG#+TCG#)
3870 IF K>30 GOTO 4880
3880 FOR I=1 TO NC
3890 IF KVALUE = 1 GOTO 3950
3900 DT1#=(1#/T#)-(1#/TB#(I))
3910 DT2#=(1#/TC#(I))-(1#/TB#(I))
3920 DLNP#=LOG(PC#(I))
3930 PS#(I)=EXP(DLNP#*DT1#/DT2#)
3940 XK1#(I)=PS#(I)/P#
3950 YY#(I)=Z#(I)*XK1#(I)
3960 NEXT I
3970 KKK=0
3980 NLOOP=1
3990 IF T#<50# GOTO 4590
4000 IF T#>1200# GOTO 4590
4010 K=K+1
4020 GOSUB 1760 'call prcons
4030 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
4040 IFASE=1
4050 GOSUB 740 'call prfuga
4060 FOR I=1 TO NC:F1#(I)=F#(I):NEXT I
4070 ZZ1#=ZZ#
4080 SUMY#=0#
4090 FOR I=1 TO NC
4100 Y2#(I)=YY#(I)
4110 SUMY#=SUMY#+YY#(I)
4120 NEXT I
4130 FOR I=1 TO NC
4140 YY#(I)=YY#(I)/SUMY#
4150 NEXT I
4160 SUMY#=0#
4170 KKK=KKK+1
4180 FOR I=1 TO NC:ZX#(I)=YY#(I):NEXT I
4190 IFASE=0
4200 GOSUB 740 'call prfuga
4210 FOR I=1 TO NC:F2#(I)=F#(I):NEXT I
4220 ZZ2#=ZZ#
4230 FOR I=1 TO NC
4240 YY#(I)=F1#(I)*YY#(I)/F2#(I)
4250 Y1#(I)=YY#(I)
4260 SUMY#=SUMY#+YY#(I)
4270 NEXT I
4280 ITEST=0
4290 FOR I=1 TO NC
4300 TEST#=ABS(Y1#(I)-Y2#(I))
4310 IF TEST#>.0001# THEN ITEST=ITEST+1
4320 YY#(I)=YY#(I)/SUMY#
4330 NEXT I
4340 IF KKK>25 GOTO 4400
4350 IF ITEST<=0 GOTO 4400
4360 FOR I=1 TO NC
4370 Y2#(I)=Y1#(I)
4380 NEXT I
4390 GOTO 4160
4400 S#(NLOOP)=SUMY#
4410 KKK=0
4420 IF (NLOOP-1)>0 GOTO 4460
4430 NLOOP=2
4440 T#=T#-.005#
4450 GOTO 4020
4460 DSDT#=(S#(2)-S#(1))/(.005#)
4470 IF (ABS(DSDT#)<.00001#) GOTO 4590
4480 DLT#=(S#(1)-1#)/DSDT#
4490 IF ABS(DLT#)<.0026# GOTO 4670
4500 IF K>50 GOTO 4650
4510 IF K<11 THEN DD#=20#
4520 IF K>=11 THEN DD#=5#
4530 IF DLT#>DD# THEN T#=T#+DD#
4540 IF DLT#>DD# GOTO 3980
4550 IF DLT#<-DD# THEN T#=T#-DD#
4560 IF -DLT#>DD# GOTO 3980
4570 T#=T#+DLT#+.0025#
4580 GOTO 3980
4590 IF K>2 GOTO 4630
4600 IF ZZ1#>=.307# THEN T#=T#-10#
4610 IF ZZ2#<=.307# THEN T#=T#+10#
4620 GOTO 3980
4630 PRINT "Calculation not converging: one-phase region or poor initial guess"
4640 GOTO 4880
4650 PRINT:PRINT "Bubble point calculation did not converge":PRINT
4660 GOTO 4880
4670 YK#=K
4680 TTEST#=(ZZ1#-ZZ2#)^2
4690 IF TTEST#>.00001# GOTO 4730
4700 IF ZZ1#>=.307# THEN T#=T#-25#/SQR(YK#)
4710 IF ZZ1#<.307# THEN T#=T#+25#/SQR(YK#)
4720 KVALUE = 0:GOTO 3870
4730 FOR I=1 TO NC
4740 XK1#(I)=YY#(I)/Z#(I)
4750 NEXT I
4760 PRINT:PRINT L1$
4770 PRINT "Bubble point temperature calculation using the "
4780 PRINT "Peng Robinson equation of state"
4790 PRINT "Pressure = ";:PRINT USING L42$; P#;:PRINT" bar"
4800 PRINT "Calculated bubble point temperature = ";:PRINT USING L42$;T#;:PRINT" K":PRINT
4810 PRINT "Component Feed Vapor K"
4820 FOR I=1 TO NC
4830 PRINT USING L3$; KK1$(I), Z#(I), YY#(I), XK1#(I)
4840 NEXT I
4850 PRINT:PRINT L4$:PRINT,:PRINT USING L5$; ZZ1#, ZZ2#: PRINT
4860 PRINT "Number of iterations = " ;K
4870 PRINT L1$: PRINT: ICONV=1
4880 PRINT "Type 1 to run again with same components and kijs"
4890 PRINT " 0 to return to main menu";:INPUT IENDFG
4900 IF IENDFG=0 GOTO 4940
4905 IF ICONV=0 GOTO 3640
4910 PRINT "Type 1 to use last K values as initial guess"
4920 INPUT " 0 to reinitialize K values ";KVALUE
4930 K = 0: GOTO 3650
4940 RETURN
4950 '========================================================================
4960 'subroutine dewptt
4970 ICOMP=1: KVALUE=0
4980 K=0
4990 SUM#=0#: ICONV=0
5000 INPUT "Enter feed pressure (bar)";P#
5010 PRINT "Enter estimated dew point temperature (K)"
5020 INPUT " or 0 for default estimate";T#
5030 IF NC=1 THEN Z#(1)=1#
5040 IF NC=1 GOTO 5130
5050 FOR I=1 TO NC
5060 PRINT "Enter feed mole fraction of ";KK1$(I);
5070 INPUT Z#(I)
5080 SUM#=SUM#+Z#(I)
5090 NEXT I
5100 FOR I=1 TO NC
5110 Z#(I)=Z#(I)/SUM#
5120 NEXT I
5130 IF T#<>0 GOTO 5210
5140 TBG#=0
5150 TCG#=0
5160 FOR I=1 TO NC
5170 TBG#=TBG#+Z#(I)*TB#(I)
5180 TCG#=TCG#+Z#(I)*TC#(I)
5190 NEXT I
5200 T#=.5#*(TBG#+TCG#)
5210 IF K>30 GOTO 6180
5220 FOR I=1 TO NC
5230 IF KVALUE = 1 GOTO 5290
5240 DT1#=(1#/T#)-(1#/TB#(I))
5250 DT2#=(1#/TC#(I))-(1#/TB#(I))
5260 DLNP#=LOG(PC#(I))
5270 PS#(I)=EXP(DLNP#*DT1#/DT2#)
5280 XK1#(I)=PS#(I)/P#
5290 XX#(I)=Z#(I)/XK1#(I)
5300 NEXT I
5310 KKK=0
5320 NLOOP=1
5330 IF T#<100 GOTO 5930
5340 IF T#>1200 GOTO 5930
5350 K=K+1
5360 GOSUB 1760 'call prcons
5370 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
5380 IFASE=0
5390 GOSUB 740 'call prfuga
5400 FOR I=1 TO NC:F1#(I)=F#(I):NEXT I
5410 ZZ1#=ZZ#
5420 SUMX#=0#
5430 FOR I=1 TO NC
5440 X2#(I)=XX#(I)
5450 SUMX#=SUMX#+XX#(I)
5460 NEXT I
5470 FOR I=1 TO NC
5480 XX#(I)=XX#(I)/SUMX#
5490 NEXT I
5500 SUMX#=0#
5510 KKK=KKK+1
5520 FOR I=1 TO NC:ZX#(I)=XX#(I):NEXT I
5530 IFASE=1
5540 GOSUB 740 'call prfuga
5550 FOR I=1 TO NC:F2#(I)=F#(I):NEXT I
5560 ZZ2#=ZZ#
5570 FOR I=1 TO NC
5580 XX#(I)=F1#(I)*XX#(I)/F2#(I)
5590 X1#(I)=XX#(I)
5600 SUMX#=SUMX#+XX#(I)
5610 NEXT I
5620 ITEST=0
5630 FOR I=1 TO NC
5640 TEST#=ABS(X1#(I)-X2#(I))
5650 IF TEST#>.0001# THEN ITEST=ITEST+1
5660 XX#(I)=XX#(I)/SUMX#
5670 NEXT I
5680 IF KKK>25 GOTO 5740
5690 IF ITEST<=0 GOTO 5740
5700 FOR I=1 TO NC
5710 X2#(I)=X1#(I)
5720 NEXT I
5730 GOTO 5500
5740 S#(NLOOP)=SUMX#
5750 KKK=0
5760 IF (NLOOP-1)>0 GOTO 5800
5770 NLOOP=2
5780 T#=T#-.005#
5790 GOTO 5360
5800 DSDT#=(S#(2)-S#(1))/.005#
5810 IF (ABS(DSDT#)<.00001#) GOTO 5930
5820 DLT#=(S#(1)-1#)/DSDT#
5830 IF K>50 GOTO 5950
5840 IF K<11 THEN DD#=50#
5850 IF K>=11 THEN DD#=10#
5860 IF DLT#>DD# THEN T#=T#+DD#
5870 IF DLT#>DD# GOTO 5320
5880 IF DLT#<-DD# THEN T#=T#-DD#
5890 IF DLT#<-DD# GOTO 5320
5900 IF ABS(DLT#)<.01# GOTO 5970
5910 T#=T#+DLT#+.0025#
5920 GOTO 5320
5930 PRINT "Calculation not converging: one-phase region or poor initial guess"
5940 GOTO 6180
5950 PRINT:PRINT "Dew point calculation did not converge": PRINT
5960 GOTO 6180
5970 YK#=K
5980 TTEST#=(ZZ1#-ZZ2#)^2
5990 IF TTEST#>.00001# GOTO 6030
6000 IF ZZ1#>=.307# THEN T#=T#-25#/SQR(YK#)
6010 IF ZZ1#<.307# THEN T#=T#+25#/SQR(YK#)
6020 KVALUE = 0: GOTO 5210
6030 FOR I=1 TO NC
6040 XK1#(I)=Z#(I)/XX#(I)
6050 NEXT I
6060 PRINT:PRINT L1$
6070 PRINT "Dew point temperature calculation using the "
6080 PRINT "Peng Robinson equation of state"
6090 PRINT "Pressure = ";:PRINT USING L42$;P#;:PRINT" bar"
6100 PRINT "Calculated dew point temperature = ";:PRINT USING L42$;T#;:PRINT" K":PRINT
6110 PRINT "Component Feed Liquid K"
6120 FOR I=1 TO NC
6130 PRINT USING L3$; KK1$(I), Z#(I), XX#(I), XK1#(I)
6140 NEXT I
6150 PRINT:PRINT L4$:PRINT,:PRINT USING L5$; ZZ1#, ZZ2#: PRINT
6160 PRINT "Number of iterations = " ;K
6170 PRINT L1$: PRINT: ICONV=1
6180 PRINT "Type 1 to run again with same components and kijs"
6190 PRINT " 0 to return to main menu";:INPUT IENDFG
6200 IF IENDFG=0 THEN RETURN
6205 IF ICONV=0 GOTO 4980
6210 PRINT "Type 1 to use last K values as initial guess"
6220 INPUT " 0 to reinitialize K values ";KVALUE
6230 GOTO 4980
6240 '========================================================================
6250 'subroutine bubptp
6260 ICOMP=1
6270 IENDFG=0
6280 KVALUE = 0
6290 PBG#=0#
6300 K=0
6310 SUM#=0#: ICONV=0
6320 INPUT "Enter feed temperature (K)";T#
6330 PRINT "Enter estimated bubble point pressure (bar)"
6340 INPUT " or 0 for default estimate";P#
6350 IF NC=1 THEN Z#(1)=1#
6360 IF NC=1 GOTO 6450
6370 FOR I=1 TO NC
6380 PRINT "Enter feed mole fraction of ";KK1$(I);
6390 INPUT Z#(I)
6400 SUM#=SUM#+Z#(I)
6410 NEXT I
6420 FOR I=1 TO NC
6430 Z#(I)=Z#(I)/SUM#
6440 NEXT I
6450 GOSUB 1760 'call prcons
6460 IF KVALUE=1 AND P#>.0001 GOTO 6600
6470 FOR I=1 TO NC
6480 IF T#>(1.1#*TC#(I)) GOTO 6540
6490 DT1#=(1#/T#)-(1#/TB#(I))
6500 DT2#=(1#/TC#(I))-(1#/TB#(I))
6510 DLNP#=LOG(PC#(I))
6520 PS#(I)=EXP(DLNP#*DT1#/DT2#)
6530 GOTO 6570
6540 TR#=T#/TC#(I)
6550 AXP#=7.224#-2.598#*LOG(TR#)-7.534#/TR#
6560 PS#(I)=PC#(I)*EXP(AXP#)
6570 PBG#=PBG#+Z#(I)*PS#(I)
6580 NEXT I
6590 IF P#<.0001# THEN P#=PBG#
6600 FOR I=1 TO NC
6610 IF KVALUE = 1 GOTO 6630
6620 XK1#(I)=PS#(I)/P#
6630 YY#(I)=Z#(I)*XK1#(I)
6640 NEXT I
6650 KKK=0
6660 NLOOP=1
6670 IF P#<0# GOTO 7270
6680 IF P#>10000# GOTO 7270
6690 K=K+1
6700 IF K>150 GOTO 7290
6710 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
6720 IFASE=1
6730 GOSUB 740 'call prfuga
6740 FOR I=1 TO NC:F1#(I)=F#(I):NEXT I
6750 ZZ1#=ZZ#
6760 SUMY#=0#
6770 FOR I=1 TO NC
6780 Y2#(I)=YY#(I)
6790 SUMY#=SUMY#+YY#(I)
6800 NEXT I
6810 FOR I=1 TO NC
6820 YY#(I)=YY#(I)/SUMY#
6830 NEXT I
6840 SUMY#=0#
6850 KKK=KKK+1
6860 FOR I=1 TO NC:ZX#(I)=YY#(I):NEXT I
6870 IFASE=0
6880 GOSUB 740 'call prfuga
6890 FOR I=1 TO NC:F2#(I)=F#(I):NEXT I
6900 ZZ2#=ZZ#
6910 FOR I=1 TO NC
6920 YY#(I)=F1#(I)*YY#(I)/F2#(I)
6930 Y1#(I)=YY#(I)
6940 SUMY#=SUMY#+YY#(I)
6950 NEXT I
6960 ITEST=0
6970 FOR I=1 TO NC
6980 TEST#=ABS(Y1#(I)-Y2#(I))
6990 IF TEST#>.0001# THEN ITEST=ITEST+1
7000 YY#(I)=YY#(I)/SUMY#
7010 NEXT I
7020 IF KKK>25 GOTO 7080
7030 IF ITEST<=0 GOTO 7080
7040 FOR I=1 TO NC
7050 Y2#(I)=Y1#(I)
7060 NEXT I
7070 GOTO 6840
7080 S#(NLOOP)=SUMY#
7090 KKK=0
7100 IF (NLOOP-1)>0 GOTO 7140
7110 NLOOP=2
7120 P#=P#-.005#
7130 GOTO 6710
7140 DSDP#=(S#(2)-S#(1))/.005#
7150 IF (ABS(DSDP#)<.000001#) GOTO 7270
7160 DLP#=(S#(1)-1#)/DSDP#
7170 IF K>50 GOTO 7290
7180 IF K<11 THEN DD#=50#
7190 IF K>=11 THEN DD#=10#
7200 IF DLP#>DD# THEN P#=P#+DD#
7210 IF DLP#>DD# GOTO 6660
7220 IF DLP#<-DD# THEN P#=P#-DD#
7230 IF DLP#<-DD# GOTO 6660
7240 IF ABS(DLP#)<.01# GOTO 7310
7250 P#=P#+DLP#+.0025#
7260 GOTO 6660
7270 PRINT "Calculation not converging: one-phase region or poor initial guess"
7280 GOTO 7520
7290 PRINT:PRINT "Bubble point calculation did not converge":PRINT
7300 GOTO 7520
7310 YK#=K
7320 TTEST#=(ZZ1#-ZZ2#)^2
7330 IF TTEST#>.00001# GOTO 7370
7340 IF ZZ1#>=.307# THEN P#=P#+25#/SQR(YK#)
7350 IF ZZ1#<.307# THEN P#=P#-25#/SQR(YK#)
7360 KVALUE = 0: GOTO 6600
7370 FOR I=1 TO NC
7380 XK1#(I)=YY#(I)/Z#(I)
7390 NEXT I
7400 PRINT:PRINT L1$
7410 PRINT "Bubble point pressure calculation using the "
7420 PRINT "Peng Robinson equation of state"
7430 PRINT "Temperature = ";:PRINT USING L42$;T#;:PRINT" K"
7440 PRINT "Calculated bubble point pressure = ";:PRINT USING L42$;P#;:PRINT" bar":PRINT
7450 PRINT "Component Feed Vapor K"
7460 FOR I=1 TO NC
7470 PRINT USING L3$; KK1$(I),Z#(I), YY#(I), XK1#(I)
7480 NEXT I
7490 PRINT:PRINT L4$:PRINT,:PRINT USING L5$; ZZ1#, ZZ2#: PRINT
7500 PRINT "Number of iterations = ";K
7510 PRINT L1$: PRINT: ICONV=1
7520 PRINT "Type 1 to run again with same components and kijs"
7530 PRINT " 0 to return to main menu";:INPUT IENDFG
7540 TOLD#=T#
7550 IF IENDFG=0 THEN RETURN
7555 IF ICONV=0 GOTO 6290
7560 PRINT "Type 1 to use last K values as initial guess"
7570 INPUT " 0 to reinitialize K values ";KVALUE
7580 GOTO 6290
7590 '=========================================================================
7600 'subroutine dewptp
7610 ICOMP=1: KVALUE=0
7620 K=0
7630 SUM#=0#: PDG#=0#: ICONV=0
7640 INPUT "Enter feed temperature (K)";T#
7650 PRINT "Enter estimated dew point pressure (bar)"
7660 INPUT " or 0 for default estimate";P#
7670 IF NC=1 THEN Z#(1)=1#
7680 IF NC=1 GOTO 7770
7690 FOR I=1 TO NC
7700 PRINT "Enter feed mole fraction of ";KK1$(I);
7710 INPUT Z#(I)
7720 SUM#=SUM#+Z#(I)
7730 NEXT I
7740 FOR I=1 TO NC
7750 Z#(I)=Z#(I)/SUM#
7760 NEXT I
7770 GOSUB 1760 'call prcons
7780 IF KVALUE=1 AND P#>.0001 GOTO 7920
7790 FOR I=1 TO NC
7800 DT1#=(1#/T#)-(1#/TB#(I))
7810 DT2#=(1#/TC#(I))-(1#/TB#(I))
7820 DLNP#=LOG(PC#(I))
7830 PS#(I)=EXP(DLNP#*DT1#/DT2#)
7840 PDG#=PDG#+.25#*Z#(I)*PS#(I)
7850 NEXT I
7860 IF P#<.0001# THEN P#=PDG#
7870 FOR I=1 TO NC
7880 IF KVALUE=1 GOTO 7900
7890 XK1#(I)=PS#(I)/P#
7900 XX#(I)=Z#(I)/XK1#(I)
7910 NEXT I
7920 KKK=0
7930 NLOOP=1
7940 IF P#<0# GOTO 8550
7950 IF P#>10000# GOTO 8550
7960 K=K+1
7970 IF K>30 GOTO 8800
7980 GOSUB 1760 'call prcons
7990 FOR I=1 TO NC:ZX#(I)=Z#(I):NEXT I
8000 IFASE=0
8010 GOSUB 740 'call prfuga
8020 FOR I=1 TO NC:F1#(I)=F#(I):NEXT I
8030 ZZ1#=ZZ#
8040 SUMX#=0#
8050 FOR I=1 TO NC
8060 X2#(I)=XX#(I)
8070 SUMX#=SUMX#+XX#(I)
8080 NEXT I
8090 FOR I=1 TO NC
8100 XX#(I)=XX#(I)/SUMX#
8110 NEXT I
8120 SUMX#=0#
8130 KKK=KKK+1
8140 FOR I=1 TO NC:ZX#(I)=XX#(I):NEXT I
8150 IFASE=1
8160 GOSUB 740 'call prfuga
8170 FOR I=1 TO NC:F2#(I)=F#(I):NEXT I
8180 ZZ2#=ZZ#
8190 FOR I=1 TO NC
8200 XX#(I)=F1#(I)*XX#(I)/F2#(I)
8210 X1#(I)=XX#(I)
8220 SUMX#=SUMX#+XX#(I)
8230 NEXT I
8240 ITEST=0
8250 FOR I=1 TO NC
8260 TEST#=ABS(X1#(I)-X2#(I))
8270 IF TEST#>.0001# THEN ITEST=ITEST+1
8280 XX#(I)=XX#(I)/SUMX#
8290 NEXT I
8300 IF KKK>25 GOTO 8360
8310 IF ITEST<=0 GOTO 8360
8320 FOR I=1 TO NC
8330 X2#(I)=X1#(I)
8340 NEXT I
8350 GOTO 8120
8360 S#(NLOOP)=SUMX#
8370 KKK=0
8380 IF (NLOOP-1)>0 GOTO 8420
8390 NLOOP=2
8400 P#=P#-.005#
8410 GOTO 7980
8420 DSDP#=(S#(2)-S#(1))/.005#
8430 IF (ABS(DSDP#)<.000001#) GOTO 8550
8440 DLP#=(S#(1)-1#)/DSDP#
8450 IF K>50 GOTO 8570
8460 IF K<11 THEN DD#=50#
8470 IF K>=11 THEN DD#=10#
8480 IF DLP#>DD# THEN P#=P#+DD#
8490 IF DLP#>DD# GOTO 7930
8500 IF DLP#<-DD# THEN P#=P#-DD#
8510 IF DLP#<-DD# GOTO 7930
8520 IF ABS(DLP#)<.01# GOTO 8590
8530 P#=P#+DLP#+.0025#
8540 GOTO 7930
8550 PRINT "Calculation not converging: one-phase region or poor initial guess"
8560 GOTO 8800
8570 PRINT:PRINT "Dew point calculation did not converge":PRINT
8580 GOTO 8800
8590 YK#=K
8600 TTEST#=(ZZ1#-ZZ2#)^2
8610 IF TTEST#>.00001# GOTO 8650
8620 IF ZZ2#>=.307# THEN P#=P#+25#/SQR(YK#)
8630 IF ZZ2#<.307# THEN P#=P#-25#/SQR(YK#)
8640 KVALUE = 0: GOTO 7870
8650 FOR I=1 TO NC
8660 XK1#(I)=Z#(I)/XX#(I)
8670 NEXT I
8680 PRINT:PRINT L1$
8690 PRINT "Dew point pressure calculation using the "
8700 PRINT "Peng Robinson equation of state"
8710 PRINT "Temperature = ";:PRINT USING L42$;T#;:PRINT" K"
8720 PRINT "Calculated dew point pressure = ";:PRINT USING L42$;P#;:PRINT" bar":PRINT
8730 PRINT "Component Feed Liquid K"
8740 FOR I=1 TO NC
8750 PRINT USING L3$; KK1$(I), Z#(I), XX#(I), XK1#(I)
8760 NEXT I
8770 PRINT:PRINT L4$:PRINT,:PRINT USING L5$; ZZ1#, ZZ2#: PRINT
8780 PRINT "Number of iterations = ";K
8790 PRINT L1$: PRINT: ICONV=1
8800 PRINT "Type 1 to run again with same components and kijs"
8810 PRINT " 0 to return to main menu";:INPUT IENDFG
8820 IF IENDFG=0 THEN RETURN
8825 IF ICONV=0 GOTO 7620
8830 PRINT "Type 1 to use last K values as initial guess"
8840 INPUT " 0 to reinitialize K values ";KVALUE
8850 GOTO 7620
8860 END
8870 '========================================================================
|
Compiler du Basic.
Répondre avec citation
Partager