Matlab中的常微分方程组答案

【问题标题】：system of ordinary differential equation in MatlabMatlab中的常微分方程组
【发布时间】：2016-09-20 08:47:32
【问题描述】：

假设我们有这个哈密顿量：

n = 10;
H = ones(n,n);

密度矩阵为：

Ro = sym('r',[n,n]);%Density matrix

运动方程为：

 H*Ro-Ro*H

上面的运动方程是方程的右边，左边是密度矩阵的时间导数。

如何在没有符号数学工具箱的情况下在 Matlab 中求解运动方程？我需要更改 n 的值。最多可达 100 个。

【问题讨论】：

“求解”是指从一些初始条件x0 将系统集成有限时间T？
是的，我是这个意思。初始条件可以是一个只有 (1,1) 元素不为零且等于一的矩阵。

标签： matlab ode

【解决方案1】：

在您的动态函数中，在向量和矩阵之间进行整形，以便使用 MATLAB 的 standard ode functions，这（据我所知）需要向量输入。请注意，符号工具箱未在此解决方案的任何地方使用。 R 可以是任意大小的 n-by-n，在您机器内存的限制范围内。

function dR = dynfun(R,H)
                                %// R    = n^2-by-1 vector
                                %// H    = n-by-n matrix
    n  = sqrt(length(R));       
    R  = reshape(R,[n,n]);      %// reshape R to n-by-n matrix
    dR = H*R-R*H;              
    dR = dR(:);                 %// reshape dR to n^2-by-1 vector


end

调用 ODE 求解器：

[tout,Rout] = ode45(@(t,R) dynfun(R,H), [0,T], R0(:));

其中T 是最终时间，R0 是n-by-n 初始条件，tout 是输出时间步长，Rout 是求解轨迹。注意由于整形，Rout 将是k-by-n^2，其中k 是时间步数。如果您想随着时间的推移获得实际矩阵，则需要重塑 Rout 的每一行。

【讨论】：