run_dist_dyn.m

% =========================Run Dynamic Distillation===================================
%
% Solve for the transient stage compositons in an ideal
% binary distillation column using ode45.
%
%  alpha		= relative volatility (1.5)
%  ns			= total number of stages (41)
%  nf			= feed stage (21)
%  feedi   	= initial feed flowrate (1)
%  zfeedi  	= initial feed composition, light comp (0.5)
%  qf      	= feed quality (1 = sat'd liqd,0 = sat'd vapor) (1)
%
% refluxi = initial reflux flowrate (2.706)
% vapori  = initial reboiler vapor flowrate (3.206)
% md      = distillate molar hold-up (5)
% mb      = bottoms molar hold-up (5)
% mt      = stage molar hold-up (0.5)
%
%
% By Giovani Tonel ([e-mail]) October 2006
  
  disp(' ');
  global DIST_PAR
  
% xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx
title1	= 'Distillation Parameters';

ind1	= {'Relative volatility (1.5)','Total number of stages (41)','Feed stage (21)',...
      'Initial feed flowrate (1)','Initial feed composition (0.5), light comp (0.5)',...
      'Feed quality (1 = sat_d liqd,0 = sat_d vapor)(1)',...
      'Initial reflux flowrate (2.706)','Initial reboiler vapor flowrate (3.206)',...
      'Distillate molar hold-up (5)','Bottoms molar hold-up (5)','Stage molar hold-up (0.5)',...
      'Magnitude step in reflux (0.02706)','Time of reflux step change (5)'};

setpar	=	{'1.5','41','21','1','0.5','1','2.706','3.206','5','5','0.5','0.02706','5'};


prompt		=	ind1;
def			=	setpar;
dlgTitle	=	title1;
lineNo		=	1;
M	=	inputdlg(prompt,dlgTitle,lineNo,def);
       
% Converte os dados para os formatos double array:


% a=char(M(1));  
% b=str2num(t0);

DIST_PAR	=	str2num(char(M));
disp(' ');      

% the initial guess is a linear interpotion between 0.01 and 0.99

% +++++++++++++++++++++++++++++++++++++++

% x0 = ones(1,41)*0.15;

x0 = 0.01:0.98/40:0.99;

x0 = x0';

% +++++++++++++++++++++++++++++++++++++++

  


%--------------------------------------------------------------------------------------
% Showing Set Parameters

lengthx0 = length(x0);

if lengthx0 == DIST_PAR(2),
   disp(' ');
   disp(['Number of initial conditions: ' num2str(lengthx0)]);
   disp(' ');

disp('       ==================================================================');

out1= 'alfa ns nf feedi zfeedi qf refluxi vapori md mb mt st_reflux t_reflux';
printmat(DIST_PAR',title1,'Values',out1)
disp(' ');
disp('       ==================================================================');

op=menu('Choose the final time (tf)', ...
   't0=0 tf=10',...
   't0=0 tf=50',...
   't0=0 tf=100',...
   't0=2 tf=300');

switch op
    case 1
       t0		= 0 ; 			% s 			 
			tf 	= 10 ; 			% s	 
    case 2
       t0		= 0 ; 			% s 			 
			tf 	= 50 ; 		% s	  
    case 3
       t0		= 0 ; 			% s 			 
       tf 	= 100 ; 		% s
       disp(' ');
       disp('This can request some seconds...')
       disp(' ');
    
 	 case 4
       t0		= 0 ; 			% s 			 
       tf 	= 300 ; 		% s
       disp(' ');
       disp('This can request very much seconds...')
       disp(' ');
  
      end


%--------------------------------------------------------------------------------------
% Solving

[t,x] = ode45('dist_dyn',0,tf,x0); 

button = questdlg('Plot Grafic ?', ...
                  'Plot','Yes','No','Yes');
               
switch button
case 'Yes',
   disp(' ');
   disp('Ploting...')
   
   
   plot(t,x); title('Transient stage compositons - Distillation Column');
   xlabel('Time'); ylabel('Compositons')
   
case 'No',
   disp('No Ploting!!!');
end

  
  %--------------------------------------------------------------------------------------
  
else
   disp(' ');
   disp('Not enough number of initial conditions!')
   
end