This is my code
clear; close; clc;
syms a1_head a2_head b hstar
%Parameter Massa
m1 = 8095; % massa train set 1 dalam kg
m2 = 8500; % massa train set 2 dalam kg
g = 10;
c_0_1 = 0.01176;
c_1_1 = 0.00077616;
c_2_1 = 4.48 ;
c_0_2 = 0.01176 ;
c_1_2 = 0.00077616;
c_2_2 = 4.48;
v_0 = 300;
hstar = 120;
a_1 = -1./m1.*(c_1_1 + 2.*c_2_1.*v_0);
a_2 = -1./m2.*(c_1_2 + 2.*c_2_2.*v_0);
a_1_head = 1-(a_1.*hstar);
a_2_head = 1-(a_2.*hstar);
b = 1;
% Model data
A = sym(zeros(4,4));
A(1,2) = a_1_head;
A(3,2) = (a_2_head) - 1; A(3,4) = a_2_head;
display(A);
B = sym(zeros(4,2));
B(1,1) = -b*hstar;
B(2,1) = b;
B(3,2) = -b*hstar ;
B(4,1) = -b; B(4,2) = b;
display(B);
% Q and R matrices for ARE
Q = sym(eye(4)); display(Q);
R = sym(zeros(2,2)); R(1,:) = [1 2]; R(2,:) = [2 3]; display(R);
% Matrix S to find
svar = sym('s',[1 16]);
S = [svar(1:4); svar(5:8); svar(9:12); svar(13:16)];
S(2,1) = svar(2);
S(2,2) = svar(1);
S(2,4) = svar(2);
S(3,1) = svar(3);
S(3,2) = svar(7);
S(3,3) = svar(1);
S(4,1) = svar(4);
S(4,2) = svar(2);
S(4,3) = svar(12);
S(4,4) = svar(1);
display(S);
% LHS of ARE: A'*S + S*A' - S*B*Rinv*B'*S
left_ARE = transpose(A)*S + S*A - S*B*inv(R)*transpose(B)*S;
display(left_ARE);
% RHS of ARE: -Q
right_ARE = -Q;
display(right_ARE);
% subs(left_ARE,{s1, s2, s3}, {1, 5, 10})
% Find S in ARE
syms s1 s2 s3 s4 s7 s12
[Sol_S] = solve(left_ARE == right_ARE,s1,s2,s3,s4,s7,s12,'Real',true,'returnconditions',true)
The result like this
Sol_S =
struct with fields:
s1: [0×1 sym]
s2: [0×1 sym]
s3: [0×1 sym]
s4: [0×1 sym]
s7: [0×1 sym]
s12: [0×1 sym]
I want to find my value s using this code
D = [s1 s2 s3 s4 s7 s12];
D = [S.s1 S.s2 S.s3 S.s4 S.s7 S.s12];
But it said
Error using sym/subsref
Too many output arguments.
Error in ARE (line 78)
D = [S.s1 S.s2 S.s3 S.s4 S.s7 S.s12];

 Respuesta aceptada

Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
The values of S are stored in the struct returned by solve(). The actual value of the symbolic variable does not change. So the correct syntax is
D = [Sol_S.s1 Sol_S.s2 Sol_S.s3 Sol_S.s4 Sol_S.s7 Sol_S.s12];
However, the result shows
Sol_S =
struct with fields:
s1: [0×1 sym]
s2: [0×1 sym]
s3: [0×1 sym]
s4: [0×1 sym]
s7: [0×1 sym]
s12: [0×1 sym]
that MATLAB is not able to find a solution. It is probably the case that a solution does not exist. Also, note that there are 16 equations and only 6 variables, so I guess the system is over-determined.

15 comentarios

