Solving simultaneous time-dependent matrix equations

12 views (last 30 days)
I am attempting to solve a set of time-dependent equations which involve matrices. I have tried to use both ode45 and ode15 which both work for equations with one-dimensional variables; as follows:
% Defining constants
A_a = 7.7790;
g1_1 = 0.0019;
g2_1 = 0.0021;
ga_1 = 0.1862;
Hv1 = 2.375;
Hv2 = 2.375;
H1 = 3.9;
CE = 0.82;
Astar1 = 0.9282;
Astar2 = 1.1602;
tspan = [0 9];
Q0 = zeros(3,1);
[t,Q] = ode45(@(t,Q) fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2), tspan, Q0)
plot(t,Q(:,1),'-o',t,Q(:,2), '.-')
function dQdt = fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2)
dQdt = zeros(3,1);
dQdt(1) = (Astar1^2)*((g1_1*Hv1 + ga_1*(H1 + CE))-((Q(3)^2)/A_a^2));
dQdt(2) = (Astar2^2)*((g2_1*Hv2 + ga_1*(CE))-((Q(3)^2)/A_a^2));
dQdt(3) = Q(1) + Q(2);
end
This solves the equation as expected. However when the constants g1_1, g2_1, and ga_1 each become a 1-by-2 matrix, I tried using a similar function (given below), however, it will not solve the code. Is there a way to use matrices in differential equations?
% Defining constants
A_a = 7.7790;
g1_1 = [0.0019 0.0076]
g2_1 = [0.0021 0.0083]
ga_1 = [0.1862 0.3725]
Hv1 = 2.375;
Hv2 = 2.375;
H1 = 3.9;
CE = 0.82;
Astar1 = 0.9282;
Astar2 = 1.1602;
tspan = [0 9];
Q0 = zeros(3,2);
[t,Q] = ode45(@(t,Q) fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2), tspan, Q0)%, options)
plot(t,Q(:,1),'-o',t,Q(:,2), '.-')
function dQdt = fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2)
dQdt = zeros(3,2);
dQdt(1) = (Astar1^2)*((g1_1*Hv1 + ga_1*(H1 + CE))-((Q(3)^2)/A_a^2));
dQdt(2) = (Astar2^2)*((g2_1*Hv2 + ga_1*(CE))-((Q(3)^2)/A_a^2));
dQdt(3) = Q(1) + Q(2);
end

Accepted Answer

Stephan
Stephan on 19 Feb 2019
Edited: Stephan on 19 Feb 2019
Hi,
try:
% Defining constants
A_a = 7.7790;
g1_1 = [0.0019 0.0076]
g2_1 = [0.0021 0.0083]
ga_1 = [0.1862 0.3725]
Hv1 = 2.375;
Hv2 = 2.375;
H1 = 3.9;
CE = 0.82;
Astar1 = 0.9282;
Astar2 = 1.1602;
tspan = [0 9];
Q0 = zeros(3,2);
[t,Q] = ode45(@(t,Q) fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2), tspan, Q0)%, options)
plot(t,Q(:,1),'-o',t,Q(:,2), '.-')
function dQdt = fun4(t,Q,A_a,g1_1,g2_1,Hv1,Hv2,ga_1,H1,CE,Astar1,Astar2)
dQdt = zeros(6,1);
dQdt(1:2) = (Astar1^2)*((g1_1*Hv1 + ga_1*(H1 + CE))-((Q(3)^2)/A_a^2));
dQdt(3:4) = (Astar2^2)*((g2_1*Hv2 + ga_1*(CE))-((Q(3)^2)/A_a^2));
dQdt(5:6) = Q(1) + Q(2);
end
You need to look if the plot has to be changed, since the solution Q now has 6 columns. One way would be:
subplot(3,1,1)
plot(t,Q(:,1),'-o',t,Q(:,2), '.-')
subplot(3,1,2)
plot(t,Q(:,3),'-o',t,Q(:,4), '.-')
subplot(3,1,3)
plot(t,Q(:,5),'-o',t,Q(:,6), '.-')
Also im not sure about the last line in the function - maybe:
dQdt(5:6) = Q(1:2) + Q(3:4);
?
Best regards
Stephan
  2 Comments
JS
JS on 19 Feb 2019
Thank you for this! This works in calculating each value of Q, however it splits Q1 into two answers. The values of g1_1, g2_1, and ga_1 are also time-dependent and were calculated earlier (in reality each is a 41-by-1 matrix rather than a 2-by-1 matrix, though I was trialling the problem with only a 2-by-1 matrix) . I was hoping to obtain time-dependent answers to Q1, Q2, and Q3 which vary in time as with ga_1, g1_1 and g2_1, and will also be of dimension 41-by-1 which can be plotted against time (time is a 41-by-1 matrix).
Perhaps this could be represented better by using tspan = [0:9:9] in the code above (so that time is also a 2-by-1 matrix).
Stephan
Stephan on 19 Feb 2019
Edited: Stephan on 19 Feb 2019
your function as you wrote it will give as many answers for Q1 as g1_1 has elements, because it follows the rules of linear algebra. You add and multiply scalar values to a 1x2 vector, which always results in a 1x2 vector.
If i understood you correctly, you should combine the previous time dependent calculation of the factors with the calculation of Q in one common function to get what you want.

Sign in to comment.

More Answers (0)

Categories

Find more on Chemistry 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!