diff --git a/CheckAndersLMI.m b/CheckAndersLMI.m deleted file mode 100644 index 325b68e0705e1b8b228429b4a51f97b8710af8ed..0000000000000000000000000000000000000000 --- a/CheckAndersLMI.m +++ /dev/null @@ -1,195 +0,0 @@ -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -% -% This code simulates distributed minimax adaptive control algorithm for -% uncertain networked systems. -% -% Copyrights Authors: 1) Venkatraman Renganathan - Lund University, Sweden. -% 2) Anders Rantzer - Lund University, Sweden. -% 3) Olle Kjellqvist - Lund University, Sweden. -% -% Email: venkatraman.renganathan@control.lth.se -% anders.rantzer@control.lth.se -% olle.kjellqvist@control.lth.se -% -% Date last updated: 20 October, 2023. -% -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% -%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - -% Make a fresh start -clear all; close all; clc; - -% set properties for plotting -set(groot,'defaultAxesTickLabelInterpreter','latex'); -set(groot,'defaulttextinterpreter','latex'); -set(groot,'defaultLegendInterpreter','latex'); - -% Flag deciding whether to generate new data or load existing data -% dataPrepFlag = 1: Generates new network data -% dataPrepFlag = 0: Loads existing network data -dataPrepFlag = 1; - -% When dataPrepFlag = 1: Generate new data for network -if(dataPrepFlag) - - disp('Generating new data for network as dataPrepFlag = 1.'); - - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - % Define the network data using graph structure - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - % Specify Number of nodes - N = 6; - % Specify Number of edges - M = N - 1; - % Create a random graph with N nodes and M edges - [adjMatrix, I, graphVariable] = GenerateRandomTreeGraph(N); - % Get number of edges on graph - numEdges = size(graphVariable.Edges.EndNodes, 1); - fprintf('Generated random graph is acyclic and has only one component \n'); - - % Plot the graph - large graph can be time consuming & may not be informative - if N <= 300 - plot(graphVariable); - else - disp('The graph is too large to display.'); - end - - % Infer state and input dimensions - [n,m] = size(I); - - % Get the degree vector - degreeVec = sum(adjMatrix, 1)'; - - % Get one of the node indeces with maximum degree for plotting - [~, maxDegreeNodeIndex] = max(degreeVec); - - % Form input matrix B with parameter b and matrix I such that B/b = I - b = 0.1; - B = b*I; - - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - % Prepare the dynamics matrix A for each Model - %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - - % Specify the number of possible plant network parameters - numModels = 2; - - % Populate all the possible plant models - for dIter = 1:numModels - % Placeholders for a matrix for each component - ai = zeros(n, 1); - % Flag for stopping while loop - dynamicsPreparationFlag = 1; - % Counter for populating ai vector - iter = 1; - % Loop through until ai^2 + bb^{T} < ai is satisfied - while dynamicsPreparationFlag - % Generate a random ai - ai(iter, 1) = rand; - % Check condition - if(ai(iter, 1)^2 - ai(iter, 1) + 2*b^2*degreeVec(iter, 1) < 0) - iter = iter + 1; - end - % if counter exceeds limit, form the A matrix & stop while loop - if(iter > n) - A = diag(ai); - dynamicsPreparationFlag = 0; - end - end - % Prepare A and B matrices - AMatrices{dIter} = A; - BMatrices{dIter} = B; - - end - - AVector = zeros(n, numModels); - for dIter = 1:numModels - AVector(:, dIter) = diag(AMatrices{dIter}); - end - Arowvectors = {}; - for i = 1:N - Arowvectors{i} = AVector(i, :); - end - Acombs = combvec(Arowvectors{:}).'; - % Find min & max a values for each node from all possible numModels values - minAVector = min(AVector, [], 2); - maxAVector = max(AVector, [], 2); - % Form the Max A matrix using maxAVector - maxAMatrix = diag(maxAVector); - % Find the network level global minimum and maximum a value - netMinA = min(minAVector); - netMaxA = max(maxAVector); - - % Infer the number of possible plants - numPlants = numModels^(N); - - % Compute Hinfty control for each possible plants - for i = 1:numPlants - APossibles{i} = diag(Acombs(i, :)); - BPossibles{i} = B; - % Compute H_infinity control as u = Kx, with K = B'*(A-I)^{-1} - KinfPossibles{i} = BPossibles{i}'*inv(APossibles{i} - eye(n)); - end - disp('Finished Computing H_infinity Control Gains'); - - % Check Anders LMI for each possible plant - gamma = 7; - P_non = (2/(1 - netMaxA))*eye(n); - P_non_inv = ((1 - netMaxA)/2)*eye(n); - P_wt = inv(P_non_inv - gamma^(-2)*eye(n)); - Q = eye(n); - R = eye(M); - happy = 1; - - for k = 1:numPlants - Ak = APossibles{k}; - if(happy == 0) - break; - end - for l = 1:numPlants - Al = APossibles{l}; - if(happy == 0) - break; - end - for p = 1:numPlants - if(happy == 0) - break; - end - if(k ~= l && l == p) - continue; - end - Kp = KinfPossibles{p}; - A_cl_kp = Ak + B*Kp; - A_cl_lp = Al + B*Kp; - AndersRHS = Q + Kp'*R*Kp + 0.25*(A_cl_kp+A_cl_lp)'*P_wt*(A_cl_kp+A_cl_lp) - 0.25*gamma^(2)*(A_cl_kp-A_cl_lp)'*(A_cl_kp-A_cl_lp); - AndersLMI = P_non - AndersRHS; - if(all(eig(AndersLMI)) >= 0) - happy = 1; - else - happy = 0; - end - end - end - end - - - % Finished Checking Anders LMI - if(happy == 1) - disp('Choice of Pkl matrices satisfy Anders LMI.'); - else - disp('Choice of Pkl matrices do not satisfy Anders LMI.'); - end - - % Save the generated network data into mat file. - disp('Saving the generated small network data.'); - %save('smallNetworkDynamicsData.mat'); -else - disp('Loading existing data for small network as dataPrepFlag = 0.'); - %load('smallNetworkDynamicsData.mat'); -end - - - -