%Calculate optimal cost for disturbance 
%Disturbance in feed varies from -20% to +20% (0.8*u_nom(6) to 1.2*u_nom(6)
%Disturbance  in compositions vary from 
%Disturbance in boilup (Vmax) vary from -20% to +20% (0.8*u_nom(2) to
%1.2*u_nom(2))
%Disturbance in the vapor split ratio Rl varies from 0.1 to 1
%The disturbance qf is always fixed and we consider that the feed is always
% a liquid feed
% Script for plotting loss using exact local method
clear all; close all; clc;

global coldata dyndata hydraulic termodata UU alfa
% COLUMN DATA %%

% Number of stages
coldata.Nt = [12 12 12 12 12 12];
% Initial holdup in reboiler and condenser
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];

hydraulic.k1 = 1/0.063;

% THERMODYNAMIC DATA %%

load Nominal_final.mat   
dyndata.x_nom = x_nom;
dyndata.u_nom = u_nom;
options = optimset('fmincon');
options = optimset(options,'TolFun',1e-8,'TolCon',1e-8,'MaxFunEvals',500000,'Display','iter',...
'Algorithm','active-set');%, 'GradObj','on','Hessian','on');%,...
% all linear equality and inequality constraints given in specs, the
% disturbances can be directly set in specs file
Specs_nullspace
%Relative volatility of Ethanol, Propanol and n-Butanol with respect to the
%heaviest component n-Butanol
alfa=[4.215374389 2.102194688 1];

% Optimization for various disturbances
%Preallocation of cost, which is the sum of impurities in product streams
%
%Joptd=zeros(1,9);                        %optimal cost for feed as disturbance 
     % for i=1:9
      % F= 0.8*u_nom(6):0.05:1.2*u_nom(6);
       %beq(1)= F(i);
      %[x_opt,fval,~,~,~,~,Hessian] = fmincon(@(x) Objective(x),[x_nom;u_nom],Ain,bin,Aeq,beq,[],[],@(x) Petlyuk_nonlcon(x),options);
      %Joptd(i)= fval;
      %end


Joptd=zeros(1,14);                       %optimal cost for boliup as disturbance
for i=1:14
    Vmax= 0.8*u_nom(2):0.049482:1.2*u_nom(2);
    beq(5)= Vmax(i);
[x_opt,fval,~,~,~,~,Hessian] = fmincon(@(x) Objective(x),[x_nom;u_nom],Ain,bin,Aeq,beq,[],[],@(x) Petlyuk_nonlcon(x),options);
Joptd(i)= fval;
end