%Minimum Singular Value Rule
%Scale the input using linear model such that each input has same effect on
%the cost J, get UscI.j= 1/(abs(Juu))^1/2, 3x1 matrix
%Load Juu.mat

clear all; close all; clc;

global coldata dyndata hydraulic termodata UU alfa

% The script for F matrix (optimal sensitivity)is script_Fmatrix

%script_Fmatrix

% instead of writing above script load F1 directly
load F1.mat
% The value of Juu as calculated using finite-diffrence approximation
% method
                    
load Juu
u_sc = [1/sqrt(Juu(1,1)); 1/sqrt(Juu(2,2));1/sqrt(Juu(3,3))];

%optimal variation in candidate controlled variables w.r.t. disturbances
%load F1
dd1= 1e-5;dd2=1e-5;dd3=1e-5;dd4=1e-5;dd5=1e-4;
delYopt = F1*[dd1; dd2; dd3; dd4; dd5];


%for controlled variable which is temperature here , the implementaion
%error 0.5 degree C added, no control error because of PI controller
nm = 0.5*ones(72,1);
CscI = abs(delYopt)+nm;


%Scaling matrix for the controlled variables and inputs
Dc = diag(CscI);
Du = diag (u_sc);

% The scaled model
load GY
G = inv(Dc)*GY*Du;
%load G
%load Juu
%Jud=inf;
%Gd1=inf;
load Wn 
load Wnd
%Branch and bound algorithm to calculate the worst minimu loss
%[B,sset,ops,ctime,flag]=b3wc(G,Gd1,Wnd,Wn,Juu,Jud,inf,20,F1)
%[B,sset,ops,ctime,flag]=pb3wc(GY,Gd1,Wnd,Wn,Juu,Jud,3,inf,10);
%[Loss,sset1,ops1,ctime,flag]=pb3wc(GY,Wnd,Wn,Juu,F1,3,20000,40);
% The branch and bound algorithm to calculate maximum of the minimum
% singular values to select the best combination of measurement which are
% three temperatures here.
%[B,sset,ops]=b3msv(Q,nc)
% G is scaled gain matrix (Ivar Halvorsen, Sigurd Skogestad, John Morud and
% V. Alstad Paper 2003,"optimal selection of controlled variables", INCER
% MSV saved as MSV and sset saved as sset, total 20 sets of measurements
[MSV,sset,ops]=b3msv(G,20);

% Calculate the loss using the value of B saved as MSV.mat(minimum singular
% values) 

Loss_MSV = norm(Juu,2)*0.5*1/norm(G,2)^2
% The end of MSV

% RGA Calculation (pairing of T10, T38, T60 with L, S and RL)
load GY
Grga = [GY(10,:);GY(38,:);GY(60,:)];
R = inv(Grga') .* Grga ;
%Relative gain array is save as R, [-0.000336,1.88e-05,1.00031712216556;0.987
%,0.012,-0.000363;0.01268,0.9872,4.6528e-05;]