if there is 16 variabel and 16 equations it can? if the S matrix
% Matrix S to find
svar = sym('s',[1 16]);
S = [svar(1:4); svar(5:8); svar(9:12); svar(13:16)];
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
because when i run 16 equation and find 16 variables when i debug it's too long to run
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
if 16 variables it too long to run
Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
Yes, It takes a long time too, I noticed it too. However, it seems that you are trying to solve the Algebraic Riccati equation. Using Symbolic toolbox is not a good strategy. You should use a numerical solver given by MATLAB: https://www.mathworks.com/help/control/ref/icare.html. However, this requires the Control system toolbox. If you don't have the toolbox, then try these packages from FEX:
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
ok thank you i will try
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
Editada: Ivan Dwi Putra el 11 de Jun. de 2020
I have read the code in the link, but I don't understand how I can get my s1-s16 value and how i can get the K matriks?
Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
Did you tried the icare() function?
Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
I think this package does not solve the type of Riccati equation you are trying to solve. I couldn't find a more comprehensive package related to the Algebraic Riccati equation on FEX. I guess the control toolbox function is the best option here.
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
in my state space formula i have A,B,W . where is W is placed?I have 4 state (x1,x3 is position), (x2,x4 is velocity) so my C is [0 0 1 0
0 0 0 1] like this right?
ok i just try icare
%Parameter Massa
m1 = 8095; % massa train set 1 dalam kg
m2 = 8500; % massa train set 2 dalam kg
g = 10;
%Parameter Gaya
f1 = 205.10^3; % dalam N
f2 = 302.10^3; % dalam N
c_0_1 = 0.01176;
c_1_1 = 0.00077616;
c_2_1 = 4.48 ;
c_0_2 = 0.01176 ;
c_1_2 = 0.00077616;
c_2_2 = 4.48;
v_0 = 300;
hstar = 120;
a_1 = -1./m1.*(c_1_1 + 2.*c_2_1.*v_0);
a_2 = -1./m2.*(c_1_2 + 2.*c_2_2.*v_0);
a_1_head = 1-(a_1.*hstar);
a_2_head = 1-(a_2.*hstar);
b = 1;
p_1 = -1./m1.*(c_0_1 - c_2_1.*(v_0).^2);
p_2 = -1./m2.*(c_0_2 - c_2_2.*(v_0).^2);
A = [0 a_1_head 0 0;
0 0 0 0;
0 (a_2_head - 1) 0 a_2_head;
0 0 0 0
];
B = [-b.*hstar 0;
b 0;
0 -b.*hstar;
-b b
];
U_t = [f1; f2;];
W = [((a_1 - 1).*v_0) - (p_1.*hstar);
0;
((a_2 - 1).*v_0) - (p_2.*hstar);
((a_1 - 1).*v_0) - (p_1.*hstar);
];
C = [0 0 1 0;
0 0 0 1;];
G = -B*B';
Q = C'*C;
[X1,K1,L1] = icare(A,[],Q,[],[],[],G)
but why the result is like this
X1 =
[]
K1 =
[]
L1 =
0
-0.3153
-1.0030
-120.0037
x1 and k1 no value
Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
I think you should put B directly and ignore the G term. Also, you can calculate R and put it as 4th argument. Try this
[X1,K1,L1] = icare(A,B,Q,R);
It is also possible that solution does not exist. For example, read here: https://pdfs.semanticscholar.org/6b69/c5dada2e3700fa7ac6ada7cf30bf33b68136.pdf
do I set my own R?
i use R
R = [1 2;
2 3;];
[X1,K1,L1] = icare(A,B,Q,R)
But still the X and K still no value
X1 =
[]
K1 =
[]
L1 =
0
0
-0.3134
-1.0087
Ameer Hamza
Ameer Hamza el 11 de Jun. de 2020
In this case, I think the equation is unsolvable for these particular matrices A, B, Q, and R.
Ivan Dwi Putra
Ivan Dwi Putra el 11 de Jun. de 2020
ok thank you for helping me

Iniciar sesión para comentar.

Más respuestas (0)

Categorías

Más información sobre Mathematics en Centro de ayuda y File Exchange.

Productos

Versión

R2020a

Etiquetas

Preguntada:

el 11 de Jun. de 2020

Comentada:

el 11 de Jun. de 2020

Community Treasure Hunt

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

Start Hunting!

Translated by