# VLE examples

From ASCEND

REQUIRE "atoms.a4l"; IMPORT "johnpye/extpy/extpy"; IMPORT "vleplot"; MODEL modifiedraoult; (* first declare the variables you are working with *) T IS_A temperature; P, P1S, P2S IS_A pressure; A IS_A factor; gamma1, gamma2 IS_A positive_factor; Pfac IS_A factor; A1, B1, C1 IS_A positive_factor; A2, B2, C2 IS_A positive_factor; x1, x2, y1, y2 IS_A fraction; (* now write down the equations that describe the vapor liquid equilibrium *) (* Antoine's equation for vapor pressure *) vp1: P1S/1000.0{Pa} = exp(A1 - B1/(T/1{K} - C1)); vp2: P2S/1000.0{Pa} = exp(A2 - B2/(T/1{K} - C2)); (* summation equations for mole fractions in vapor and liquid phases *) sumx: x1 + x2 = 1.0; sumy: y1 + y2 = 1.0; (* the modified form of Raoult's Law - with a simple model for the activity coefficients *) vle1: x1*gamma1*P1S = y1*P; vle2: x2*gamma2*P2S = y2*P; (* equations for the activity coefficients *) activity1: gamma1 = exp(A*x2*x2); activity2: gamma2 = exp(A*x1*x1); eqA: A = 2.771 - 0.00523*T/1.0{K}; Pfac*100000.0 {Pa} = P; METHODS METHOD default_self; A1 := 16.59158;B1 := 3643.31;C1 := 33.424; A2 := 14.25326;B2 := 2665.54;C2 := 53.424; END default_self; METHOD specify; FIX A1, B1, C1, A2, B2, C2; FIX T, x1; END specify; METHOD values; T := 318.15 {K}; x1 := 0.25; END values; METHOD on_load; RUN reset; RUN values; RUN default_self; END on_load; METHOD fancyplot; EXTERNAL vleplot(SELF); END fancyplot; END modifiedraoult;

Here is the vleplot.py (most of the work by John)

# python script for calculating Pressure versus y from T and x # binary system, modified raoult's law, simple activity coefficient model import extpy from pylab import * from solverreporter import * def vleplot(self): browser = extpy.getbrowser() ioff() figure() # # I have chosen three temperatures # for T in [320,340,370]: self.T.setRealValue(T) # # collect the data for plotting in two sets of arrays (one for X, one for Y) # I have two sets here - one for P versus y and other for P versus x # XX1 = [] PP1 = [] XX2 = [] PP2 = [] # # change x1 from 0 to 1.0 # there has to be a space between "in" and "[" # for x1 in [0.01,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.40,0.45,0.50,0.55,0.60,0.65,0.7,0.75,0.8,0.85,0.9,0.95,0.99]: self.x1.setRealValue(x1) # # send the pair of values T x1 to the solver # and append the Pressure and y1 (from the solver) to the arrays # the x's are also appended try: browser.sim.solve(browser.solver,SimpleSolverReporter(browser,message="T = %f, x1 = %f" % (T,x1))) XX1.append(float(self.x1)) PP1.append(float(self.P)) XX2.append(float(self.y1)) PP2.append(float(self.P)) except: browser.reporter.reportError('Failed to solve for x1 = %f' % x1) continue ## plot the data plot(XX1,PP1) plot(XX2,PP2) hold(1) ## legend() ion() show() extpy.registermethod(vleplot) #the above method can be called using "EXTERNAL vleplot(self)" in ASCEND. # if you want to see the azeotrope clearly, restrict the calculation to one # temperature