kill(all)$

analytical method : binary distillation
peter.vlasschaert@gmail.com,05/04/2017
part1: enriching-section 'operating line',revision1

alpha :constant relative volatility

assume(alpha >= 1);

p1:hyperbole       : equilbrium equation.
p2:linear equation : operating line.

p1:y=(alpha*x)/(1+(alpha-1)*x);
p2:y=m*x+b;

removes (p1,p2),solve for x

p3:ratsimp(rhs(p1)-rhs(p2))=0;
p4:(-1)*num(part(p3,1))=0;

p5:coeff(part(p4,1),x^2);
p6:coeff(part(p4,1),x);
p7:ev(part(p4,1),x=0);

solve : quadratic equation ,two roots.

p8:float(solve([p4], [x])[1]);
p9:float(solve([p4], [x])[2]);

first root : assume (alpha >=1)
           : x=(0 = < first root =< 1)

first  value :interval 'x=0' in p4

p10:ev(p4,x=0);

second value :interval 'x=1' in p4

p11:factor(ratsimp(ev(p4,x=1)));

reason : alpha > 0

p12:part(p11,1,2)=0;

rem : A=(alpha-1)*m          -> p5
      B=m+(alpha-1)*b-alpha  -> p6
      C=b                    -> p7
 A*x^2+B*x+C = 0

same: p11 : A+B+C=0

p13:factor(ratsimp(p5+p6+p7))=0;

which root  p8 or p9 stay 
interval :        0 < x = 'root' < 1 (negative root,p8)
may be other root,x > 1 (positive root,p9)'for enriching-section'

x=1,y(1)=m*x+b -> y=m+b, because  y < 1 -> m+b < 1
y= mol fraction of vapor  , 0 = < x = < 1 (p1)
x= mol fraction of liquid , 0 = < y = < 1 (p1)
y = x       :  azeotrope                   ('p1',alpha = 1)
reason :
y <> x = 1  :  distillation => y < 1 => m+b < 1,(alpha > 1)

rem: y=m*x+b , general equation
-------------------------------
x=0 =>  b>0
x=1 =>  m+b < 1
rem:
why : ' b > 0 '
see : above 'given'
-------------------------
 B=m+(alpha-1)*b-alpha
 A=(alpha-1)*m > 0
 C= b
-------------------------
 B+A+C < 0
 B+A-A+C-C  < -A-C
(A+B+C)-A-C < -A-C
 B          < -A-C -> -A-C <0  , A>0,C=b -> b > 0
example : appendix 'example'
exit : root x > 1 ( root 'p8' )
( A+B+C=alpha*(m+b-1) < 0) because ' m+b < 1 'and' alpha > 0'

p14:subst(1, x, p2);

assume(b>0,part(p14,2)<1);

m:0.501;
b:0.43;
alpha:5.5;

negativ root,quadratic equation (p8)

assume(b>0,part(p14,2)<1);

p15: first  intersection 'green line' and 'blue line'
p16: second intersection 'green line' and 'blue line'

p15:ev(part(p8,2));
p16:ev(part(p9,2));

y=x ,45 degree line.'ref line'

p15:y=x;

find : intersection between 'green line ' and 'red line'
green line : y=m*x+b  ' enriching-section operating line.
red   line : y=x      ' above in the tower y=x=xd'

p16:solve([p2, p15], [x,y])$
p17:xd=float(p16[1][1]);
p18:yd=float(p16[1][2]);

plot2d([part(p1,2),part(p15,2),part(p2,2)],
     [x,-.5,1.3],[y,-.5,1.3],
 [plot_format, gnuplot])$