欧拉格式 改进欧拉格式 后退欧拉格式
MATLAB改进欧拉法与四阶龙格-库塔求解一阶常微分方程
姓名:樊元君学号:02 日期:一、实验目的掌握MATLAB语言、C/C++语言编写计算程序的方法、掌握改进欧拉法与四阶龙格-库塔求解一阶常微分方程的初值问题。
掌握使用MATLAB程序求解常微分方程问题的方法。
:二、实验内容1、分别写出改进欧拉法与四阶龙格-库塔求解的算法,编写程序上机调试出结果,要求所编程序适用于任何一阶常微分方程的数值解问题,即能解决这一类问题,而不是某一个问题。
实验中以下列数据验证程序的正确性。
求,步长h=。
*2、实验注意事项的精确解为,通过调整步长,观察结果的精度的变化^)三、程序流程图:●改进欧拉格式流程图:~|●四阶龙格库塔流程图:]四、源程序:●改进后欧拉格式程序源代码:function [] = GJOL(h,x0,y0,X,Y)format longh=input('h=');…x0=input('x0=');y0=input('y0=');disp('输入的范围是:');X=input('X=');Y=input('Y=');n=round((Y-X)/h);\i=1;x1=0;yp=0;yc=0;for i=1:1:nx1=x0+h;yp=y0+h*(-x0*(y0)^2);%yp=y0+h*(y0-2*x0/y0);%·yc=y0+h*(-x1*(yp)^2);%yc=y0+h*(yp-2*x1/yp);%y1=(yp+yc)/2;x0=x1;y0=y1;y=2/(1+x0^2);%y=sqrt(1+2*x0);%fprintf('结果=%.3f,%.8f,%.8f\n',x1,y1,y);:endend●四阶龙格库塔程序源代码:function [] = LGKT(h,x0,y0,X,Y)。
format longh=input('h=');x0=input('x0=');y0=input('y0=');disp('输入的范围是:');"X=input('X=');Y=input('Y=');n=round((Y-X)/h);i=1;x1=0;k1=0;k2=0;k3=0;k4=0;for i=1:1:n~x1=x0+h;k1=-x0*y0^2;%k1=y0-2*x0/y0;%k2=(-(x0+h/2)*(y0+h/2*k1)^2);%k2=(y0+h/2*k1)-2*(x0+h/2)/(y0+h/2*k1);% k3=(-(x0+h/2)*(y0+h/2*k2)^2);%k3=(y0+h/2*k2)-2*(x0+h/2)/(y0+h/2*k2);% k4=(-(x1)*(y0+h*k3)^2);%k4=(y0+h*k3)-2*(x1)/(y0+h*k3);%…y1=y0+h/6*(k1+2*k2+2*k3+k4);%y1=y0+h/6*(k1+2*k2+2*k3+k4);%x0=x1;y0=y1;y=2/(1+x0^2);%y=sqrt(1+2*x0);%fprintf('结果=%.3f,%.7f,%.7f\n',x1,y1,y);end·end*五、运行结果:改进欧拉格式结果:;}四阶龙格库塔结果:步长分别为:和时,不同结果显示验证了步长减少,对于精度的提高起到很大作用,有效数字位数明显增加。
四阶龙格库塔和向前欧拉方法和隐式欧拉法以及改进欧拉法教材
实验九欧拉方法信息与计算科学金融崔振威201002034031一、实验目的:1、掌握欧拉算法设计及程序实现二、实验内容:1、p364-9.2.4、p386-9.5.6三、实验要求:主程序:欧拉方法(前项):function [x,y]=DEEuler(f,a,b,y0,n);%f:一阶常微分方程的一般表达式的右端函数%a:自变量的取值下限%b:自变量的取值上限%y0:函数的初值%n:积分的步数if nargin<5,n=50;endh=(b-a)/n;x(1)=a;y(1)=y0;for i=1:nx(i+1)=x(i)+h;y(i+1)=y(i)+h*feval(f,x(i),y(i));endformat short欧拉方法(后项):function [x,y]=BAEuler(f,a,b,y0,n);%f:一阶常微分方程的一般表达式的右端函数%a:自变量的取值下限%b:自变量的取值上限%y0:函数的初值%n:积分的步数if nargin<5,n=50;endh=(b-a)/n;x(1)=a;y(1)=y0;for i=1:nx(i+1)=x(i)+h;y1=y(i)+h*feval(f,x(i),y(i));y(i+1)=y(i)+h*feval(f,x(i+1),y1);endformat short梯形算法:function [I,step,h2] = CombineTraprl(f,a,b,eps)%f 被积函数%a,b 积分上下限%eps 精度%I 积分结果%step 积分的子区间数if(nargin ==3)eps=1.0e-4;endn=1;h=(b-a)/2;I1=0;I2=(subs(sym(f),findsym(sym(f)),a)+subs(sym(f),findsym(sym(f)),b))/h;while abs(I2-I1)>epsn=n+1;h=(b-a)/n;I1=I2;I2=0;for i=0:n-1x=a+h*i;x1=x+h;I2=I2+(h/2)*(subs(sym(f),findsym(sym(f)),x)+subs(sym(f),findsym(sym(f)),x1));endendI=I2;step=n;h2=(b-a)/n;改进欧拉方法:function [x,y]=MoEuler(f,a,b,y0,n);%f:一阶常微分方程的一般表达式的右端函数%a:自变量的取值下限%b:自变量的取值上限%y0:函数的初值%n:积分的步数if nargin<5,n=50;endh=(b-a)/n;x(1)=a;y(1)=y0; for i=1:nx(i+1)=x(i)+h;y1=y(i)+h*feval(f,x(i),y(i)); y2=y(i)+h*feval(f,x(i+1),y1); y(i+1)=(y1+y2)/2; endformat short四阶龙格-库塔法:function y = DELGKT4_lungkuta(f, h,a,b,y0,varvec) %f:一阶常微分方程的一般表达式的右端函数 %h:积分步长%a :自变量的取值下限 %b:自变量的取值上限%varvec :常微分方程的变量组 format long; N = (b-a)/h;y = zeros(N+1,1); y(1) = y0; x = a:h:b;var = findsym(f); for i=2:N+1K1 = Funval(f,varvec,[x(i-1) y(i-1)]);%Funval 为程序所需要的函数 K2 = Funval(f,varvec,[x(i-1)+h/2 y(i-1)+K1*h/2]); K3 = Funval(f,varvec,[x(i-1)+h/2 y(i-1)+K2*h/2]); K4 = Funval(f,varvec,[x(i-1)+h y(i-1)+h*K3]); y(i) = y(i-1)+h*(K1+2*K2+2*K3+K4)/6; endformat short;p364-9.2.4欧拉方法(前项):1、22)(,1)0(,22'+-+-==-=-t t e t y y y t y t解:执行20步时:编写函数文件doty.m 如下: function f=doty(x,y); f=x^2-y;在Matlab 命令窗口输入:>> [x1,y1]=DEEuler('doty',0,2,1,20)可得到结果:x1 =0 0.1000 0.2000 0.3000 0.4000 0.5000 0.6000 0.70000.8000 0.9000 1.0000 1.1000 1.2000 1.3000 1.4000 1.50001.6000 1.7000 1.8000 1.90002.0000y1 =1.0000 0.9000 0.8110 0.7339 0.6695 0.6186 0.5817 0.55950.5526 0.5613 0.5862 0.6276 0.6858 0.7612 0.8541 0.96471.0932 1.2399 1.4049 1.5884 1.7906在Matlab 命令窗口输入:>> y3=-exp(-x1)+x1.^2-2*x1+2求得解析解:y3 =1.0000 0.9052 0.8213 0.7492 0.6897 0.6435 0.6112 0.59340.5907 0.6034 0.6321 0.6771 0.7388 0.8175 0.9134 1.02691.1581 1.3073 1.4747 1.6604 1.8647输入:>> plot(x1,y1,'o');hold on>> plot(x1,y3,'*');hold on可得近似值与解析解的图像比较:执行40步时:在Matlab 命令窗口输入:>> [x1,y1]=DEEuler('doty',0,2,1,40)可得到结果:x1 =0 0.0500 0.1000 0.1500 0.2000 0.2500 0.3000 0.3500 0.4000 0.4500 0.5000 0.5500 0.6000 0.6500 0.7000 0.75000.8000 0.8500 0.9000 0.9500 1.0000 1.0500 1.1000 1.15001.2000 1.2500 1.3000 1.3500 1.4000 1.4500 1.5000 1.55001.6000 1.6500 1.7000 1.7500 1.8000 1.8500 1.9000 1.95002.0000y1 =1.0000 0.9500 0.9026 0.8580 0.8162 0.7774 0.7417 0.7091 0.6798 0.6538 0.6312 0.6121 0.5967 0.5848 0.5767 0.5724 0.5719 0.5753 0.5826 0.5940 0.6094 0.6290 0.6526 0.68050.7126 0.7490 0.7897 0.8347 0.8841 0.9379 0.9961 1.05881.1260 1.1977 1.2739 1.3547 1.4401 1.5301 1.6247 1.7240 1.8279在Matlab 命令窗口输入:>> y3=-exp(-x1)+x1.^2-2*x1+2求得解析解:y3 =1.0000 0.9513 0.9052 0.8618 0.8213 0.7837 0.7492 0.7178 0.6897 0.6649 0.6435 0.6256 0.6112 0.6005 0.5934 0.5901 0.5907 0.5951 0.6034 0.6158 0.6321 0.6526 0.6771 0.70590.7388 0.7760 0.8175 0.8633 0.9134 0.9679 1.0269 1.09031.1581 1.2305 1.3073 1.3887 1.4747 1.5653 1.6604 1.7602 1.8647输入:>> plot(x1,y1,'o');hold on>> plot(x1,y3,'*');hold on可得近似值与解析解的图像比较:从上面结果可以看出,执行40步时近似解的值要接近于解析解,误差更小,结果更精确。
第五章:常微分方程数值解法第一节欧拉法
改进欧拉法公式
改进欧拉法公式
拓扑学中欧拉公式的证明
用数学归纳法证明
( 1)当 r= 2时,由表明 1,这两个区域可以想象为以赤道为边界的两个半球面,
赤道上存有两个“顶点” 将赤道分为两条“边界”,即为 r= 2,v= 2,e= 2;于是 r+
v- e= 2,欧拉定理设立.。
( 2)设r= m(m≥ 2)时欧拉定理成立,下面证明 r= m+ 1时欧拉定理也成立。
由表明 2,我们在 r= m+ 1的地图上自由选择一个区域 x ,则 x 必存有与它如此相连的区域 y ,使在换成 x 和 y 之间的唯一一条边界后,地图上只有 m 个区域了;在
换成 x 和 y 之间的边界后,若原该边界两端的顶点现在都还是 3条或 3条以上边界
的顶点,则该顶点留存,同时其他的边界数维持不变;若原该边界一端的或两端的顶
点现在沦为 2条边界的顶点,则换成该顶点 ,该顶点两边的两条边界便沦为一条边界。
于就是 ,在换成 x 和 y之间的唯一一条边界时只有三种情况:
①减少一个区域和一条边界;
②增加一个区域、一个顶点和两条边界;
③减少一个区域、两个顶点和三条边界;
即为在换成 x 和 y 之间的边界时 ,不论何种情况都必定存有“增加的区域数 + 增
加的顶点数 = 增加的边界数”我们将上述过程反过来 (即将 x 和 y之间换成的边界又
照曝光原样图画上 ) ,就又沦为 r= m+ 1的地图了,在这一过程中必然就是“减少的
区域数 + 减少的顶点数 = 减少的边界数”。
因此,若r= m (m≥2)时欧拉定理成立,则 r= m+ 1时欧拉定理也成立.。
由 ( 1)和 ( 2)所述 ,对于任何正整数r≥2,欧拉定理设立。
复变函数中欧拉公式的证明。
欧拉格式 改进欧拉格式 后退欧拉格式
functionf=doty(x,y)
f=y-2*x/y;
在matlab命令窗口输入:
[x,y]=euler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
(2)、后退的Euler格式
e3=[0 0.0005 0.0009 0.0013 0.0017 0.0022 0.0027 0.0033 0.0040 0.0048 0.0058];
plot(x,abs(e1),'*');
hold on
plot(x,abs(e2),'b');
hold on
plot(x,abs(e3),'r-');
5、实验总结
通过编程实现了常微分方程初值问题数值解法中的欧拉方法及其后退、改进的算法,并比较了其数值解与精确解之间的误差。可以看出后退的欧拉方法得到的数值解精确度较差,而改进的欧拉方法得到的结果则相对较好。
6、教师评语及评分
hold on
4、实验结果与分析
(1)、Euler格式
计算常微分方程的结果为:
>> [x,y]=euler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
x y y1
0 1.0000 1.0000
0.1000 1.1000 1.0954
Euler格式和RungeKutta格式.ppt
讨论常微分方程(ordinary differential equation)的定解问题,这类问题
的最简单型是如下一阶方程的初值问题:
y f ( x, y), a x b
y( x0 )
y0 ,
y .
(*)
由于数值解是找精确解 y( x) 的近似值,因此,总假设方程的解 y( x)
可用拟合方法求该组数据 ( xi , yi ) 的近似曲线
➢得到数值解有两个基本途径:
⑴ 把近似解表示成有限个独立函数之和,例如:截断的幂(Taylor)
级数或正交函数展开式中的前几项. 涉及到计算高阶导,尽管可用“数 值微分”技术,但得到的公式太长、太复杂!通常比较适用于手算.
⑵ 离散化方法(也称为差分方法),它提供了用当前节点上的或前几
存在且唯一,并具有充分的光滑性!
定理 设函数f(x,y)在区域 D : a x b, y R 上连续,且在区域D内关于y
满足Lipschitz (李普希兹)条件,即存在正数L,使得对所有的 x [a, b]
和任何 y1, y均2有
f ( x, y1 ) f ( x, y2 ) L y1 y2 ,
ba n
②从初始条件 y( x0 ) y0 出发,按照节点的排序,依次逐个计算 yi 的
值,称为步进法. 一般有两种类型:单步法、多步法.
➢为了考察数值解是否具有使用价值,必须解决的基本问题:
①当步长h取得充分小时,所得的近似解yn能否以足够的精度逼近初值
问题的精确解y(xn).这就是收敛问题。即当h 0时, yn y(xn) ?
再从 ( x1, y1) 出发,以 f ( x1, y1) 为斜率作直线 y y1 f ( x1, y1)( x x1)
偏微分方程数值解例题答案
yyy[y例11110.1[1(101)]0.9,10.1[0.9(10.10.9)]0.9019,1(0.90.9019)0.900952p c y y y ì=-´´+´=ïï=-´´+´=íïï=+=î20.900950.1[0.90095(10.10.90095)]0.80274,0.900950.1[0.80274(10.20.80274)]0.80779,1(0.802740.80779)0.805262p c y y yì=-´´+´=ïï=-´´+´=íïï=+=î 这样继续计算下去,其结果列于表9.1. 表9.1 Euler 方法方法改进的Euler 方法方法准确值准确值n xn yny)(n x y0.1 0.9000000 0.9009500 0.9006235 0.2 0.8019000 0.8052632 0.8046311 0.3 0.7088491 0.7153279 0.7144298 0.4 0.6228902 0.6325651 0.6314529 0.5 0.5450815 0.5576153 0.5563460 0.6 0.4757177 0.4905510 0.4891800 0.7 0.4145675 0.4310681 0.4296445 0.8 0.3610801 0.3786397 0.3772045 0.9 0.3145418 0.3326278 0.3312129 1.0 0.2741833 0.2923593 0.2909884 从表9.1可以看出,Euler 方法的计算结果只有2位有效数字,而改进的Euler 方法确有3位有效数字,这表明改进的Euler 方法的精度比Euler 方法高. 例2 试用Euler 方法、改进的Euler 方法及四阶经典R-K 方法在不同步长下计算初值问题ïîïíì=££+-=1)0(,10),1(d d y x xy y xy 在0.2、0.4、0.8、1.0处的近似值,并比较它们的数值结果. 解 对上述三种方法,每执行一步所需计算)1(),(xy y y x f +-=的次数分别为1、2、4。
用改进欧拉法和梯形法解初值问题
在数值分析领域,求解初值问题是一个非常重要的课题。
常见的数值方法有欧拉法和梯形法,但它们都存在一定的局限性。
本文将深入探讨如何改进欧拉法和梯形法,以更高效、精确地解决初值问题。
让我们回顾一下初值问题的基本概念。
初值问题是指,在已知微分方程初值条件的情况下,求解微分方程的数值解。
通常情况下,微分方程很难直接求解,因此需要借助数值方法来逼近微分方程的解。
欧拉法和梯形法是最常见的数值方法之一,它们都是通过离散化微分方程来逼近微分方程的解。
接下来,让我们讨论欧拉法的改进方法。
欧拉法是一种一阶精度的数值方法,它的局限性在于步长选择较大时,数值解会出现较大的误差。
为了改进欧拉法,可以考虑使用改进的欧拉法,如改进的Euler-Cauchy法,它通过使用更复杂的插值公式来提高精度。
也可以考虑使用四阶Runge-Kutta法等更高阶的方法来提高数值解的精度和稳定性。
另让我们分析梯形法的改进方法。
梯形法是一种二阶精度的数值方法,它对微分方程积分过程进行了更加精确的近似。
但梯形法在处理刚性微分方程时会出现数值不稳定的情况。
为了改进梯形法,可以考虑使用隐式的Runge-Kutta法,它可以更好地处理刚性微分方程,提高数值解的稳定性和精度。
总结来说,改进欧拉法和梯形法是为了解决数值解精度不高、数值不稳定等问题。
在实际应用中,根据不同的初值问题特性,选择合适的数值方法很关键。
只有综合考虑数值解精度、稳定性等因素,才能更好地求解初值问题。
从我个人的观点来看,改进欧拉法和梯形法是数值分析领域的一个重要研究方向。
随着科学技术的不断发展,解决初值问题需要更高效、精确的数值方法。
改进欧拉法和梯形法对于提高数值解的精度和稳定性具有重要意义。
在本文中,我们深入探讨了如何改进欧拉法和梯形法,以更高效、精确地解决初值问题。
通过对不同的改进方法进行分析,我们可以更好地理解数值方法的优缺点,为实际问题的求解提供参考。
希望本文能够对初值问题的数值解求解提供一定的帮助,也希望读者可以从中获得一定的启发和思考。
欧拉及改进的欧拉法求解常微分方程
生物信息技术0801徐聪U200812594
#include<stdio.h>
#include<math.h>
void f1(double *y,double *x,double *yy)
{
y[0]=2.0;
x[0]=0.0;
yy[0]=2.0;
f1(y,x,yy);
}
{
y[0]=2.0;
x[0]=0.0;
yy[0]=2.0;
for(int =1;i<=9;i++)
{
x[i]=x[i-1]+0.2;
y0[i]=y[i-1]+0.2*(y[i-1]-x[i-1]);
y[i]=y[i-1]+0.1*(y[i-1]-x[i-1]+y0[i-1]-x[i-1]);
yy[i]=x[i]+1+exp(x[i]);
printf("若x=%f,计算值是%f,真实值是%f,截断误差是%f\n ",x[i],y[i],yy[i],y[i]-yy[i]);
}
};
void f4(double *y,double *x,double *yy,double *y0)
{
y[0]=1.0;
x[0]=0.0;
yy[0]=1.0;
}
};
void f2(double *y,double *x,double *yy)
{
y[0]=1.0;
x[0]=0.0;
yy[0]=1.0;
for(int i=1;i<=9;i++)
第六章_常微分方程初值问题的数值解法_习题课
h2 h3 y ( x n ) y ( x n ) O(h 4 ) 2 6 而且 y ( x n ) f ( x n , y ( x n )) , y ( x n 1 ) f ( x n 1 , y ( x n 1 )) ,对 y ( x n 1 ) 也在 x n 处作 Talor 展开, y ( x n 1 ) y ( x n ) hy ( x n )
湖北民族学院理学院《数值计算方法》教学辅导材料
陈以平编写
h2 h3 y ( x n ) y ( x n ) O(h 4 ) 2 6 h h h2 h3 y ( x n ) y ( x n ) y ( x n ) y ( x n ) y ( x n ) O(h 4 ) 2 2 2 12 h3 y ( x n ) O(h 4 ) O(h 3 ) 12 h3 所以,梯形公式是 2 阶方法,其截断误差的主项是 y ( x n ) 。 12 y ( x n ) hy ( x n )
y k (0.9 0.1y k sin x k ) 0.1( y k 1 y k 1 sin x k 1 )
2
当 k=0,x0=1, y0=1 时,x1=1.2,有 y y (. . y sin x ) (. sin ) .
y f ( x, y ) 3.求解初值问题 欧拉法的局部截断误差是( y ( x ) y 改进欧拉法的局部截断误差是( ); 四阶龙格-库塔法的局部截断误差是( ). (A)O(h2) (B)O(h3) (C)O(h4) (D)O(h5)
4. 改进欧拉法的平均形式公式是( ) y p y k hf ( x k , y k ) y p y k hf ( x k , y k ) (B) y c y k hf ( x k , y p ) .(A) y c y k hf ( x k , y p ) y k ( y p y c ) y k ( y p y c ) y p y k hf ( x k , y k ) y p y k hf ( x k , y k ) (C) y c y k hf ( x k , y p ) (D) y c y k hf ( x k , y p ) y k h ( y p y c ) y k ( y p y c ) (D) 答案:
改进的欧拉公式求微分方程 解释说明以及概述
改进的欧拉公式求微分方程解释说明以及概述1. 引言1.1 概述欧拉公式是数学上一项重要且经典的公式,它将复数、三角函数和指数函数之间建立了一个重要的联系。
自从欧拉在18世纪首次提出这个公式以来,它在各个领域中得到了广泛的应用,尤其是在微分方程的求解过程中起到了关键作用。
1.2 文章结构本文主要围绕改进的欧拉公式在微分方程求解中的应用展开研究。
文章分为五个部分。
首先,引言部分对文章进行了概述,并介绍了文章的结构。
接下来,在第二部分中我们回顾了传统欧拉公式及其在微分方程中的应用情况,并探讨了目前存在的局限性和挑战。
然后,在第三部分中我们详细地介绍了改进的欧拉公式求解微分方程方法,包括新型欧拉公式的推导和定义、改进方法所具有的优势以及针对具体微分方程问题的求解实例和结果分析。
第四部分着重对理论验证与实验结果进行对比与分析,包括理论模型与改进方法之间差异说明、实验设计及数据收集处理方法介绍以及实验结果对比分析与结论得出。
最后,第五部分总结本文的研究工作,并提出了改进欧拉公式研究领域未来发展方向的建议和期待值得探讨的部分。
1.3 目的本文旨在通过改进欧拉公式求解微分方程方法,提高传统欧拉公式在微分方程求解中的应用效果,克服其局限性,并验证改进方法在不同情景下的适用性。
通过理论推导、实验验证和结果分析,将有效地展示改进欧拉公式的优势和改善效果,并从中得出结论,为后续研究提供参考和启示。
同时,文章还将就未来改进欧拉公式研究领域的发展方向进行讨论和展望,为相关领域的学者提供思路和借鉴。
2. 欧拉公式及其应用2.1 欧拉公式的定义欧拉公式,也被称为欧拉恒等式,是数学中一个重要的关系式。
它由数学家莱昂哈德·欧拉在18世纪提出,并以他的名字命名。
这个公式可以被写为:e^(i * π) + 1 = 0其中,e表示自然对数的底数(约等于2.718),i是虚数单位(满足i^2 = -1),π代表圆周率。
2.2 欧拉公式在微分方程中的应用欧拉公式在微分方程领域有着广泛的应用。
_改进欧拉方法范文
_改进欧拉方法范文改进欧拉方法的方法有很多,并且可以从多个方面进行改进。
在下面的文本中,将介绍三种常见的改进方法:改进的欧拉方法、修正的欧拉方法和改进的欧拉方法。
改进的欧拉方法(Improved Euler Method)是一种将欧拉方法进行改进的方法,它通过对函数的斜率进行线性插值来提高计算精度。
具体来说,改进的欧拉方法使用欧拉方法所计算出的斜率与下一个时间步上使用欧拉方法所计算出的斜率的平均值来计算下一个时间步的值。
改进的欧拉方法的迭代公式如下:y(i+1)=y(i)+(1/2)*h*(f(t(i),y(i))+f(t(i+1),y(i)+h*f(t(i),y(i ))))其中,i表示当前时间步,i+1表示下一个时间步,h表示时间步长,t(i)表示当前时间,t(i+1)表示下一个时间,y(i)表示在当前时间步处的函数值,y(i+1)表示下一个时间步处的函数值,f(t,y)表示在时间t处函数的斜率。
修正的欧拉方法(Modified Euler Method)是在改进的欧拉方法的基础上进行改进的方法,它通过对两个时间步使用欧拉方法所计算出的斜率的平均值来计算下一个时间步的值。
修正的欧拉方法的迭代公式如下:y(i+1)=y(i)+h*f(t(i)+1/2*h,y(i)+1/2*h*f(t(i),y(i)))其中,i表示当前时间步,i+1表示下一个时间步,h表示时间步长,t(i)表示当前时间,y(i)表示在当前时间步处的函数值,y(i+1)表示下一个时间步处的函数值,f(t,y)表示在时间t处函数的斜率。
改进的欧拉方法(Heun's method)是一种通过将两个时间步的斜率进行加权平均来计算下一个时间步的值的方法。
改进的欧拉方法的迭代公式如下:y(i+1)=y(i)+(1/2)*h*(f(t(i),y(i))+f(t(i+1),y(i)+h*f(t(i),y(i ))))其中,i表示当前时间步,i+1表示下一个时间步,h表示时间步长,t(i)表示当前时间,t(i+1)表示下一个时间,y(i)表示在当前时间步处的函数值,y(i+1)表示下一个时间步处的函数值,f(t,y)表示在时间t处函数的斜率。
matlab改进欧拉法求解转子运动方程
matlab改进欧拉法求解转子运动方程欧拉法是一种常用的数值计算方法,适用于求解微分方程问题。
在转子运动方程中,我们通常通过欧拉法来模拟转子的运动状态,以了解其运动轨迹和相关参数的变化。
转子运动方程描述了转子在给定的初始条件下的运动状态。
它包含了转子的质量、惯性矩阵、刚度矩阵以及施加在其上的外力。
通过数值方法来求解这些微分方程,我们可以获得转子的角位移、角速度等相关信息。
下面我们将介绍如何使用改进的欧拉法来求解转子的运动方程。
首先,我们需要将转子的运动方程转化为适合数值计算的形式。
通常,在转子运动方程中,我们可以使用角位移和角速度来描述转子的状态。
通过求解转子的运动方程,我们可以得到转子在不同时间点的角位移和角速度。
接下来,我们需要将时间进行离散化,以便在数值计算中使用。
我们可以将整个运动过程划分为许多小的时间步长,然后使用欧拉法逐步求解转子的运动状态。
在标准的欧拉法中,我们使用当前时间点的状态来估计下一个时间点的状态。
但是,这种方法存在精度不高的问题,特别是对于非线性和高阶微分方程。
因此,我们引入改进的欧拉法。
改进的欧拉法使用当前时间点和下一个时间点的状态来进行估计。
具体步骤如下:1. 使用欧拉法计算当前时间点的估计状态。
2. 使用当前时间点的估计状态和下一个时间点的状态进行插值,得到改进的估计状态。
3. 使用改进的估计状态来更新下一个时间点的状态。
通过使用改进的欧拉法,我们可以更准确地估计转子的运动状态,并减小数值计算误差。
为了求解转子的运动方程,我们需要选择合适的时间步长和初始条件。
时间步长应足够小,以保证数值计算的精度,但也不能过小,以充分考虑计算效率。
初始条件应符合实际物理情况,并且应该与欧拉法的初始状态相匹配。
最后,我们使用Matlab编程来实现转子运动方程的求解。
Matlab是一种强大的数值计算工具,提供了许多函数和工具箱,方便我们进行数值计算和可视化。
通过使用改进的欧拉法求解转子运动方程,我们可以更准确地了解转子的运动状态和相关参数的变化。
欧拉公式的改进
假设 yi 1 = y( xi 1 ), yi = y( xi ) ,则可以导出 Ri = y( xi +1 ) yi +1 = O(h3 ) 即中点公式具有 2 阶精度。
§1 Euler’s Method
方 法 显式欧拉 隐式欧拉 梯形公式 中点公式
简单 稳定性最好 精度提高 精度提高, 显式
d f ( x, y) dx 首先希望能确定系数 1、2、p,使得到的算法格式有 2阶 dy 精度,即在 yi = y( xi ) 的前提假设下,使得 = f x ( x, y) + f y ( x, y) dx Ri = y( xi +1 ) yi +1 = O(h3 ) = f x ( x, y) + f y ( x, y) f ( x, y) y( x ) =
§1 Euler’s Method
改进欧拉法 /* modified Euler’s method */
Step 1: 先用显式欧拉公式作预测,算出 y i +1 = y i + h f ( x i , y i ) Step 2: 再将 yi +1 代入隐式梯形公式的右边作校正,得到
y i +1
§2 Runge-Kutta Method Sample Judge Program #include <stdio.h> #include <math.h> #include "98115001_15.h" double f0(double x, double y) { return (yx*x+1.0); } void main( ) { FILE *outfile = fopen("out.txt", "w"); int n; double a, b, y0; a = 0.0; b = 2.0; y0 = 0.50; n = 10; RK4(f0, a, b, y0, n, outfile); fprintf(outfile, "\n"); fclose(outfile); } Sample Output ( represents a space)
计算方法 第6章 常微分方程数值解
已知Euler格式 yn1 yn hf ( xn , yn )
h2 y( xn1 ) yn1 2 y''( xn )
即Euler格式具有一阶精度
如果令
y( xn1 ) y( xn1 ) 2h
y'( xn )
f ( xn , yn )
并假定 y( xn1 ) yn1, y( xn ) yn
常微分方程数值解
常微分方程的数值解法
§1 引 言 §2 欧拉方法 §3 龙格-库塔方法
2
§1 引 言
在工程和科学技术的实际问题中,常需要解常微 分方程。但常微分方程组中往往只有少数较简单和典 型的常微分方程(例如线性常系数常微分方程等)可 求出其解析解。对于变系数常微分方程的解析求解就 比较困难,而一般的非线性常微分方程就更不用说了。 在大多数情况下,常微分方程只能用近似法求解。这 种近似解法可分为两大类:一类是近似解析法,如级 数解法、逐次逼近法等;另一类则是数值解法,它给 出方程在一些离散点上的近似解。
yn
2
xn yn
令 h 0.1 将 x0 0, y0 1 代入Euler格式
步进计算结果见P106表5.1
第五章:常微分方程数值解
Euler值
y 1 2x
第五章:常微分方程数值解
Euler格式的误差分析
pn1
事实上Euler格式的每一步都存在误差,为了方便讨论y算( x)
d2x
dt 2 x(t
0
c
m )
x x
0
0 (t
t) 0
x(t ) x
0
0
5
数学物理方程论文
常微分及椭圆型偏微分方程的数值算法彭毅九江学院理学院 A0821Email:*******************摘要:对于一些不能求解解析解的常微分方程和偏微分方程进行精确求解是非常困难的,本文探讨了应用欧拉法,求解该类常微分方程,通过Matlab的平台,执行欧拉法各步骤。
欧拉法简单地取切线的端点作为下一步的起点进行计算,当步数增多时,误差会因积累而越来越大。
为提高精度,在欧拉格式的基础上进行改进。
采用区间两端的函数值的平均值作为直线方程的斜率,改进的中和欧拉法的精度为二阶。
而在求解偏微分方程数值解的过程中应用最常用的有限元方法求解。
本文以泊松方程为例,在基于变分问题的近似求解中,选择合适的线性基函数,在函数空间中求其弱解,对其区间有限单元化,单元越小(网络越细)则离散域的近似程度越好,计算结果也越精确,但计算量及误差都将增大,对单元构造一个适合的近似解,即推导有限单元的列式,其中包括选择合理的单元坐标系,建立单元基函数,以某种方法给出单元各状态变量的离散关系,从而形成刚度矩阵,将单元总装形成离散域的总矩阵方程,联立方程组求解和结果,再应用数学软件(Matlab)与精确解进行比较,很好的阐述了该方法是一种具有收敛快、精度高、简便有效的通用方法,在工程中具有广阔的应用前景。
关键字:Matlab 欧拉法有限元法初值问题1.引言我们知道微分方程就是联系着自变量、未知函数及其导数的关系式,如果在微分方程中,自变量的个数只有一个,我们称这种微分方程为常微分方程;而含有未知函数偏导数的等式叫偏微分方程。
一般情况下,一个偏微分方程可以写成:0),,,,,,,,,,(= yy xy xx y x u u u u u u y x f 其中,f 是自变量x ,y , 和未知函数u 及其偏导数 ,,,,,yy xy xx y x u u u u u 的已知函数。
解空间的有限维子空间N V 通常由在每一个单元上是自变量x 的多项式,在整个区间[]b a ,上连续,在a x =时取值为零的全体函数所构成,我们称N V 为函数空间。
二阶欧拉格式时间离散化
二阶欧拉格式时间离散化
二阶欧拉格式是一种用于时间离散化的方法,用于将连续时间系统转换为离散时间系统。
在连续时间系统中,二阶欧拉格式可以用于近似求解微分方程。
这种方法可以用于解决线性时不变系统的稳定性问题。
此外,二阶欧拉格式还可以用于离散PID控制器的设计。
在离散时间系统中,二阶欧拉格式也可以用于状态转移矩阵和多变量系统的时间响应分析。
同时,离散时间调制视频介绍了离散时间条件下通过复指数载波信号和正弦载波信号进行调制的情况,与连续时间信号的调制进行对比,提出二者的异同点。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
>> x=0:0.1:1;
e1=[0 0.0046 0.0086 0.0125 0.0166 0.0209 0.0257 0.0311 0.0373 0.0445 0.0527];
e2=[0 -0.0046 -0.0089 -0.0132 -0.0179 -0.0232 -0.0293 -0.0364 -0.0449 -0.0551 -0.0674];
1.0000 1.7848 1.7321
(2)、后退的Euler格式
计算常微分方程的结果为:
>> [x,y]=backeuler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
x y y1
0 1.0000 1.0000
0.8000 1.6165 1.6125
0.9000 1.6782 1.6733
1.0000 1.7379 1.7321
(4)、三者的误差图
注:‘*’表示Euler格式的数值解与准确解的误差,‘b’表示后退的Euler格式的数值解与准确解的误差,‘r-’表示改进的Euler格式的数值解与准确解的误差。
5、实验总结
通过编程实现了常微分方程初值问题数值解法中的欧拉方法及其后退、改进的算法,并比较了其数值解与精确解之间的误差。可以看出后退的欧拉方法得到的数值解精确度较差,而改进的欧拉方法得到的结果则相对较好。
6、教师评语及评分
(4)、根据结果,绘制Euler格式、后退的Euler格式和改进的Euler格式三者的误差图
根据表1可得:
Euler格式的数值解与准确解之间的误差e1:
e1=y-y1=[0 0.0046 0.0086 0.0125 0.0166 0.0209 0.0257 0.0311 0.0373 0.04450.0527];
1.4832
0.7
1.5803
1.5128
1.5525
1.5492
0.8
1.6498
1.5676
1.6165
1.6125
0.9
1.7178
1.6183
1.6782
1.6733
1.0
1.7848
1.6647
1.7379
1.7321
3、详细设计
(1)、Euler格式
欧拉格式的matlab实现程序(euler.m)如下:
微分方程数值解实验报告
实验序号:1日期:2013年10月15日
班级
10应数A班
姓名
黄达
学号
201005050148
实验名称
Euler格式
实验所用软件及版本
Maltab2008
1、实验目的
进一步理解Euler格式、后退的Euler格式、改进的Euler格式的设计思路和算法流程,培养动手实践能力和分析能力。
function[x,y]=neweuler(fun,x0,xfinal,y0,n)
ifnargin<5,n=10;
end
h=(xfinal-x0)/n;
x(1)=x0;y(1)=y0;
fori=1:n
x(i+1)=x(i)+h;
z0=y(i)+h*feval(fun,x(i),y(i));
y(i+1)=y(i)+h/2*(feval(fun,x(i),y(i))+feval(fun,x(i+1),z0));
end
编写一个doty.m文件,如下:
functionf=doty(x,y)
f=y-2*x/y;
在matlab命令窗口输入:
[x,y]=neweuler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
f=y-2*x/y;
在matlab命令窗口输入:
[x,y]=backeuler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
(3)、改进的Euler格式
改进的Euler格式的matlab实现程序(neweuler.m)如下:
后退的Euler格式的matlab实现程序(backeuler.m)如下:
function[x,y]=backeuler(fun,x0,xfinal,y0,n)
ifnargin<5,n=10;
end
h=(xfinal-x0)/n;
x(1)=x0;y(1)=y0;
fori=1:n
x(i+1)=x(i)+h;
1.1000
1.0909
1.0959
1.0954
0.2
1.1918
1.1743
1.1841
1.1832
0.3
1.2774
1.2517
1.2662
1.2649
0.4
1.3582
1.3237
1.3434
1.3416
0.5
1.4351
1.3910
1.4164
1.4142
0.6
1.5090
1.4540
1.4860
0.1000 1.0909 1.0954
0.2000 1.1743 1.1832
0.3000 1.2517 1.2649
0.4000 1.3237 1.3416
0.5000 1.3910 1.4142
0.6000 1.4540 1.4832
0.7000 1.5128 1.5492
0.8000 1.5676 1.6125
2、实验内容
编写Euler格式、后退的Euler格式、改进的Euler格式的程序代码,用于计算下列常微分方程,将计算结果列于表1,并绘制三者的误差图,给出相应的结论。
准确解 ,步长取 。
表1. Euler格式的数值结果(保留5位有效数字)
Euler格式
后退Euler
格式
改进的Euler格式
准确解
0.1
e3=[0 0.0005 0.0009 0.0013 0.0017 0.0022 0.0027 0.0033 0.0040 0.0048 0.0058];
plot(x,abs(e1),'*');
hold on
plot(x,abs(e2),'b');
hold on
plot(x,abs(e3),'r-');
0.2000 1.1918 1.1832
0.3000 1.2774 1.2649
0.4000 1.3582 1.3416
0.5000 1.4351 1.4142
0.6000 1.5090 1.4832
0.7000 1.5803 1.5492
0.8000 1.6498 1.6125
0.9000 1.7178 1.6733
z0=y(i)+h*feval(fun,x(i),y(i));
fork=1:3
z1=y(i)+h*feval(fun,x(i+1),z0);
ifabs(z1-z0)<1e-3
break;
end
z0=z1;
end
y(i+1)=z1;
end
编写一个doty.m文件,如下:
functionf=doty(x,y)
后退的Euler格式的数值解与准确解之间的误差e2:
e2=y-y1=[0 -0.0046 -0.0089 -0.0132 -0.0179 -0.0232 -0.0293 -0.0364 -0.0449 -0.0551-0.0674];
改进的Euler格式的数值解与准确解之间的误差e2:
e3=y-y1=[0 0.0005 0.0009 0.0013 0.0017 0.0022 0.0027 0.0033 0.0040 0.0048 0.0058].
0.9000 1.6183 1.6733
1.0000 1.6647 1.7321
(3)、改进的Euler格式
计算常微分方程的结果为:
>> [x,y]=neweuler('doty',0,1,1,10);
y1=sqrt(1+2.*x); % 对应的准确解
disp(' x y y1')
disp([x',y',y1'])
function[x,y]=euler(fun,x0,xfinal,y0,n)
ifnargin<5,n=10;
end
h=(xfinal-x0)/n;
x(1)=x0;y(1)=y0;
fori=1:n
x(i+1)=x(i)+h;
y(i+1)=y(i)+h*feval(fun,x(i),y(i));
end
x y y1
0 1.0000 1.0000
0.1000 1.0959 1.0954
0.2000 1.1841 1.1832
0.3000 1.2662 1.2649
0.4000 1.3434 1.3416
0.5000 1.4164 1.4142
0.6000 1.4860 1.4832