function varargout = multicomA(t,varargin)

% ***************** Program written by Stathis Skouras in March, 2001 *****************
% ***************** This version: March 21, 2008 by Antonio Araujo ********************
% 
% multicomA- This is a nonlinear model of a continuous distillation column with
%            NT-1 theoretical stages including a reboiler (stage 1) plus a
%            total condenser ("stage" NT). 
%            Model assumptions: 
%            NC components (multicomponent mixture); Component NC is the heavy component
%            constant relative volatilities; 
%            no vapor holdup; 
%            one feed and two products;
%            constant molar flows (same vapor flow on all stages); 
%            Liquid flow dynamics modelled by Franci's Weir Formula.
%            total condenser.  
%
%            The model is based on column A in Skogestad and Postlethwaite
%            (1996). The model has NC*NT states.
%
% Inputs:    t    - time in [min].
%            X    - State, the first NT states are compositions of light
%                   component A, the next NT states are compositions of component B, etc,
%                   compositions of the heavy component are not reported directly  
%                   reboiler/bottom stage is stage (1) and condenser is stage (41). 
%                   The last NT states are liquid holdups in each stage. 
%            U(1) - reflux L,
%            U(2) - boilup V,
%            U(3) - top or distillate product flow D,
%            U(4) - bottom product flow B,
%            U(5) - feed rate F,
%            U(6) - feed liquid fraction, qF.
%             A   - feed compositions, zF.
%         U=[U A] - Disturbances vector
% Outputs:   xprime - vector with time derivative of all the states 
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

parameters
switch nargin
    case 2 % Optimization
        X  = t;
        U  = X(NC*NT+1:end);
        Dd = varargin{1};
    case 3 % Steady state
        X  = t;
        U  = varargin{1};
        Dd = varargin{2};
    case 4 % Dynamic
        X  = varargin{1};
        U  = varargin{2};
        Dd = varargin{3};
end

% Splitting the states
x = X(1:(NC-1)*NT);                    % Liquid compositions from btm to top
M = X((NC-1)*NT+1:NC*NT);              % Liquid hold up from btm to top

% Rearrange elements of composition vector (x) for later use
Iu = reshape(1:(NC-1)*NT,NT,NC-1);
x  = x(Iu);                          

% Inputs
LT = U(1);                             % Reflux
VB = U(2);                             % Boilup

% Disturbances
F  = Dd(1);                             % Feedrate
qF = Dd(2);                             % Feed liquid fraction
zF = Dd(3:end);                         % Feed compositions

% The distillation model
flag = 0; % This flag signals the inclusion of equations for liquid levels. If flag~=0, Simulink is being used and no liquid level is calculated in the code.
model

% Output
xprime = [dxdt; dMdt];
x = x(:); c = [-D + 0; -B + 0; -LT + 0; -VB + VBmin; VB - VBmax]; ceq = [xprime; -x(2*NT) + 0.995];

switch nargin
    case 2 % Optimization
        varargout = {c,ceq};
    case 3 % Steady state
        varargout = {[],xprime};
    case 4 % Dynamic
        varargout = {xprime};
end

