%initial steady state simulation
%%Three Product Petlyuk Column
%Constant Relative Voltality
% Initializes the model and runs
% The otputs are compositions (A, B and C) of  (D, S, R, D1 and B1)
% Plot profiles give the composition profiles at steady-state
% Plot temperature profile gives temperature profile at steady-state.
%%
clear all; close all; clc;
global coldata dyndata hydraulic termodata UU alfa tend NI i
%% COLUMN DATA %%

% Number of stages in 6 sections of Petlyuk 
coldata.Nt = [12 12 12 12 12 12];
% Initial holdup stages
coldata.M0 = [1 1 1 1 1 1 1 1 1 1 1 1];
%Boiling points of ethanol, propanol and n-butanol
coldata.Tboil = [78.37 97 117.4];
termodata.P = 1.0;
hydraulic.k1 = 1/0.063;
%% THERMODYNAMIC DATA %%  (Realtive volatilities alfa)
alfa=[4.215374389	2.102194688	1] %ethanol propanol butanol

%% load nominal x u
load Nominal_final.mat        % The nominal value obtained from optimization for ...
%minimizing the sum of impurities in the product streams (D, S and R)
UU=[];                        %
i=0;%to initiate Petlyuk_openloop make change in CS by i+1
indmat_Yes_S1
% Make the nominal values of x, u and Temp global
dyndata.x_nom = x_nom; 
dyndata.u_nom=u_nom;          
dyndata.temp_nom=temp_nom;
%% Run this file to get initial values with LV Configuartion, uncomment
% line number 
NI=5; % Number of integerators
options = odeset('MaxStep',10);
%tend=20000;
% at steady state reduce the end time to 5000
tend = 5000;
[t,x]=ode15s(@Simulation_Petlyuk,[0 tend],[dyndata.x_nom;zeros(NI,1)],options);
%% xf states to plot profiles.
xf=x(end,1:225)';
u_nom=UU(end,2:12)';
Plot_Profiles
dxdt = norm(Petlyuk_NLmodel(0, xf, u_nom));
% Print the values of compositions in D (distillate), S (side stream), R (residue)...
%and D1 and R1, (top and bottom product from the prefractionator both  end)
dis_composition                      % saved as init_compSS.mat
%Temperature as measurement Y
[sys Temp]= Petlyuk_NLmodel(t,xf,u_nom);
Yss = Temp;                          % Steady-state Temperature, saved as Yss.mat 
%plot_temp_profile                   % Temperature Profile  in the column,
%use command window to make temp profile


