function sfdist(block)

nu = 6;
ny = 2;
NT = 4;

% block definitions
block.NumInputPorts  = nu;
block.NumOutputPorts = ny;

% Setup port properties to be inherited or dynamic
block.SetPreCompInpPortInfoToDynamic;
block.SetPreCompOutPortInfoToDynamic;


for i = 1:nu
    block.InputPort(i).Dimensions = 1;
    block.InputPort(i).DirectFeedthrough = true;
end

for i = 1:ny
    block.OutputPort(i).Dimensions = 1;
end

block.OutputPort(1).Dimensions = NT;
block.OutputPort(2).Dimensions = NT;


block.SampleTimes = [0 0];

% Setup Dwork
block.NumContStates = NT*2;

block.RegBlockMethod('InitializeConditions',    @InitializeConditions);
block.RegBlockMethod('Outputs',                 @Outputs);
block.RegBlockMethod('Derivatives',             @Derivatives);
block.RegBlockMethod('SetInputPortSamplingMode',@SetInpPortFrameData);


function SetInpPortFrameData(block, idx, fd)
  
block.InputPort(idx).SamplingMode = fd;
for i = 1:block.NumOutputPorts
    block.OutputPort(i).SamplingMode = fd;
end



function InitializeConditions(block)

load steadystate % y0
block.ContStates.Data = y0;

function Outputs(block)
NT=4; NF = 2;

% State variables
x = block.ContStates.Data(1:NT);
M = block.ContStates.Data(NT+1:end);

% Disturbances
F   = block.InputPort(1).Data;
zF  = block.InputPort(2).Data;

% Manipulation
LT  = block.InputPort(3).Data;
VB  = block.InputPort(4).Data;
D   = block.InputPort(5).Data;
B   = block.InputPort(6).Data;

% Output variables
block.OutputPort(1).Data = x;
block.OutputPort(2).Data = M;

function Derivatives(block)

  NT=4; NF = 2;
% State variables
  x = block.ContStates.Data(1:NT)';
  M = block.ContStates.Data(NT+1:end)';

% Disturbances
  F   = block.InputPort(1).Data;
  zF  = block.InputPort(2).Data;

% Manipulation
  LT  = block.InputPort(3).Data;
  VB  = block.InputPort(4).Data;  
  D   = block.InputPort(5).Data;
  B   = block.InputPort(6).Data;

% Model
  
  i  = 1:NT-1;
  V = VB*ones(1,4);
  L(2:NF)      = LT+F;         % Liquid flows below feed (Kmol/min)
  L(NF+1:NT)   = LT*ones(1,2); % Condenser's liquid flow (Kmol/min)
   % VLE
     alpha = 4.78;
     y(i) = alpha*x(i)./(1+(alpha-1)*x(i));
     y(NT) = x(NT); 

   
   


   i  = 2:NT-1;

    % Holdups inside column will be constant because of assumption of flow V and L given above
     dMdt(i)   = L(i+1) - L(i) + V(i-1) - V(i); 
     dMdt(NF)  = dMdt(NF) + F;      
     
     dMxdt(i)  = L(i+1).*x(i+1) - L(i).*x(i) + V(i-1).*y(i-1) - V(i).*y(i);
     dMxdt(NF) = dMxdt(NF) + F*zF;
     
      dMdt(NT)  = V(NT-1) - LT -D; %Condenser holdup
     dMxdt(NT) = V(NT-1)*y(NT-1) - LT*x(NT) -D*x(NT);
     
     dMdt(1)   = L(2) - V(1) -B; %Reboiler holdup
     dMxdt(1)  = L(2)*x(2) - V(1)*y(1) -B*x(1);

     i=1:NT;
     dxdt(i) = (dMxdt(i)-x(i).*dMdt(i))./M(i);
     
     
  block.Derivatives.Data = [dxdt';dMdt'];