function Y = WILSON(T); 

% This function uses the Wilson model at given values of temperature (T), pressure (P) and liquid composition (X)
% Returns the vector of vapour compositions (Y)

global X NC P R A V W 

% Antoine equation where Psat [mmHg] and temperature T [K]
for i=1:NC
  log10Ps = A(i,1) - A(i,2) / (T - 273.15 + A(i,3));
  Ps(i)   = 10^log10Ps;
end

% Starting Wilson model
% Calculate the matrix of Wilson coefficients LAMBDA_ij
for i=1:NC
  for j=1:NC
    if i==j
     LAMBDA(i,j)=1;
    else
     LAMBDA(i,j)=(V(j)/V(i))*exp(-W(i,j)/(R*T));
    end
  end
end

for i=1:NC
  for k=1:NC
    RATIO(k)=LAMBDA(k,i)/(X'*LAMBDA(k,:)');
  end
LNGAMMA(i)=-log(X'*LAMBDA(i,:)')+1-X'*RATIO';          % logarithm of the activity coefficients
end
GAMMA = exp(LNGAMMA)';                                  % activity coefficients
% Finishing Wilson model

% VLE formula, P*Y(i)=X(i)*GAMMA(i)*Ps(i) 
Y=Ps'.*X.*GAMMA./(P*760);                              % Vector of vapour compositions