MATLAB code error, LQR

10 views (last 30 days)
Mostafa Sallam
Mostafa Sallam on 1 Jun 2023
Answered: Nathan Hardenberg on 30 Jun 2023
Hello: for the follwoing code:
% Part 1: Define the system parameters and initial conditions
% Define the system parameters
eta = ones(4,1); N = 2*ones(4,1); A = ones(4,1); Z_m = 10*ones(4,1); m = ones(4,1); H = 0.5*ones(4,1);
% Define the torque disturbance
T_d = @(t) [0.05*sin(t); 0.05*cos(t); 0.05*sin(1.3*t); 0.1*sin(t)];
% Define the communication topology
L = [1 0 0 -1;
-1 1 0 0;
0 -1 1 0;
0 0 -1 1];
% Define the initial conditions
x0 = [0; 1; 0.5;
1; 1; 0;
2; 0.1; 1;
3; 0; 1];
% Part 2: Define the LQR controller and simulate the system
% Define the state-space matrices for the system
A_sys = zeros(12,12);
B_sys = zeros(12,4);
C_sys = eye(4,12);
D_sys = zeros(4,4);
for i = 1:4
% Define the state matrices for the ith link
A_i = [0 1 0; -(N(i)/eta(i))*sin(x0(i,1)) -(A(i)/eta(i)) 0; 0 0 0];
B_i = [0; 1/eta(i); 1/(m(i)*eta(i))];
% Update the system matrices with the ith link matrices
A_sys((3*i-2):(3*i),(3*i-2):(3*i)) = A_i;
B_sys((3*i-2):(3*i),i) = B_i;
end
% Define the LQR controller
Q = eye(12);
R = 0.1*eye(4);
[K,~,~] = lqr(A_sys,B_sys,Q,R);
% Define the reference signal
yd = 2;
% Define the simulation time span and time step
tspan = [0 15];
dt = 0.01;
% Define the simulation function
f = @(t,x) multi_agent_system(t,x,K,L,N,eta,A,Z_m,m,H,T_d,yd);
% Simulate the system using theode45 solver and plot the results
% Simulate the system using ode45
[t, x] = ode45(f, tspan, x0);
% Plot the results
figure;
for i = 1:4
subplot(4,1,i)
plot(t, x(:,(3i-2)), 'r-', 'LineWidth', 2); hold on;
plot(t, ydones(size(t)), 'b--', 'LineWidth', 2); hold off;
xlabel('Time (s)');
ylabel(['Position of Link ' num2str(i)]);
legend(['Link ' num2str(i) ' Position'], 'Reference Signal');
end
% Part 3: Define the system dynamics function**
function dxdt = multi_agent_system(t,x,K,L,N,eta,A,Z_m,m,H,T_d,yd)
% Define the state variables
w = x(:,1:3:end);
dotw = x(:,2:3:end);
T = x(:,3:3:end);
% Define the control input
u = -K*x';
% Define the time derivatives of the state variables
d_w = dotw;
d_dotw = -((N./eta).*sin(w) + (A./eta).*dotw - T_d(t)./eta) + u;
d_T = (1./(m.*eta)).*u - (Z_m./(m.*eta)).*dotw -(H./m).*T;
% Define the consensus term
consensus_term = L*x;
% Define the time derivatives of the state variables for each agent
dxdt = zeros(size(x));
for i = 1:4
dxdt(:,(3i-2):(3i)) = [d_w(:,i), d_dotw(:,i), d_T(:,i)] + consensus_term(:,(3i-2):(3i));
end
% Define the reference signal tracking error
error = x(:,1:3:end) - yd;
% Add the reference signal tracking error to the first agent's torque
dxdt(:,3) = dxdt(:,3) - error./eta;
end
it give me the following error:
Error using lqr (line 42)
Cannot compute the stabilizing Riccati solution S for the LQR design.
This could be because:
* R is singular,
* [Q N;N' R] needs to be positive definite,
* The E matrix in the state equation is singular.
Error in New_Model7 (line 44)
[K,~,~] = lqr(A_sys,B_sys,Q,R);
if anyone can help me?

Answers (1)

Nathan Hardenberg
Nathan Hardenberg on 30 Jun 2023
I do not know the exact working behind the "Riccati Solution", but I have had this error before. For me it helped to adjust Q and/or R values slightly.
In your example you could for example try this:
R = 0.1*eye(4);
R(1,1) = 0.09; R(2,2) = 0.99; R(3,3) = 0.101; R(4,4) = 0.102
R = 4×4
0.0900 0 0 0 0 0.9900 0 0 0 0 0.1010 0 0 0 0 0.1020
This works for you example. If you have more Info abour you system you can also adjust Q and R a bit more informed. But for a first shot this should be good enough
Also: There is another error in your code afterwards, but this does not seem to have anything to do with the lqr method

Tags

Products


Release

R2020a

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!