How to transform a second order ODE to the format xddot=A.X+B.U

4 ビュー (過去 30 日間)
Ítalo Almeida
Ítalo Almeida 2023 年 1 月 30 日
コメント済み: Sam Chak 2023 年 2 月 8 日
I have a system of ODEs of the type that are currently simbolic functions and I want to find a way to obtain it's matrices automatically cause the system is too big to do it by hand.
  1 件のコメント
Ítalo Almeida
Ítalo Almeida 2023 年 2 月 8 日
編集済み: Ítalo Almeida 2023 年 2 月 8 日
the solution i found was to use equationsToMatrix(equations,list). I created a 'list' with all the simbolic variables x and u and it gave me the matrix equivalent to then I just had to split the matrix into two. As i said in the answers bellow, I don't know any built in matlab function to solve the diferential riccati equation, so I'll have to use the system matrices(Riccati diferential equations's parameters) to solve it with ode45.

サインインしてコメントする。

回答 (3 件)

Star Strider
Star Strider 2023 年 1 月 30 日
See if using odeToVectorField and matlabFunction will do what you want.
Those should work if you want to solve them using the MATLAB ordinary differential equation integrators.
If you want to use them with the Control System Toolbox, that will require a different approach.
  2 件のコメント
Ítalo Almeida
Ítalo Almeida 2023 年 2 月 6 日
My intention was to simulate with the matlab integrator, yes. But I needed the matrices to calculate the control law. I took a look in the Control System Toolbox MATLAB page and think this might be the way to go. My end goal is to find a way to calculate the LQR optimal control of the plant.
Star Strider
Star Strider 2023 年 2 月 6 日
編集済み: Star Strider 2023 年 2 月 7 日
The linear control systems assume that all the derivatives are first-degree. The odeToVectorField function will put them in that form, and the second ‘Subs’ output will give the row order of the variables.
Adding:
sympref('AbbreviateOutput',false);
after the initial syms call might also be helpful.
EDIT — (7 Feb 2023 at 11:54)
Another option (that I do not often use and so forgot about) is the odeFunction function. That might be more appropriate here.
.

サインインしてコメントする。


Paul
Paul 2023 年 2 月 6 日
編集済み: Paul 2023 年 2 月 6 日
Hi Ítalo,
Let p = x and q = xdot. Then the standard, first order state space model is
[pdot;qdot] = [zeros(n) eye(n); A zeros(n)] * [p ; q] + [0*B; B] * u;
For example, let A = 2 and B = 3. The output of the system is x.
A = 2; B = 3;
hss = ss([0 1;A 0],[0*B;B],[1 0],0)
hss = A = x1 x2 x1 0 1 x2 2 0 B = u1 x1 0 x2 3 C = x1 x2 y1 1 0 D = u1 y1 0 Continuous-time state-space model.
htf = tf(hss)
htf = 3 -------------------- s^2 - 2.22e-16 s - 2 Continuous-time transfer function.
We see from the transfer function that the input/ouput ode is: ydddot = 2*y + 3*u, exactly as it should be because y = x.
You might also be interestes in mechss and related functionality to represent your model. You can always convert it to a full ss model for input to lqr if lqr doesn't accept a sparse second order model.

Sam Chak
Sam Chak 2023 年 2 月 7 日
The odeToVectorField() should work. But since you already have the linear 2nd-order system in this ODE form , then you should be able to readily transform it to the State-space form through some basic matrix operations:
, where and
Then, the optimal gain can be obtained from the lqr() function.
Check this example if this is something you are look for. You can also try applying the sparse 2nd-order state-space model as recommended by @Paul.
A = magic(4)
A = 4×4
16 2 3 13 5 11 10 8 9 7 6 12 4 14 15 1
B = diag([2, 3, 5, 7])
B = 4×4
2 0 0 0 0 3 0 0 0 0 5 0 0 0 0 7
n = size(A)
n = 1×2
4 4
Aa = [zeros(n) eye(n);
A zeros(n)]
Aa = 8×8
0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 1 16 2 3 13 0 0 0 0 5 11 10 8 0 0 0 0 9 7 6 12 0 0 0 0 4 14 15 1 0 0 0 0
Bb = [zeros(n); B]
Bb = 8×4
0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 2 0 0 0 0 3 0 0 0 0 5 0 0 0 0 7
K = lqr(Aa, Bb, eye(2*n), eye(n))
K = 4×8
10.4604 0.7265 1.4724 6.7601 2.8207 -0.0576 0.0713 0.5318 1.1589 6.4533 4.7276 2.7826 -0.0863 2.1027 0.4705 0.2195 2.6327 4.3779 4.8063 3.6314 0.1784 0.7842 1.5818 0.3147 8.5926 4.1670 4.4758 7.5283 1.8613 0.5123 0.4405 1.6495
  2 件のコメント
Ítalo Almeida
Ítalo Almeida 2023 年 2 月 8 日
the lqr matlab function does not work for me beacause I have an specific finite end time to regulate the states and also I am trying to get the states to a non zero value, so it is the tracking version of the lqr. Because of that I need to solve the diferential riccati equation and not the algebraic one. my ideia was to solve it using the ode45 solver, thats why i needed the coeficient matrices of the system.
Sam Chak
Sam Chak 2023 年 2 月 8 日
Sounds like a finite horizon problem. You can check if you can apply one of the QP solvers for solving the Riccati Differential Equations. You can find some info here:

サインインしてコメントする。

カテゴリ

Help Center および File ExchangeOrdinary Differential Equations についてさらに検索

製品


リリース

R2021b

Community Treasure Hunt

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

Start Hunting!

Translated by