变步长龙格库塔法matlab代码,MATLAB实例源码教程：龙格库塔法求解微分方程组源代码实例.doc...

MATLAB实例源码教程：龙格库塔法求解微分方程组源代码实例

题目：

用经典Runge-Kutta方法求下列一阶微分方程组的近似解：

y1 = 3y1 + 2y2 − (2x2 + 1)e2x, y1(0) = 1 e2x表示exp(2*x)

y2 = 4y1 + y2 + (x2 + 2x − 4)e2x, y2(0) = 1

y3 = 2y1 − y2 − xe3x, y3(0) = 1

y4 = y1 + x2ex, y4(0) = 1

y5 = y2 − e2x, y5(0) = 1

其中初值条件y0为一个五维数组，包含了这五个方程在区间[0,1]左端点0的值，并假设N＝10为区间等分数

程序：

function rk(A,x,h,y0)

%A为字符串函数的元胞数组

%x为x轴上的各点

%h为点距

%y0为初值

i=1;

y(i,:)=y0;

m=length(A);

b=x(length(x));

while x(i)> y1=3*a(2)+2*a(3)-(2*a(1)^2+1)*exp(2*a(1)); 解释：a(1)为x，a(2)、a(3)……为y1，y2……

>> y2=4*a(2)+a(3)+(a(1)^2+2*a(1)-4)*exp(2*a(1));

>> y3=2*a(2)-a(3)-a(1)*exp(3*a(1));

>> y4=a(2)+a(1)^2*exp(a(1));

>> y5=a(3)-exp(2*a(1));

>> A=[{y1},{y2},{y3},{y4},{y5}]

A =

[1x38 char] [1x41 char] [1x28 char] [1x21 char] a(3)-exp(2*a(1))

>> rk(A,0:0.1:1,0.1,ones(1,5))

y =

1.0000 1.0000 1.0000 1.0000 1.0000

1.4692 1.1649 1.1312 1.1226 0.9963

2.1246 1.5112 1.3309 1.3031 0.9930

3.0680 2.1507 1.6105 1.5679 1.0078

4.4629 3.2638 1.9873 1.9573 1.0718

6.5725 5.1403 2.4884 2.5335 1.2376

9.8237 8.2476 3.1587 3.3940 1.5930

14.9117 13.3402 4.0745 4.6925 2.2837

22.9721 21.6384 5.3681 6.6746 3.5497

35.8640 35.1212 7.2699 9.7346 5.7843

56.6365 57.0045 10.1818 14.5095 9.6315