由于数值分析临时布置了有关于MATLAB一个大作业,才发现身边有相当一部分的同学对matlab相当不熟悉,于是做了一个关于这三个方法的小文章做一点抛砖引玉的作用。
首先是关于这三个方法的总结,欧拉法属于显式的算法,即可以直接由已知的数据推出所需要的数据,但是贴合的程度较低。改进欧拉法增加了一步计算量,使曲线能够更加的贴近于原数据,而四阶龙格法则是贴合程度最高,但相对计算量最大的一个算法。
各个算法的公式如下:
一阶微分方程:
欧拉法:
文章来源地址https://www.toymoban.com/news/detail-778705.html
改进欧拉法:
四阶龙格法:
下面直接上MATLAB代码来实现操作,具体操作的步骤和原因看代码的注释,
%func1.m文件 常微分方程为y' = 2*y/(x+1)
function z = func1(x, y)
z = 2*y/(x+1); %欧拉法函数定义
%forwardEuler1.m文件 向前欧拉法
function [x, y] = forwardEuler(x0, y0, x1, h) %定义函数
n = floor((x1 - x0)/h); %得到所要建立的向量长度
x = zeros(n + 1, 1); %建立x和y的向量
y = zeros(n + 1, 1);
x(1) = x0; %赋值:由于matlab是从1开始而不是从0开始
y(1) = y0;
% 带入公式,用for循环得出向量各个值
for i = 1 : n
x(i + 1) = x(i) + h;
y(i + 1) = y(i) + h * (2 * y(i) / (x(i)+1));
end
%命令1
x0 = 0; %x起始值
x1 = 1; %x终止值
h = 0.1; %步长
y0 = 1; %y起始值
[x, y] = forwardEuler(x0, y0, x1, h); %函数代入
hold on %将多个图均显示出来
plot(x, y,'r','LineWidth', 2) %画出图像,红色,线的宽度为2
l = x0 : 0.1 : x1; %取标准值
lu = (x+1).*(x+1);
plot(l, lu, 'g','LineWidth', 2)
legend('欧拉方法','理论值')
%func2.m文件 常微分方程为y' = 2*y/(x+1)
function z = func2(x, y)
z = 2*y/(x+1); %改进欧拉法函数定义
%Euler2.m文件 未改进的欧拉法
function [x, y] = Euler(x0, x1, y0, h)
n = floor((x1 - x0)/h); % 求n值
x = zeros(n + 1, 1); % 定义x与y向量
y = zeros(n + 1, 1);
x(1, 1) = x0;
y(1, 1) = y0;
for i = 2 : n+1
x(i, 1) = x(i - 1, 1) + h;
y(i, 1) = y(i - 1, 1) + h * (2 * y(i - 1, 1) / (x(i - 1, 1) + 1));
end
%improvedEuler2.m文件 改进后的欧拉法
function [x, y] = improvedEuler(x0, x1, y0, h)
n = floor((x1 - x0)/h);
x = zeros(n + 1, 1);
y = zeros(n + 1, 1);
x(1) = x0;
y(1) = y0;
% 改进欧拉法公式代入,由于有不同的公式则先得出z1,才能进行z2的运算
for i = 1 : n
x(i + 1) = x(i) + h;
z1 = 2 * y(i) / (x(i)+1);
z2 = 2 * (y(i) + h * z1) / (x(i) + h + 1);
y(i + 1) = y(i) + h / 2 * (z1 + z2);
end
%命令2
x0 = 0;
x1 = 1;
y0 = 1;
h = 0.1;
[x, y] = useEuler(x0, x1, y0, h);
yy = (1+x).*(1+x); %精确解
[x, q] = improvedEuler(x0, x1, y0, h);
hold on
plot(x, y, 'k');
plot(x, yy, 'r');
plot(x, q, 'b');
legend('传统欧拉法','理论值','改进欧拉法')
%func3.m文件 常微分方程为y' = 2*y/(x+1)
function z = func3(x, y)
z = 2*y/(x+1); %四阶龙格法函数定义
%runge3.m文件 四阶龙格法方法展示
function [x, y] = runge(x0, x1, y0, h)
n = (x1 - x0) / h;
x = zeros(n + 1);
y = zeros(n + 1);
x(1) = x0;
y(1) = y0;
%按照四阶龙格公式得出下列式子
for i = 1:n
x(i + 1) = x(i) + h;
k1 = func(x(i), y(i));
k2 = func(x(i) + 0.5*h, y(i) + k1*h/2);
k3 = func(x(i) + 0.5*h, y(i) + k2*h/2);
k4 = func(x(i)+ h, y(i) + k3*h);
y(i + 1) = y(i) + h*(k1 + 2*k2 + 2*k3 + k4)/6;
end
end
%命令3
x0 = 0;
x1 = 1;
y0 = 1;
h = 0.1;
[x, y] = runge(x0, x1, y0, h);
yy = (1+x).*(1+x); %精确解
hold on
plot(x, yy, 'r');
plot(x, y, 'b.');
%由于四阶龙格法与理论值的贴合度太高,所以一个采用了线段形式,一个采用了点的形式
%其中点为四阶龙格,线是标准值
其中根据1.2.3分别进行不同的分文件选取,即后缀含有1的代表欧拉法展示,后缀含有2的代表改进欧拉法的展示,后缀带3的代表四阶龙格法展示,大概会产生十个分文件,分别在命令的文件中进行运行即可得到所需的下列三个图像。
文章来源:https://www.toymoban.com/news/detail-778705.html
到了这里,关于matlab分别用欧拉法,改进欧拉法,和四阶龙格法对同一个函数进行拟合的文章就介绍完了。如果您还想了解更多内容,请在右上角搜索TOY模板网以前的文章或继续浏览下面的相关文章,希望大家以后多多支持TOY模板网!