function [Ne,Ns]=smoker(alpha,R,xd,xb,zf,q)
% This function uses analytical solution for binary system
% alpha - relative volatility
% R     - reflux ratio
% xd    - distillate comp. mole fraction
% xb    - bottom comp. mole fraction
% zf    - feed composition
% q     - feed thermal condition

%Solve for (xi,yi) and store in vector xi

a(1,1) = R/(R+1); a(1,2) = -1; b(1) = -xd/(R+1);
a(2,1) = q/(q-1); a(2,2) = -1; b(2) = zf/(q-1);
xi=a\b';

%Enriching section

m=R/(R+1); b= xd/(R+1);

%Compute coefficients of quadratic.

c(1) = m*(alpha-1);
c(2) = m + b*(alpha-1) -alpha;
c(3) = b;

%Determine the roots K1, K2 and number of trays
temp = roots(c); K1 =min(temp); K2=max(temp);
Ne=log(((xd-K1)/(xi(1)-K1))*((K2-xi(1))/(K2-xd)))/ ...
   log((1+(alpha-1)*K2)/(1+(alpha-1)*K1));

%Striping section

m=(R*zf+q*xd-(R+q)*xb)/((R+1)*zf+(q-1)*xd-(R+q)*xb);
b= (1-m)*xb;

%Compute coefficients of quadratic.

c(1) = m*(alpha-1);
c(2) = m + b*(alpha-1) -alpha;
c(3) = b;

%Determine the roots K1, K2 and number of trays
temp = roots(c); K1 =min(temp); K2=max(temp);
Ns=log(((xi(1)-K1)/(xb-K1))*((K2-xb)/(K2-xi(1))))/ ...
   log((1+(alpha-1)*K2)/(1+(alpha-1)*K1));