%TEST_ACCURACY_AREA

clear all
clc

global h
h=actxserver('UnisimDesign.Application');

load x0
load ua
load disturbances

x0 = x0_area ;
lb = [4240 ; 500 ; 101/3600 ; 40/3600 ; 0/3600 ;  6/3600 ; -159.15 ; -105 ] ;
ub = [4440 ; 537 ; 103/3600 ; 45/3600 ; 4/3600 ; 10/3600 ; -158.00 ;  -85 ] ;

% scale_x = diag(2.^(fix(log2(ub - lb)))) ;
scale_x = eye(8) ;
x0 = scale_x\x0 ;
lb = scale_x\lb ;
ub = scale_x\ub ;

parameters.scale_x = scale_x ;

x0_red = x0(1:6) ;
T0 = x0(7:8);

T0(1) = T0(1)*1.00 ;
T0(2) = T0(2)*1.00 ;


% parameters.solver = 'fsolve';
% options = optimset('Display','iter','MaxIter',100,'TolFun',2e-4,'Algorithm','levenberg-marquardt') ;
% [T,fval,exitflag] = fsolve(@recycles,T0,options,x0_red,d,parameters) ;
% % [T,resnorm,res_final,exitflag] = lsqnonlin(@recycles,T0,[],[],options,x0_red,d,parameters) ;
% x_full = scale_x*[x0_red ; T];
% T_final = x_full(7:8) 
% res_final = recycles(T,x0_red,d,parameters) 
% exitflag

% parameters.solver = 'fmincon' ;
% options_II = optimset('Display','iter','Algorithm','active-set','MaxIter',50,'TolCon',1e-3);
% [T,fval,exitflag] = fmincon(@dummyJ,T0,[],[],[],[],lb(7:8),ub(7:8),@recycles,options_II,x0_red,d,parameters) ;
% x_full = scale_x*[x0_red ; T];
% T_final = x_full(7:8) 
% [c res_final] = recycles(T,x0_red,d,parameters) 
% exitflag

% Have to check the gradients here.

% jac = zeros(2) ;
% res0 = recycles(T0,x0_red,d,parameters) ;
% for i=1:2
%     T= T0 ;
%     T(i) = T(i)*1.01 ;
%     res = recycles(T,x0_red,d,parameters) ;
%     jac(:,i) = (res - res0)/(T(i) - T(i)/1.01) ;
% end

% T1 = linspace(lb(7),ub(7),10) ;
% T2 = linspace(lb(8),ub(8),10) ;
% resmat = zeros(10) ;
% for i=1:10
%     disp(i)
%     for j=1:10
%         T = [T1(i) ; T2(j)];
%         res = recycles(T,x0_red,d,parameters) ;
%         normres= res'*res ;
%         resmat(i,j) = normres ;
%     end
% end
% 
% figure(1),surf(T1,T2,resmat)

%Test how often we can converge to a solution when changing x using
%different solvers. First use fsolve:

parameters.solver = 'fsolve' ;
fsolve_options = optimset('Display','off','MaxIter',50,'TolFun',1e-2) ; % ,'Algorithm','levenberg-marquardt'
areafail = 0 ;
for i=1:6
    disp(['Changing ' num2str(i) 'th element of x'])
    x = x0_red ;
    xpert = linspace(lb(i),ub(i),25) ;
    for j=1:25
        x(i) = xpert(j) ;
        [T,fval,exitflag] = fsolve(@recycles,T0,fsolve_options,x,d,parameters) ;
        if exitflag > 0
            disp('Solution found')
        elseif exitflag == 0
            disp('Maximum iterations used')
            areafail = areafail + 1 ;
        elseif exitflag < 0
            disp('Did not converge')
            areafail = areafail + 1 ;
        end
    end
end

% fmincon_options = optimset('Algorithm','interior-point','TolCon',1e-1,'Display','off','MaxIter',50,'MaxFunEvals',500);
% parameters.solver = 'fmincon' ;
% areafail = 0 ;
% for i=1:6
%     disp(['Changing ' num2str(i) 'th element of x'])
%     x = x0_red ;
%     xpert = linspace(lb(i),ub(i),25) ;
%     for j=1:25
%         x(i) = xpert(j) ;
%         [T,fval,exitflag] = fmincon(@dummyJ,T0,[],[],[],[],lb(7:8),ub(7:8),@recycles,fmincon_options,x,d,parameters) ;
%         if exitflag > 0
%             disp('Solution found')
%         elseif exitflag == 0
%             disp('Maximum iterations used')
%             areafail = areafail + 1 ;
%         elseif exitflag < 0
%             disp('Did not converge')
%             areafail = areafail + 1 ;
%         end
%     end
% end

% lsqnonlin_options = optimset('Display','off','MaxIter',100,'TolFun',1e-4,'Algorithm','levenberg-marquardt') ;
% parameters.solver = 'lsqnonlin' ;
% areafail = 0 ;
% for i=1:6
%     disp(['Changing ' num2str(i) 'th element of x'])
%     x = x0_red ;
%     xpert = linspace(lb(i),ub(i),25) ;
%     for j=1:25
%         x(i) = xpert(j) ;
%         [T,resnorm,residual,exitflag] = lsqnonlin(@recycles,T0,[],[],lsqnonlin_options,x,d,parameters) ;
%         if exitflag > 0
%             disp('Solution found')
%         elseif exitflag == 0
%             disp('Maximum iterations used')
%             areafail = areafail + 1 ;
%         elseif exitflag < 0
%             disp('Did not converge')
%             areafail = areafail + 1 ;
%         end
%     end
% end


Reset_MCHE_area(x0_area,d,h) ;