清空筛选条件
清空筛选条件

ode45 with matrix 1st order ode

1 次查看（过去 30 天）
Haya M
Haya M 2020-11-14
评论： Haya M 2020-11-15
Hi everyone, I'm trying to solve a first order ode in a matrix form using ode45:
where and
$Q=\begin{pmatrix} \sin(x) & 0 \\
0 & \cos(x) \\
\end{pmatrix}$ on .
Here is my code:
clear all
z = 0, %parameter
n = 2;
T0 = eye(n,n)
xspan = [0 5*pi];
opts = odeset('RelTol',1e-8,'AbsTol',1e-8);
[Tl] = ode45(@(x,T) odeTL(x,T,z,n),xspan,T0,opts);
[q,~] = qr(Tl);
Tl = q;
T0 = Tl;
[Tl] = ode45(@(x,T) odeTL(x,T,z,n),xspan,T0,opts);
x1 = 5*pi;
T = deval(Tl,x1);
function dTdx = odeTL(x,T,z,n)
Q = [sin(x) 0;0 cos(x)];
V = Q-z*eye(n,n);
W = T+eye(n,n);
R = T-eye(n,n);
Omg = eye(n,n)-0.5*R'*V*W;
dTdx = T*Omg;
end
As I run the code it said 'Matrix dimensions must agree' and I dont really see how the ode45 works for the matrix case?

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

采纳的回答

Alan Stevens
Alan Stevens 2020-11-15
Does this help
z = -1; %parameter Needs to be negative or T1 blows up.
n = 2;
T0 = eye(n,n);
xspan = [0 5*pi];
opts = odeset('RelTol',1e-8,'AbsTol',1e-8);
[x,Tl] = ode15s(@(x,T) odeTL(x,T,z,n),xspan,T0,opts);
plot(x,Tl),grid
function dTdx = odeTL(x,T,z,n)
T = reshape(T,2,2); % T comes in as a column vector
I = eye(n,n);
Q = [sin(x) 0;0 cos(x)];
V = Q-z*I;
W = T+I;
R = T-I;
O = I-0.5*R'*V*W;
dTdx = T*O;
dTdx = dTdx(:); % dTdx must return as a column vector
end
  3 个评论
显示 1更早的评论
Alan Stevens
Alan Stevens 2020-11-15
编辑：Alan Stevens 2020-11-15
Include these lines just before (or just after) the plot command
Tlat5pi = reshape(Tl(end,:),n,n);
disp(Tlat5pi)
The values, of course, will depend on the value you use for z. I chose z = -1 arbitrarily.
Incidentally, I should really have put
T = reshape(T,n,n);
in the loop. However, if you change n to anything other than 2, you will need to specify an equivalent size matrix for Q.
Haya M
Haya M 2020-11-15
yes, I got it! Thank you very much Alan.

请先登录，再进行评论。

更多回答（0 个）

类别

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

标签

产品

Community Treasure Hunt

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

Start Hunting!

Translated by