﻿功能:一阶常微分方程组初值问题求解(4阶龙格库塔法-4阶米尔恩法-汉明法)

格式:
[y1,y2,...,yn,error]=ODESolveHamming(f,fy,t,y0,t0,Error0,n)
[y1,y2,...,yn,error]=ODESolveHamming(f,fy,t,y0,t0,Error0)
[y1,y2,...,yn,error]=ODESolveHamming(f,fy,t,y0,t0)
[y1,y2,...,yn,error]=ODESolveHamming(f,fy,t,y0)
f:用符号变量存储的每个微分表达式,每个表达式以逗号分隔。注意默认的自变量名称为"t"
fy:符号变量存储的对应f当中每个表达式的变量名称
t:为矩阵变量或者数值,表示要求的此时刻的值
y0:为矩阵变量或者数值,表示边界条件值,这个y0对应fy当中的个数
t0:数值,对应t0时的y0,此值默认为0
Error0:变步长控制的相对误差,默认1E-12
n:步长等分最大深度,默认为20

y1,y2,...,yn:对应fy里的变量,表示返回值
error:返回的残差平方和

参考:
[1]王超能.数值分析简明教程(第二版)[M].高等教育出版社,北京,2004:105
[2]现代应用数学手册.计算与数值分析卷[M].清华大学出版社,2007:607

例子:

/*
已知:
y1'=3*y1+2*y2
y2'=4*y1+y2
y1(0)=0
y2(0)=1
求:t=0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1时刻的y1,y2的值
*/

f="3*y1+2*y2,4*y1+y2";
y="y1,y2";
y0=[0 1];
t0=0;
t=0.1:0.1:1;
[y1,y2,er]=ODESolveHamming(f,y,t,y0,t0)//回车后得到如下的结果
y1 =
[ 0.24796127072687    0.63318369140495    1.24695694591787    2.23957868424636    3.85865443480733    6.51224177167353    10.8729555553195    18.0496070374712    29.8701872636631    49.3484265606370 ]
y2 =
[ 1.15279868876949    1.45191444448315    1.98777516660111    2.90989873028219    4.46518509452053    7.06105340776891    11.3695408591110    18.4989360015886    30.2767569234039    49.7163060018085 ]
er =
[ 5.0797269776E-13    6.0854718878E-16    1.0895572055E-13    1.1868359439E-15    2.9911871401E-14    6.0931710620E-13    3.5322628285E-15    5.1678628184E-14    5.4160563150E-13    4.5586247718E-14 ]
