how to solve matrix differential equation

3 views (last 30 days)
riccardo carrieri
riccardo carrieri on 19 Aug 2019
Edited: David Wilson on 20 Aug 2019
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.

Answers (1)

David Wilson
David Wilson on 20 Aug 2019
Edited: David Wilson on 20 Aug 2019
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 Comments
riccardo carrieri
riccardo carrieri on 20 Aug 2019
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 on 20 Aug 2019
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.

Sign in to comment.

Categories

Find more on Programming in Help Center and File Exchange

Community Treasure Hunt

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

Start Hunting!