非线性ode及在实际工程中的应用（第一次大作业）

题目：考察具有粘性阻尼的单摆绘制其在区域中的相图。

 

首先编写一段 MAPLE 子程序 pendphaseplot 如下:

pendphaseplot := proc (c, k, ui, uf, n1, vi, vf, n2, ti, tf, n, vu, vv)

local s, i, a, j, b, sys, fcns, q, p; s := {};

定义变量

for i from 0 to n1 do a := ui+i*(uf-ui)/n1;

a是u变量变化的步长

for j from 0 to n2 do b := vi+j*(vf-vi)/n2;

b是v变量变化的步长

sys := diff(u(t), t) = v(t), diff(v(t), t) = -k*sin(u(t))-c*v(t);

定义系统

fcns := {u(t), v(t)};

定义方程（带入dsolve中用的，就是个过程变量）

q := dsolve({u(0) = a, v(0) = b, sys}, fcns, type = numeric, method = rkf45);

带入不同的初值求得相轨迹

p := plots[odeplot](q, [u(t), v(t)], ti .. tf, numpoints = n, view = [vu, vv], scaling = CONSTRAINED);

将上一步求得的相轨迹作出相应的图像（通过数值解作图）

s := s union p

union是取并集的操作，其实就是用s保存p

end do

end do;

plots[display]([op(s)])

op的含义就是提取整个s用于作图，从这里感觉maple编程是一个比较笨重的语言，当然了，用op也可以提取固定的一些项，详细看帮助即可

end proc

至此整个程序结束

 

命令行：

无阻尼情况：

pendphaseplot(0, 2, -4, 4, 10, -3, 3, 6, -5, 5, 100, -10 .. 10, -5 .. 5)

有阻尼的情况：

pendphaseplot(.5, 2, -4, 4, 8, -3, 3, 6, -5, 5, 100, -10 .. 10, -6 .. 6)

求解过程到这里结束，原始的题目和代码出自胡海岩编写的《应用非线性动力学》，由于****不允许上传课本资料，如果需要可以添加我的微信私下联系：