%%%fmincon%%%

clear all
clc

h = actxserver('UnisimDesign.Application');     %creating active x server, a link between matlab and unisim
uniCase = h.ActiveDocument;                     %entering document in unisim
f = uniCase.Flowsheet;                          %entering flowsheet in unisim
sol = uniCase.Solver;                           %entering solver in unisim

%Make unisim converge before optimization
sol.CanSolve = 0;
f.MaterialStreams.Item('MR1').MassFlowValue = 580;
f.MaterialStreams.Item('1').TemperatureValue = -156.3;
f.MaterialStreams.Item('MR3').PressureValue = 400.9;
f.MaterialStreams.Item('MR5').PressureValue = 2886;
f.MaterialStreams.Item('MR1').PressureValue = 2886;

sol.CanSolve = 1;

global LB UB 

LB = [400 -165  300 0.4];                         %lower and upper bounds 
UB = [1000 -150 900 1.0];

x0_nons = zeros(4,1);
x0 = zeros(4,1);
X = zeros(4,1);

x0_nons(1) = 627.5;                                 %initial molar flow MR before HEX
x0_nons(2) = -155.7;                              %initial temperature in MR before valve
x0_nons(3) = 581.2;                           %initial pressure in MR after valve
x0_nons(4) = 0.63;

scaled = nonscaled_3(x0_nons);                    %scaled variables
x0 = scaled;                             
                             
lb = [0 0 0 0];                                 % scaled lower and upper bounds scaled
ub = [1 1 1 1];


%Optimization with nonlinear constraints
options = optimset('Algorithm','active-set','TolCon',1e-4,'TolFun',1e-4,'Display','iter','MaxIter',125,'MaxFunEvals',1e4,'FinDiffType','central'); % Options for fmincon

[x,fval,exitflag,lambda,grad,hessian] = fmincon(@objective_3,x0,[],[],[],[],lb,ub,@mycon_3,options);

for i = 1:length(x)
    X(i) = (x(i)*(UB(i)-LB(i)))+LB(i);       %scale back variables and function value
end

FVAL = fval*1e6;

% %If no feasible solution, try different alogrithm
% while exitflag == 5
%      x0_nons(1) = X(1)+0.1;
%      x0_nons(2) = X(2);
%      x0_nons(3) = X(3);
%      x0_nons(4) = X(4);
% 
%      scaled = nonscaled_3(x0_nons);                                           
%      x0 = scaled;
%     
%     options = optimset('Algorithm','active-set','TolCon',1e-4,'TolFun',1e-4,'Display','iter','MaxIter',125,'MaxFunEvals',1e4,'FinDiffType','central'); % Options for fmincon
%     [x,fval,exitflag,lambda,grad,hessian] = fmincon(@objective_3,x0,[],[],[],[],lb,ub,@mycon_3,options);
% 
% 
%     for i = 1:length(x)
%     X(i) = (x(i)*(UB(i)-LB(i)))+LB(i);       %scale back variables and function value
%     end
%     FVAL = fval*1e6;
% end

% %Checking for better points
%     
%       x0_nons(1) = X(1)-0.1;
%       x0_nons(2) = X(2)-0.1;
%       x0_nons(3) = X(3)-0.1;
%       x0_nons(4) = X(4)-0.1
% 
%       scaled = nonscaled_3(x0_nons);                                           
%       x0 = scaled;   
%        
%     options = optimset('Algorithm','interior-point','TolCon',1e-4,'TolFun',1e-4,'Display','iter','MaxIter',125,'MaxFunEvals',1e4,'FinDiffType','central'); % Options for fmincon
% 
%     [x,fval,exitflag] = fmincon(@objective_3,x0,[],[],[],[],lb,ub,@mycon_3,options);
% 
%     X = zeros(4,1);
% 
%     for i = 1:length(x)
%         X(i) = (x(i)*(UB(i)-LB(i)))+LB(i);       %scale back variables and function value
%     end
% 
%     FVAL = fval*1e6;




