计算方法上机作业

微分方程求解

二、目的和意义

在科学和工程技术中，很多数学模型是从常微分方程中出现的，因此，常微分方程求解在科学计算中占有重要地位。常微分方程求解主要有两大类方法，即单步方法和线性多步方法。

单步方法容易编写程序，且收敛性好，缺点是计算结果的精度依赖于调用函数f(x,y)的次数。线性多步方法具有计算函数值次数少，结果精度高的特点，并且容易估计误差。

二、计算公式 1.改进欧拉方法

2.经典四级四阶Rung-Kutta方法

四、结构程序设计代码

2.经典四级四阶Rung-Kutta方法

(1).范德堡方程函数

function F=vandepol(t,y)

global mu; % 全局变量

F=[y(2);mu*(1-y(1)^2)*y(2)-y(1)]; end

(2).主函数

global mu; % 定义全局变量，以实现参数?在MATLAB的基本工作空间和函数的专用 % 空间之间数据的传递 mu=1.0;

tspan=[0,30]; % 积分区间 y0=[0;1]; % 初值向量

options=odeset('RelTol',1e-03,'AbsTol',[1e-06,1e-06]); %设置参数 % RelTol 相对容许误差，是个标量；

% AbsTol 绝对容许误差，标量或矢量。若是矢量则指定y的每个分量(y1,y2)的 % % 绝对容差；若为标量则指定所有分量的绝对容差

[t Y]=ode45('vandepol',tspan,y0,options); % vandepol是自编ode函数名 % 输出矩阵Y的第一列是范德堡方程的解y(t).

subplot(2,1,1); % 将图形窗口分成 2x1 的矩阵，第1窗口画y1(t) plot(t,Y(:,1),'-o'); % Y的第1列是解y(t) xlabel('Time (second)','color','r');

title(['Van de Pol equation for mu=',num2str(mu,5)],'color','r', 'FontSize',12); % 加标题并在title中插入变量mu的值,并设置颜色，字号 hold on % 添加新的图，将两张图形画在一起

plot(t,Y(:,2),'-r'); % Y的第2列是y2，画y2(t) xlabel('Time (second)','color','r');

title(['Van de Pol equation for mu=',num2str(mu,5)],'color','r', 'FontSize',12);

legend('y(t)','y''(t)'); %加图例：y(t)和y’(t), 注意单引号「'」重覆使用代表单引号’

hold off % 关闭图形添加功能

subplot(2,1,2); % 第2个窗口画平面相图 plot(Y(:,1),Y(:,2)); % 画相图 xlabel('y1','color','r'); ylabel('y2','color','r');

title('相平面图','color','r','FontSize',12);

四、结果及其讨论

2.经典四级四阶Rung-Kutta方法

对于第二种情况，要将mu和tspan分别改为1000和[0,3000]，得出结果：

3.讨论

梯形方法比欧拉方法精度高，但梯形方法是一种隐式方法，每一步计算都要进行迭代，计算量较大，故在实际中将二者结合起来，形成改进的欧拉公式。改进的欧拉方法比欧拉方法精度高。

一般讲，p级R-K方法所达到的最高阶却不一定是p阶。四级以下的R-K方法其最高阶数与计算f的次数一致。对一般p阶R-K公式中，当p>4时随计算f的次数增加，方法的阶数并不一定增加，四级四阶R-K公式是应用最广泛的公式。