用Mathematica求解矩阵微分方程

Mus*_*ush 2 math wolfram-mathematica

我需要在Mathematica中解决这个等式:

d/DX v(x)= .v(x)

这里v列向量 {v1(x),v2(x),v3(x),v4(x)}, A4×4矩阵.我想解决v1, v2, v3, v4任何初始条件下的函数.x的范围是0到1000.

如何使用这种类型的微分方程编写Mathematica代码NDSolve

Sim*_*mon 7

所以,如果你有一些可怕的矩阵

A =  RandomReal[0.1, {4, 4}]; (* A horrible matrix *)
Run Code Online (Sandbox Code Playgroud)

我们做反对称(因此解决方案是振荡的)

A = A - Transpose@A;
Run Code Online (Sandbox Code Playgroud)

定义函数的向量及其初始条件

v[x_] := {v1[x], v2[x], v3[x], v4[x]};

init = v[0] == RandomReal[1, 4]
Run Code Online (Sandbox Code Playgroud)

然后NDSolve命令看起来像

sol = NDSolve[LogicalExpand[v'[x] == A.v[x] && init], 
        {v1, v2, v3, v4}, {x, 0, 1000}]
Run Code Online (Sandbox Code Playgroud)

并且可以绘制解决方案

Plot[Evaluate[v[x] /. sol], {x, 0, 1000}]
Run Code Online (Sandbox Code Playgroud)

da plot


注意,上述微分方程是具有常系数的线性一阶方程,因此简单地使用矩阵指数求解.但是,如果矩阵A是函数x,则分析解决方案变得困难,但数字代码保持不变.

例如,尝试:

A = RandomReal[1/10, {4, 4}] - Exp[-RandomReal[1/100, {4, 4}] x^2];
A = A - Transpose@A;
Run Code Online (Sandbox Code Playgroud)

哪个可以产生像

更多