how to solve matrix differential equation

5 次查看（过去 30 天）
riccardo carrieri
riccardo carrieri 2019-8-19
编辑： David Wilson 2019-8-20
hi,
I' m trying to solve this differential equation for education purposes but I'm having difficulty to solve it.
does someone know a simple way to solve it?
? ∗ ?̈ = ?? − (?? + ??) ∗ ?
I attached a pdf where you can find all the information about the program.
thank you very much.

请先登录，再回答此问题。

回答（1 个）

David Wilson
David Wilson 2019-8-20
编辑：David Wilson 2019-8-20
Here's my attempt:
M = [
97227/21875, 0, 82269/109375, 0, 67311/43750, 0, -97227/218750, 0; ...
0, 97227/21875, 0, 82269/109375, 0, 67311/43750, 0, -97227/218750;
82269/109375, 0, 89748/546875, 0, 97227/218750, 0, -67311/546875, 0;
0, 82269/109375, 0, 89748/546875, 0, 97227/218750, 0, -67311/546875;
67311/43750, 0, 97227/218750, 0, 97227/21875, 0, -82269/109375, 0;
0, 67311/43750, 0, 97227/218750, 0, 97227/21875, 0, -82269/109375;
-97227/218750, 0, -67311/546875, 0, -82269/109375, 0, 89748/546875, 0;
0, -97227/218750, 0, -67311/546875, 0, -82269/109375, 0, 89748/546875]
Qv = [ 0; -20651534751534355/422212465065984; 0; -4130306950306871/422212465065984; 0;
-20651534751534355/422212465065984; 0; 4130306950306871/422212465065984];
I couldn't be bothered using your matrices, so I just used some random, symmetric ones. You'll need to type in yours ...
Kt = randn(8,8); Kt = Kt'*Kt; % force symmetric
Kl = randn(8,8); Kl = Kl'*Kl; % force symmetric
L = 5; % ???
e0=[0; 0; 1; 0; L; 0; 0; 1];
Now I'm going to
(1) Convert the 2nd order ODE to Cauchy form, and
(2) Generate the appropriate Mass matrix, (although you can easily invert this since the condition number is reasonable.)
Once all set up, just throw to ode45 to solve it. Once done, extract the first 8 elements & plot them.
%-------------------------------------------------
K = (Kt+Kl);
Ma = eye(16); % augmented mass matrix.
Ma(9:end, 9:end) = M; % Bottom right corner is "M"
z0 = [e0; zeros(8,1)];
odefun = @(t,z) [z(9:16,1); ...
Qv - K*z(1:8,1)];
tspan = [0,5]; % who knows what this should be ???
optns = odeset('Mass',Ma); % Include Mass matrix & keep on LHS.
[t,z] = ode45(odefun, tspan, z0, optns );
plot(t,z(:,[1:8]));
  2 个评论
riccardo carrieri
riccardo carrieri 2019-8-20
thank you very much for your quick reply, but the code don't work
this is the code:
M = eye(16); % augmented mass matrix.
M(9:end, 9:end) = Ma
z0 = [e0; zeros(8,1)];
odefun = @(t,z) [z(9:16,1);inv(Ma)*(Qv - K*z(1:8,1))];
tspan = [0,2]; % who knows what this should be ???
% try explicit method
options = odeset('Mass',M);
[t,z] = ode45(odefun, tspan, z0,options);
plot(t,z(:,[1:8]));
this is the error:
Error using odearguments (line 113)
Inputs must be floats, namely single or double.
Error in ode45 (line 115)
odearguments(FcnHandlesUsed, solver_name, ode, tspan, y0, options, varargin);
Error in prova_del_code (line 135)
[t,z] = ode45(odefun, tspan, z0,options);
do you know why?
David Wilson
David Wilson 2019-8-20
First up, you shouldn't use the inv(Ma) as well as the "Mass" matrix. Keep the mass matrix on the LHS, and drop the inv(Ma) part.
If you are going to use the inv(Ma), then drop the mass matrix part, however this is not generally recommended nor efficient.
You error seems to suggest that one (or more) ofthe arguments to ode45is the wrong type. Check each of the arguments.

请先登录，再进行评论。

类别

在 Help Center 和 File Exchange 中查找有关 Programming 的更多信息

标签

Community Treasure Hunt

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

Start Hunting!

Translated by