> ansl:=dsolve({ODE,initvals},{u(t),v(t)},type=numeric,method=gear);
ansl := proc(x_gear) ... end proc
> ansl(10,0);
[t???10.,u(t)???.989893921726687442,v(t)???.979787842765888594]
> p1:=plots[odeplot] (ansl,[t,u(t)],0..20,color=red): p2:=plots[odeplot] (ansl,[t,v(t)],0..20,color=blue):
plots[display] ({p1,p2}, title=\
2.5.3 经典数值方法
Maple中常微分方程数值解法中有一类被称作是“经典”(classical)方法. 当然, 称其为经典方法不是因为它们常用或是精度高, 而是因为它们的形式简单, 经常被用于计算方法课上的教学内容. 它们是一些常见的固定步
- 130 -
长方法, 在dsolve中用参数method=classical[方法名称], 如果不特别指出, 将默认采用向前欧拉法. 主要有： foreuler：向前欧拉法(默认) hunform：Heun公式法(梯形方法, 改进欧拉法) imply：改进多项式法 rk2：二阶龙格库塔法 rk3：三阶龙格库塔法 rk4：四阶龙格库塔法 adambash：Adams-Bashford方法(预测法) abmoulton：Adams-Bashford-Moulton方法(预测法) 下面给出微分方程数值方法的参数表： 参数参数类参数用参数用法 名 型 途 initi浮点数指定初 al 的一维值向量 数组 num正整数 指定向 ber 量个数 outp'proced指定生Procedurelis：单个ut urelist'(成单个函数, 返回有序表 默认) 函数或Listprocedure：函数- 131 -
或多个函的有序表 'listproc数的有edure' 序表 proc子程序用子程参数1：未知函数的edur名 序形式个数 e 指定第参数2：自变量 一尖常参数3：函数向量 微分方参数4：导函数向量 程组的右边部分 start 浮点数 自变量 起始值 start布尔量指定数对dverk78不适用 init (默认值积分FALSE) 是否总是从起始值开始 valu浮点数指定需如果给定, 结果是e 向量(一要输出一个2?2的矩阵. 元维数组) 函数值素[1,1]是一个向量, 的自变含自变量名和函数量数值名称; 元素[2,1]是点 一个数值矩阵, 其中第一列value的输入相同, 其他列中是相应的函数值 另外, 还有一些特殊的附加参数： maxfun：整数类型, 用于最大的函数值数量, 默认值50000, 为负数时表示无限制 corrections：正整数类型, 指定每步修正值数量, 在abmoulton中使用, - 132 -
建议值≤4 stepsize：浮点数值, 指定步长 下面看一个简单的例子： > ODE:=diff(y(x),x)=y(x)-2*x/y(x); ?2xODE := y(x)???y(x)??? ?xy(x)> initvals:=y(0)=1; initvals := y(0)???1 >
sol1:=dsolve({ODE,initvals},y(x),numeric,method=classical,stepsize=0.1,start=0);
sol1 := proc(x_classical) ... end proc
而其解析解为： > sol2:=dsolve({diff(y(x),x)=y(x)-2*x/y(x), y(0)=1}, y(x)); sol2 := y(x)???2x???1 将两者图形同时绘制在同一坐标系中比较, 可以发现, 在经过一段时
- 133 -
间后, 欧拉法的数值结果会产生较大误差.
> plot({rhs(sol2),'rhs(sol1(x)[2])'},x=0..2);
求解微方程, 无论使用什么方法或者加入什么选项, 求解完成后必须利用相关数学知识进行逻辑判断, 绝对不对简单迷信Maple给出的结果, 否则很有可能得到一个对于方程本身也许还看得过去, 但在数学或者物理意义上不合理的解.
2.6摄动法求解常微分方程
由于微分方程求解的复杂性, 一般微分方程常常不能求得精确解析解, 需要借助其它方法求得近似解或数值解, 或者两种方法兼而有之. 摄动法是重要的近似求解方法.
- 134 -