实验八一阶常微分方程初值问题Matlab实现

不实名感谢某位大兄弟的贡献

【作业内容】: 用Euler法、改进的欧拉法、四阶龙格-库塔法求如下一阶常微分方程
{y′=y−2xyy(0)=1\left\{\begin{array}{l} {y}^{\prime}={y}-\frac{2 {x}}{{y}} \\ {y}({0})={1} \end{array}\right. {y′=y−y2xy(0)=1

x∈[0,2]\mathbf{x} \in[0,2] x∈[0,2]

初值问题的解（小数点后保留8位小数），取h = 0.05.

【作业要求】：

将各种方法求出的解连同精确解使用同一个表格进行数据比较；
画出精确解（见教材）的图形，并且将其它3种方法得到的近似解画出图形,进行比较，见如上的例子7.1。
解释各种求解方法。

【题目分析及所采用的算法】

按题目要求，需要我们计算一阶常微分方程问题，分别使用Euler法、改进的欧拉法、四阶龙格-库塔法来进行迭代计算。

Euler法：
{yi+1=yi+hf(xi,yi)y0=y(x0),i=0,1,…,n\left\{\begin{array}{l}{y}_{\mathrm{i}+1}={y}_{\mathrm{i}}+{h} {f}\left({x}_{\mathrm{i}}, {y}_{\mathrm{i}}\right) \\ {y}_{0}={y}\left({x}_{0}\right)\end{array}, {i}={0}, {1}, \ldots, {n}\right. {yi+1=yi+hf(xi,yi)y0=y(x0),i=0,1,…,n
hhh为步长。依次从初值y0=y(x0)y_0 = y(x_0)y0=y(x0)开始迭代。
改进的欧拉法(将梯形公式和Euler法结合，精度更高):
yp=yi+hf(xi,yi)yc=yi+hf(xi+1,yp)yi+1=12(yp+yc)i=0,1,2…,n−1\begin{array}{l}{y}_{{p}}={y}_{{i}}+{h f}\left({x}_{{i}}, {y}_{{i}}\right) \\ {y}_{{c}}={y}_{{i}}+{h f}\left({x}_{{i}+{1}}, {y}_{{p}}\right) \\ {y}_{{i}+{1}}=\frac{{1}}{{2}}\left({y}_{{p}}+{y}_{{c}}\right) \\ {i}={0}, {1}, {2} \ldots, {n}-{1}\end{array} yp=yi+hf(xi,yi)yc=yi+hf(xi+1,yp)yi+1=21(yp+yc)i=0,1,2…,n−1
hhh为步长。依次从初值y0=y(x0)y_0 = y(x_0)y0=y(x0)开始迭代。
四阶龙格-库塔法:
{k1=f(xi,yi)k2=f(xi+h2,yi+h2k1)k3=f(xi+h2,yi+h2k2)k4=f(xi+h,yi+hk3)yi+1=yi+h6(k1+2k2+2k3+k4)i=0,1,…,n\left\{\begin{array}{l}{k}_{1}={f}\left({x}_{i}, {y}_{i}\right) \\ {k}_{2}={f}\left({x}_{i}+\frac{{h}}{2}, {y}_{i}+\frac{{h}}{2} {k}_{1}\right) \\ {k}_{3}={f}\left({x}_{i}+\frac{{h}}{2}, {y}_{i}+\frac{{h}}{2} {k}_{2}\right) \\ {k}_{4}={f}\left({x}_{i}+{h}, {y}_{i}+{h} {k}_{3}\right) \\ {y}_{{i}+{1}}={y}_{{i}}+\frac{{h}}{{6}}\left({k}_{{1}}+{2} {k}_{{2}}+{2} {k}_{3}+{k}_{4}\right) \quad {i}={0}, {1}, \ldots, {n}\end{array}\right. ⎩⎪⎪⎪⎪⎨⎪⎪⎪⎪⎧k1=f(xi,yi)k2=f(xi+2h,yi+2hk1)k3=f(xi+2h,yi+2hk2)k4=f(xi+h,yi+hk3)yi+1=yi+6h(k1+2k2+2k3+k4)i=0,1,…,n
hhh为步长。依次从初值y0=y(x0)y_0 = y(x_0)y0=y(x0)开始迭代。

【主要代码（要求程序代码有详细的注释）】

代码含有相关注释。

%Ordinary_differential.m
format long;
% 步长
h = 1/20;
% 右边界
xf = 2;
% 设置步长
x = 0 : h : xf;
% 已知函数
Func = @(x, y)(y - ((2 * x) / y));%Euler
% 初始化y数组
y_E = zeros(1, length(x));
y_E(1) = 1;
% 复现公式 & 迭代计算
for i = 1:(length(x)-1)y_E(i + 1) = y_E(i) + Func(x(i), y_E(i)) * h;
end%Improve Euler
% 初始化y数组
y_IE = zeros(1, length(x));
y_IE(1) = 1;
% 复现公式 & 迭代计算
for i = 1:(length(x) - 1)p = y_IE(i) + h * Func(x(i), y_IE(i));c = y_IE(i) + h * Func(x(i + 1), p);y_IE(i + 1) = (1/2) * (p + c);
end%Rung Kutta
% 初始化y数组
y_RK = zeros(1, length(x));
y_RK(1) = 1;
% 复现公式 & 迭代计算
for i=1:(length(x) - 1) k1 = Func(x(i), y_RK(i));k2 = Func(x(i) + 0.5 * h, y_RK(i) + 0.5 * h * k1);k3 = Func((x(i) + 0.5 * h), (y_RK(i) + 0.5 * h * k2));k4 = Func((x(i) + h), (y_RK(i) + k3 * h));y_RK(i + 1) = y_RK(i) + (1 / 6) * (k1 + 2 * k2 + 2 * k3 + k4) * h;
end%Exact
% 使用ode45计算精确函数
x0 = 0;
y0 = 1;
yspan = [x0 xf];
[x_ode45, y_ode45] = ode45(Func,yspan,y0);% 绘图
subplot(221)
plot(x_ode45, y_ode45, '-');
xlabel('x');
ylabel('y');
legend('Exact');subplot(222)
plot(x, y_E, '-');
xlabel('x');
ylabel('y');
legend('Euler');subplot(223)
plot(x, y_IE, '-');
xlabel('x');
ylabel('y');
legend('Improve Euler');subplot(224)
plot(x, y_RK, '-');
xlabel('x');
ylabel('y');
legend('Rung Kutta');% table
res = [round(transpose(x), 8), round(transpose(y_E), 8), round(transpose(y_IE), 8), round(transpose(y_RK), 8), round(y_ode45, 8)];
table = array2table(res, 'VariableNames', {'x','Euler','Improve Euler', 'Rung Kutta', 'Exact'});
table;

【结果比较与（或）分析】

绘图：

明显的看出来欧拉法并没有其他两种方法得来的曲线光滑，更不精确。

可以从数值上进一步判断其精确程度：

>> tabletable =41×5 tablex        Euler       Improve Euler    Rung Kutta      Exact   ____    __________    _____________    __________    __________0             1              1               1             10.05          1.05     1.04886905      1.04880886    1.048810180.1     1.0977381     1.09556112      1.09544514    1.095445380.15    1.14351536     1.14034459      1.14017546    1.140174580.2    1.18757368     1.18343698        1.183216    1.183216060.25    1.23011131     1.22501749      1.22474493    1.224745260.3    1.27129351     1.26523585      1.26491113    1.264911110.35    1.31126017     1.30421877      1.30384056    1.303840290.4     1.3501313     1.34207454      1.34164088    1.341640950.45    1.38801111      1.3788967      1.37840498    1.378405090.5    1.42499118     1.41476663      1.41421368    1.414213620.55     1.4611528     1.44975574      1.44913781    1.449137690.6    1.49656893     1.48392711      1.48323985    1.483239910.65    1.53130567     1.51733679      1.51657526    1.516575290.7    1.56542352     1.55003492      1.54919353    1.549193450.75    1.59897836     1.58206653      1.58113904    1.581138970.8    1.63202233     1.61347233      1.61245178    1.612451820.85    1.66460451     1.64428925      1.64316793    1.643167930.9    1.69677156     1.67455096      1.67332034    1.673320260.95    1.72856823      1.7042883      1.70293895    1.702938891    1.76003786     1.73352962      1.73205115    1.732051171.05    1.79122279     1.76230109      1.76068206    1.760682041.1    1.82216476     1.79062697      1.78885479    1.788854711.15    1.85290524     1.81852981      1.81659066     1.81659061.2     1.8834858     1.84603071      1.84390938    1.843909391.25    1.91394844     1.87314942      1.87082923    1.870829191.3    1.94433584     1.89990456      1.89736719    1.897367091.35    1.97469176     1.92631373      1.92353905    1.923538981.4    2.00506125     1.95239363      1.94935958    1.949359561.45    2.03549101      1.9781602      1.97484254    1.974842471.5    2.06602967      2.0036287      2.00000085    2.000000731.55    2.09672813      2.0288138       2.0248466     2.02484651.6    2.12763984     2.05372973      2.04939117    2.049391111.65    2.15882113     2.07839027      2.07364525    2.073645151.7    2.19033159     2.10280892      2.09761891    2.097618771.75    2.22223435     2.12699893      2.12132168    2.121321541.8     2.2545965     2.15097339      2.14476252    2.144762411.85    2.28748942     2.17474529      2.16794994    2.167949781.9     2.3209892     2.19832764      2.19089198    2.190891781.95    2.35517701     2.22173348      2.21359628    2.213596082    2.39013954     2.24497601      2.23607008     2.2360699>>

可以看出来。精确程度：欧拉法 < 改进的欧拉法 < 龙格库塔法。

实验八一阶常微分方程初值问题Matlab实现相关推荐

实验八.方程根的MATLAB求解
实验内容实验步骤.过程 1.(1)二分法第一步:在 MATLAB 软件,建立一个实现二分法的 MATLAB 函数文件su1.m如下: function x=su1(fname,a,b,e) ...
matlab欧拉法截断误差,一阶常微分方程欧拉法与梯形公式局部截断误差与p阶精度Range.PPT...
一阶常微分方程欧拉法与梯形公式局部截断误差与p阶精度Range * 一阶常微分方程欧拉法与梯形公式局部截断误差与p阶精度 Range-Kutta公式常微分方程MATLAB求解 <数值分析& ...
计算物理学（数值分析）上机实验答案5、常微分方程初值问题的数值解法
实验五.常微分方程初值问题的数值解法常微分方程的求解问题在实践中经常遇到,因此研究常微分方程的数值解法就很有必要.欧拉方法是最简单.最基本的方法,利用差商代替微商,就可得到一系列欧拉公式.这些公 ...
计算方法（六）：常微分方程初值问题的数值解法
文章目录常微分方程初值问题的数值解法欧拉(Euler)方法与改进欧拉方法欧拉方法欧拉公式的局部截断误差与精度分析改进欧拉方法龙格-库塔(Runge-Kutta)法构造原理经典龙格-库塔 ...
一阶欧拉近似matlab,MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程.doc
MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程姓名:樊元君学号:2012200902 日期:2012.11.06 一.实验目的掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题 ...
matlab初值微分方程,常微分方程初值问题的MATLAB解法
大连圣亚海洋世界官网-2021年2月7日发(作者:转身之后还是你) 用 Matlab 求常微分方程 < br>(ODE) 的初值问题 (IVP) 本节考虑一阶常微分方程  u   f ...
matlab实验8数据分析与多项式计算,hashidamatlab实验八数据处理与多项式计算.doc
实验八电子二班张秀云一.实验目的 [据处理与多项式计算 1.掌握数据统计和分析的方法 2.掌握数值插值与曲线拟合的方法及其应用 3.掌握多项式的常用运算二.实验内容 1.利用MATLAB提供的ra ...
【数值分析】用matlab解决插值问题、常微分方程初值问题
一. 二. 三. 四. 一. 问题描述:已知 sin(0.32)=0.314567, sin(0.34)=0.333487, sin(0.36)=0.352274, sin(0.38)=0.37092 ...
matlab实验八,matlab实验八
实验八数据处理与多项式计算一.实验目的 1. 掌握数据统计和分析的方法. 2. 掌握数值插值与曲线拟合的方法及其应用. 3. 掌握多项式的常用运算. 二.实验内容 1. 利用MATLAB提供的ran ...

实验八一阶常微分方程初值问题Matlab实现

实验八一阶常微分方程初值问题Matlab实现相关推荐

最新文章

热门文章

实验八 一阶常微分方程初值问题Matlab实现

实验八 一阶常微分方程初值问题Matlab实现相关推荐

最新文章

热门文章

实验八一阶常微分方程初值问题Matlab实现

实验八一阶常微分方程初值问题Matlab实现相关推荐