常微分方程数值解法-欧拉方法(1)
课题报告
题目：常微分方程的数值解法-欧拉方法
院 (系)：理学院
专业：数学与信息专业
指导教师：岳宗敏
组员：艾佳欢(组长) 邓云娜
柏茜钟岩刘磊
2015 年 5 月 11 日
常微分方程数值解法-欧拉方法
摘要：从常微分方程数值解的基本概念入手，了解最基本的数值解法--欧拉方法。并利用
欧拉方法显式隐式的特点探究如何求解微分方程，以及欧拉方法的误差分析及校正。 关键词：数值解，欧拉方法，误差，校正
ABSTRACT: From the basic concept of numerical solution of ordinary differential equations, and understand the most basic numerical solution of euler method. And by using euler explicitly implicit characteristics and explore how to solve differential equations and the error analysis and correction of euler method.
KEYWORDS ：arithmetic solution,Euler's method,error,revise
1.初值问题数值解基本概念
初值问题的数值解法，是通过微分方程离散化而给出解在某些节点上的近似值。 在[]b a ,上引入节点{}),,1(,:1100n k x x h b x x x a x k k k n n
k k =-==<<<=-=称
为步长。在多数情况下，采用等步长，即),1,0(,n k kh a x n
a
b h k =+=-=
。记准确解为)(x y ,记)(k x y 的近似值为k y ,记),(k k y x f 为k f .一阶常微分方程的初值问题
??
?=∈=')
2.1()()()
1.1(),())(,()(0 x y a y b a x x y x f x y ， 若f 在{}
?+∞≤≤=y b x a D ,内连续，且满足Lip 条件：0≥?L 使
2121),(),(y y L y x f y x f -≤-，则初值问题的连续可微解)(x y 在],[b a 上唯一存在，
称解)(x y 在节点i x 处的近似值)(i i x y y =为其数值解，该方法称为数值方法。
2. Euler 方法
在工程计算中许多实际问题的数学模型都可以用常微分方程来描述。除少常系数线性微分方程和少数特殊的微分方程可用解析方法求解外，大多数常微分方程难以求得其精确解。因此研究常微分方程的数值解法具有重要的意义。本课题研究的是关于常微分方程初值问题的最简单的数值解法，单步法中的一种---Euler 方法。
2.1 显式Euler 方法
设节点为b x x x a n =<<= 10。初值问题的显式Euler 方法为
???-=+==+1,,1,0,1
0n k f h y y a
y k k k k （2.1）
其中
),(,1k k k k k k y x f f x x h =-=+
2.1.1显式Euler 方法的导出
方法1：Taylor 展开法：
将)(1+k x y 在k x x =点进行Taylor 展开,得
[]12
1,,!
2)())(,()()(++∈''+
+=k k k k k k k k k k x x h y x y x f h x y x y ξξ 忽略2
k h 这一阶项,分别用()k k k k k y x f f y y ,,,1=+近似)(k x y ，)(1+k x y 和))(,(k k x y x f ，
得k k k k f h y y +==+1。结合初值条件α=)0(y 即得（2.1）。 方法 2：向前差分近似微分法：
用向前差分
k
k k h x y x y )
()(1-+近似微分)('k x y ，得
))(,()()(1k k k k k x y x f h x y x y ≈-+ 将近似号改作等号，用k k k f y y ,,1+近似)(k x y ，)(1+k x y ))(,(k k x y x f ，并结合初值条件即得（2.1）。
方法3：左矩数值积分法：
将（1.1）两边从k x 到1+k x 积分得
dx x y x f x y x y k k
x x k k ?
+=-+1
))(,()()(1
用k y ，1+k y 近似)(k x y 、)(1+k x y ，数值积分采用左矩公式得),(1k k k k k y x f h y y =-+，从而亦得（2.1）
Euler 方法有几何意义，如下图，式（1.1）,(1.2)的解曲线)(x y 过点),(000y x P ，
且具斜率0f 。从0P 出发以0f 为斜率作直线
段，交1x x =于),(111y x P ，显然0001f h y y +=。式（1.1）过),(111y x P 的解曲线具有斜率1f 从1P 出发以1f 为斜率作直线要交2x x =于),(222y x P ，余类推。这样我们得到一条折线n P P P 10，它在点k P
的右侧具有斜率k f ，与（1.1）过k P 的解曲线相切。我们取折线n P P P 10，作为（1.1）、（1.2）解曲线)(x y y =的近似曲线，所以Euler 方法又称折线
法。
例1取h=0.1，利用Euler 公式求解
??
??
?=≤≤-='
1)0()
10(,2y x y x y y 解：欧拉公式的具体形式为:
()
????
?
?-+=+=+n n n n n n n y x y y y x hf y y n 21.0,1
其中)10,...,1,0(1.0==+=n n nh a x n ，已知1=n y ，由此式可得:
1.11.01200001=+=????
?
?-+=y x y h y y
191818.11.12.01.11.01.1211112=??? ??
-+=???? ?
?-+=y x y h y y ... ... ....
依次计算可得 3y ，4y ，5y ，6y ，7y ，8y ，9y ，10y
2.2 隐式 Euler 方法和梯形方法
若将)(k x y 在1+k x 展开
11()()(++-=k x k k x f h x y x y 、,)(!
!21
))(21k k k h y x y η''"++1+≤≤k k k x x η
忽略2h 项，用1,+k k y y 和),(111+++=k k k y x f f 分别近似)(k x y ，)(1+k x y 及))(,(11++k k x y x f 可以得另一计算公式
1,,1,0),,(111-=+=+++n k y x f h y y k k k k k （2.2） （2．2）式称为隐式Enler 方法。隐式Euler 方法也可以利用向后差分近似微分或用右矩数值求积公式来建立。
隐式 Euler 方法(2.5)给出了1+k y 要满足的方程，要通过解方程才能得到1+k y 。 在显式和稳式Euler 方法中，忽略的项都是2h 项，为了得到更高精确度的方法，我们可将
1211,)(21
))(,()()(+++≤≤''+
+=k k k k k k k k k k x x h y x y x f h x y x y ξξ 12111,)(21
))(,()()(++++≤≤''-+=k k k k k k k k k k x x h y x y x f h x y x y ηη
取平均得
)]()([4
))](,())(,([2)()(2
111k k k k k k k k k k y y h x y x f x y x f h x y x y ηξ''-''+++=+++
当)(x y 三次连续可微时,)()()(k k k h O y y =''-''ηξ。忽略)(3
k h O 项，用1,+k k y y 分别近似
)(),(1+k k x y x y ，得
)],(),([2
111+++++
=k k k k k
k k y x f y x f h y y （2.3）
（2.3） 称为梯形方法。取这个名称的原因是利用梯形求积公式
?
+=1
))(,([2
)(,(k k
x x k k k
x y x f h dx x y x f ξ=-+++x x y x f D h x y x f x k k k ))(,(12))](,(2311 其中x D 表示关于x 的全微分，忽略数值求积余项也可建立（2.3）。
梯形方法也是隐式方法，要通过解（2.3）来得到1+k y 。
显式Euler 方法取),(),,(k k k k k y x f h y x ==φφ隐式Euler 方法取
()()
111,,,+++==k k k k k y x f h y x φφ
梯形方法取
2
1
),(21),,(1+=
=+k k k k k y x f h y x φφ,(1+k x f )1+k y 当),(y x f 在D 上满足基本条件，),(y x f 关于y 的Lipschitz 常数为L 时，只要
)2.2(,1 为例，当,2 )],(),([2 1y x f y x f h y k k k k k +++ 在∞<<∞-y 上关于y 满足Lipschitz 条件，且Lipschitz 常数为,12 从而由压缩不动点定理得方程 )],(),([2 1y x f y x f h y y k k k k k +++ = 有唯一点不动,1+yk 而且从任意) 0(1+k y 出发，迭代 )],,(),([2 ) (11)1(1i k k k k k k i k y x f y x f h y y +-++++ = 1,,1,0-=n i 都收敛到1+k y ，在实际计算中总希望有较好的) 0(1+k y ，用较少的迭代步，取得有足够精度的 1+k y 例2取h=0.1，利用隐式Euler 公式求解 ?? ? ? ?=≤≤-=' 1)0()10(,2y x y x y y 解：本题中的()y x y y x f 2,- =，故欧拉格式为 .0,1,0,)1(),(01==+-=+=+y n hx y h y x hf y y n n n n n n 隐式欧拉格式为 ,0,,1,0,),(011111==+-=+=+++++y n hx hy y y x hf y y n n n n n n n 整理得 ,0,,1,0),(11 011==++= ++y n hx y h y n n n 取h=0.1，计算结果如下表： 2．3 预估 – 校正Euler 方法 在实际计算中，),()(11i k k y x f ++的计算量比较大，往往取)1() (1≥+m y m k 作为1+k y 来用。我们称)(1m k y +为) 0(1+k y 的m 次迭代改进。最常用的方法之一是先用显式Euler 方法所得的1+k y 为)0(1+k y ，再用梯形方法改进一次 ?? ?? ???-=++=+=++++1,,1,0)],,(),([2) ,(1111n k y x f y x f h y y y x f h y y k k k k k k k k k k k x （2.4） 方法（2.4）称为预估-校正Euler 方法，或改进Euler 方法。预估-校正Euler 方法 还可写成 ))],(,(),([2 11k k k k k k k k k k y x f h y x f y x f h y y +++=++ （2.5） 或 2112 2k h k h y y k k k k ++ =+ （2.6） ),(1k k y x f k = ，),(112k h y x f k k k k +=+ 例3取h=0.1，用改进欧拉格式求下列初值问题： ?? ? ? ?=≤≤-=' 1)0()10(,2y x y x y y 并与精确解进行比较。 解：(1)理论计算 改进欧拉方法的具体形式为： ???? ? ?-+=n n n n p y x y h y y 2 ??? ? ??- +=+p n p n c y x y h y y 12 () 2 1c p n y y y += + 具体计算过程如下: ()1.1011.0121.00000=-?+=???? ? ?-?+=y x y y y p 091818.11.12.01.11.0121.010=??? ??-?+=???? ? ?-?+=p p c y x y y y ()()095909.1091818.11.12 121 1=+=+= c p y y y ... ... .... 由此依次可计算出2y ，3y ，4y ，5y ，6y ，7y ，8y ，9y ，10y 其部分结果见下表 (2)改进欧拉格式MATLAB 程序 f u nction [x, y]=maeuler(dyfun,xspan,y0,h) x=xspan(1):h:xspan(2); y(1)=y0; for n=1:(length(x0)-1) k1=feval(dyfun,x(n),y(n)); y(n+1)=y(n)+h*k1; k2=feval(dyfun,x(n+1),y(n+1)); y(n+1)=y(n)+h*(k1+k2)/2; end x=x';y=y'; 在MA TLAB 命令窗口执行： >> clear;dyfun=inline('x+y'); >> [x,y]=maeuler(dyfun, [0, 0.5], 1, 0.1); [x';y'] ans = 0 0.2000 0.4000 0.6000 0.8000 1.0000 1.0000 1.1841 1.3434 1.4860 1.6165 1.7379 >> y1=2.*exp(x)-x-1; y' %精确解 ans = 1.0000 1.1832 1.3416 1.4832 1.6125 1.7321 >> (y1-y)' %误差 ans = 0 0.0001 0.0017 0.0027 0.0040 0.0058 3.报告总结 欧拉方法是解初值问题最简单的数值方法，在这次课题报告中我们用欧拉方法，隐式欧拉方法及改进欧拉式对同一道习题进行了解答，并对比得出了关于这三种方法的一些结论。 3.1稳定性 稳定性在微分方程的数值解法中是一个非常重要的问题。因为微分方程初值问题的数值方法是用差分格式进行计算的，而在差分方程的求解过程中，存在着各种计算误差，这些计算误差(如舍入误差)引起的扰动，在传播过程中，可能会大量积累，对计算结果的准确性将产生影响。这就需要我们讨论稳定性问题。通过以上习题的计算我们给出了以上几种方法的 由此可知，隐式欧拉法的绝对稳定性比同阶的显式法好，且隐式欧拉法是绝对稳定的(无条件稳定的)，对任何0>h 。 3.2收敛性 常微分方程初值问题的求解，是将微分方程转化为差分方程来求解，并用计算值i y 来近似替代()i x y ，这种替代是否合理，还须看分割区间[]i i x x ,1-的长度h 越来越小时，即 01→-=-i i x x h 时，()i i x y y →是否成立。若成立，则称该方法是收敛的，否则称为不收 敛的。 经过上述习题的计算，我们可知欧拉方法，隐式欧拉方法，梯形方法，改进的欧拉方法的精度p 分别是1阶，1阶，2阶，2阶，关于y 满足Lipschitz 条件且0y 是准确的。因而以上几种方法都满足收敛性定理(收敛性定理:若单步法()h y x h y y n n n n ,,1?+=+具有p 阶精度，且增量函数()h y x ,,?关于y 满足 Lipschitz 条件，即 ()()() y y L h y x h y x -≤-???,,,,，且初值准确)所以以上求解常微分方程的方法都是收敛 的。 另外基于对以上几种方法的探讨，我们又分析了它们的收敛速度。 收敛速度反映了一种方法计算问题的快慢，收敛速度快说明计算效率高，很快就能得到精确的结果。收敛速度的快慢与精度p 的值有关，p 的值越大，收敛速度越快；p 的值越 由上表可以看出四种方法中，改进欧拉方法的收敛速度最快，效果最好。 最后通过对以上几种方法的学习探讨我们总结了这几种方法各自的优缺点，具体如下欧拉方法仅是常微分方程的数值解法的一小部分，另外还有龙格-库塔方法，亚当姆斯方法，这些都是很重要的方法。因此在以后的学习中，我们要想更加全面彻底地了解常微分方程就要不断地学习新的知识，不断地充实自己。 P77 31.利用改进欧拉方法计算下列初值问题，并画出近似解的草图：dy + =t = t y y ≤ ≤ ,2 ;5.0 0,3 )0( )1(= ,1 ? dt 代码： %改进欧拉法 function Euler(t0,y0,inv,h) n=round(inv(2)-inv(1))/h; t(1)=t0; y(1)=y0; for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i)); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i))+ fun(t(i+1),y1(i+1))) end plot(t,y,'*r') function y=fun(t,y); y=y+1; 调用：Euler(0,3,[0,2],0.5) 得到解析解：hold on; y=dsolve('Dy=y+1','(y(0)=3)','t'); ezplot(y,[0,2]) 图像： dy y =t - t y ;2.0 t = ≤ )0( 0,5.0 ,4 )2(2= ≤ ? ,2 dt 代码： function Euler1(t0,y0,inv,h) n=round(inv(2)-inv(1))/h; t(1)=t0; y(1)=y0; for i=1:n y1(i+1)=y(i)+h*fun(t(i),y(i)); t(i+1)=t(i)+h; y(i+1)=y(i)+1/2*h*(fun(t(i),y(i))+ fun(t(i+1),y1(i+1))) end plot(t,y,'*r') function y=fun(t,y); y=y^2-4*t; 调用： Euler1(0,0.5,[0,2],0.2) 图像： 求欧拉回路的Fleury算法 一、实验内容： 判断图G是否存在欧拉回路，若存在，输出其中一条欧拉回路。否则，显示无回路。 二、实验环境：vc++ 三、实验过程与结果 1.问题简介：通过图（无向图或有向图）中所有边一次且仅一次行遍所有 顶点的回路称为欧拉回路。具有欧拉回路的图称为欧拉图 2.算法思想（框图）： （1）任取v0∈V(G)，令P0=v0. （2）设Pi=v0e1v1e2…e i v i已经行遍，按下面方法来从E(G)-{e1,e2,…,e i}中选取e i+1： （a）e i+1与v i相关联； （b）除非无别的边可供行遍，否则e i+1不应该为 G i=G-{e1,e2,…,e i}中的桥。 （3）当（2）不能再进行时，算法停止。 可以证明，当算法停止时所得简单回路Pm=v0e1v1e2…e m v m(v m=v0)为G中一条欧拉回路。 3.数据输入： 边数5，点数6 相关联的点1 2 1 3 2 5 4 2 3 2 4 5 4.运行结果： 存在欧拉回路 1，3，2，4，5，2，1 5.分析总结： Fleury算法是求欧拉图的十分有效的算法，在执行过程中需要用到类似于图的深度优先遍历，因为该算法就是需要将已找到的路径不断的扩展下去，直到将所有边扩展进路径。 四、完整源程序 #include 毕业论文 题目：常微分方程的Euler解法 及其程序设计 学院：数学与信息科学学院 专业：数学与应用数学 毕业年限： 2011年6月 学生： 学号： 指导教师： 常微分方程的Euler解法及其程序设计 摘要本文总结了常微分方程的Euler解法，对各种格式给出了误差估计，设计了这些格式的计算程序． 关键词常微分方程；Euler解法；误差分析；程序设计 Euler Method of Ordinary Differential Equation and Its Programming Abstract Euler method of ordinary differential equation is summarized，the error of each format is analyzed and its programming is designed in this paper. Keywords Ordinary differential equation; Euler method; Error analysis; Programming 科学技术中常常需要求解常微分方程的定解问题，这类问题最简单的形式，即为微分方程 (,)dy f x y dx = （1） 的初值问题 00(,),(). dy f x y dx y x y ?=???=? （2） 定理 （存在与唯一性定理）如果方程（1）的右端函数(,)f x y 在闭矩形域 000000:,R x a x x a y b y y b -≤≤+-≤≤+ 上满足如下条件： （1）在R 上连续； （2）在R 上关于变量y 满足利普希茨（Lipschitz ）条件，即存在常数L ，使 对于R 上任何一对点(,)x y 和(,)x y 有不等式： (,)(,)f x y f x y L y y -≤-, 则初值问题（2）在区间00000x h x x h -≤≤+上存在唯一解 00(),()y y x y x y ==, 其中0(,)min(,),max (,)x y R b h a M f x y M ∈==． 根据存在与唯一性定理，只要(,)f x y 关于y 满足Lipschitz 条件 (,)(,)f x y f x y L y y -≤-， 即可保证其解()y y x =存在并唯一． 然而解析方法只能用来求解少数较简单和典型的常微分方程，例如线性常系数微分方程等，对于变系数常微分方程的解析求解就比较困难，而一般的非线性常微分方程就更不用说，因此，在大多数情况下，实际问题中归结出来的微分方程主要靠数值解法求解． 一、实验内容： 判断图G是否存在欧拉回路，若存在，输出其中一条欧拉回路。否则，显示无回路。 二、实验过程与结果 1.问题简介：通过图（无向图或有向图）中所有边一次且仅一次行遍所有顶点的回路称为欧拉回路。 具有欧拉回路的图称为欧拉图 2.算法思想（框图）： （1）任取v0∈V(G)，令P0=v0. （2）设Pi=v0e1v1e2…eivi已经行遍，按下面方法来从E(G)-{e1,e2,…,ei}中选取ei+1： （a）ei+1与vi相关联； （b）除非无别的边可供行遍，否则ei+1不应该为Gi=G-{e1,e2,…,ei}中的桥。 （3）当（2）不能再进行时，算法停止。 可以证明，当算法停止时所得简单回路Pm=v0e1v1e2…emvm(vm=v0)为G中一条欧拉回路。 3.数据输入： 边数5，点数6 相关联的点1 2 1 3 2 5 4 2 3 2 4 5 4.运行结果： 存在欧拉回路1，3，2，4，5，2，1 5.分析总结： Fleury算法是求欧拉图的十分有效的算法，在执行过程中需要用到类似于图的深度优先遍历，因为该算法就是需要将已找到的路径不断的扩展下去，直到将所有边扩展进路径。 三、完整源程序 #include <> #include <> #include <> struct stack { int top , node[81]; } T,F,A; //顶点的堆栈 int M[81][81]; //图的邻接矩阵int n; int degree[81]; bool brigde(int i,int j) { int flag[81],t,s; for(s=1;s<=n;s++) flag[s]=0; if(degree[i]==1) return false; else { M[i][j]=0;M[j][i]=0; =1; [1]=i; flag[i]=1; t=i; while>0) { for(s=1;s<=n;s++) { if(degree[s]>0){ if(M[t][s]==1) if(flag[s]==0) { ++; []=s; flag[s]=1; t=s; break; } } } if(s>n){ ; t=[]; } } for(s=1;s<=n;s++) { if(degree[s]>0) if(flag[s]==0) { M[i][j]=M[i][j]=1; return true; break; } } if(s>n) return false; } } void Fleury(int x) //Fleury算法{ int i,b=0; if<=n+1){ ++;[]=x; for(i=1;i<=n;i++) if(M[x][i]==1) if(brigde(x,i)==false) { b=1; break; } if(b==1) { M[x][i]=M[i][x]=0; degree[x]--; degree[i]--; Fleury(i); } } 第8章 常微分方程边值问题的数值解法 引 言 第7章介绍了求解常微分方程初值问题的常用的数值方法；本章将介绍常微分方程的边值问题的数值方法。 只含边界条件(boundary-value condition)作为定解条件的常微分方程求解问题称为常微分方程的边值问题(boundary-value problem). 为简明起见，我们以二阶边值问题为 则边值问题(8.1.1)有唯一解。 推论 若线性边值问题 ()()()()()(),, (),()y x p x y x q x y x f x a x b y a y b αβ'''=++≤≤?? ==? (8.1.2) 满足 (1) (),()p x q x 和()f x 在[,]a b 上连续； (2) 在[,]a b 上, ()0q x >, 则边值问题(8.1.1)有唯一解。 求边值问题的近似解，有三类基本方法： (1) 差分法(difference method)，也就是用差商代替微分方程及边界条件中的导数，最终化为代数方程求解； (2) 有限元法(finite element method)； (3) 把边值问题转化为初值问题，然后用求初值问题的方法求解。 差分法 8.2.1 一类特殊类型二阶线性常微分方程的边值问题的差分法 设二阶线性常微分方程的边值问题为 (8.2.1)(8.2.2) ()()()(),,(),(), y x q x y x f x a x b y a y b αβ''-=< 生物信息技术0801 徐聪U200812594 #include . 欧拉路径（欧拉回路）相关定义： 若图G中存在这样一条路径，使得它恰通过G中每条边一次,则称该路径为欧拉路径。若该路径是一个圈，则称为欧拉回路。 具有欧拉回路的图称为欧拉图（简称E图）。具有欧拉路径但不具有欧拉回路的图称为半欧拉图。 判断图是否为欧拉图： 无向图为欧拉图，当且仅当为连通图且所有顶点的度为偶数。 无向图为半欧拉图，当且仅当为连通图且除了两个顶点的度为奇数之外，其它所有顶点的度为偶数。 有向图为欧拉图，当且仅当其基图（忽略有向图所有边的方向，得到的无向图称为该有向图的基图）连通，且所有顶点的入度等于出度。 有向图为半欧拉图，当且仅当其基图连通，且有且仅有一个顶点的入度比出度大1、有且仅有一个顶点的入度比出度小1，其它所有顶点的入度等于出度。 欧拉回路（路径）的算法： 有向图： 第一步，判断是否存在欧拉路径（欧拉回路），如果不存在，算法结束。 第二步，如果存在欧拉路径，从入度比出度小1的点开始BFS；如果存在欧拉回路，从任意一点开始。 第三步，设DFS到点u，遍历u的出边e(u,v)。 第四步，如果e未被标记，转到第五步，否则转到第三步。 第五步，标记e(u,v)，DFS点v。 第六步，边e(u,v)入栈。 第七步，完成DFS后，从栈顶顺序输出边构成一个欧拉路径（欧拉回路）。 无向图： 第一步，判断是否存在欧拉路径（欧拉回路），如果不存在，算法结束。 第二步，如果存在欧拉路径，从度为奇数的点开始BFS；如果存在欧拉回路，从任意一点开始。 第三步，设DFS到点u，遍历u的出边e(u,v)。 第四步，如果e未被标记，转到第五步，否则转到第三步。 第五步，标记e(u,v)及它的反向边，DFS点v。 第六步，边e(u,v)入栈。 第七步，完成DFS后，从栈顶顺序输出边构成一个欧拉路径（欧拉回路）。 '. 欧拉回路的求解 左图是一个井田图，由于2、3、5、8、9、12、14、15几个点都是奇数连线，故不存在欧拉回路。而右图增加几条连线后，该图就存在欧拉回路。 a)假设点1和点2 之间的连线消失， b)建立数学模型把右图的拓扑关系（并考虑a中连线消失的因素）表达出来 c)理解fleury算法，并计算一条欧拉迹，使得该欧拉迹从点1出发，经过b中的每一条边， 最终达到点2 d)使用plot命令把该欧拉迹显示出来。这个动画过程可以用一个for循环语句实现，如下。 其中pos是个2x16的矩阵，2行分别代表x/y轴坐标，每一列表示每个点的坐标，共16个点；另外，T是个2xN的矩阵，每一列表示一条边从T(1，i) 点到T(2，i)点。fleury 算法的目的就是要产生这样一个T 矩阵。 for i … draw_arrow(pos(:,T(1,i))',pos(:,T(2,i))',0.5) pause; end 以下是这段动画的其中几个截图 clear all hold off A=zeros(16); for i=1:16 %A(i,i)=1/2; if i+1<=16 && mod(i,4)~=0 A(i,i+1)=1; end if i+4<=16 A(i,i+4)=1; end end A(2,5)=1; A(3,8)=1; A(9,14)=1; A(12,15)=1; A=A+A'; [T3 c3]=fleury3(A); pos3(1:2,1)=0; for i=1:16 if i==1, pos3(1:2,i)=0; else if mod(i-1,4)~=0 pos3(1,i)=pos3(1,i-1)+1; pos3(2,i)=pos3(2,i-1); else pos3(1,i)=pos3(1,i-4); pos3(2,i)=pos3(2,i-4)-1; end end end figure (1), title('3rd Question') for i=1:16 for j=i:16 if A(i,j)==1, plot([pos3(1,i),pos3(1,j)],[pos3(2,i),pos3(2,j)]); hold on, end end end for i=2:size(T3,2) draw_arrow(pos3(:,T3(1,i))',pos3(:,T3(2,i))',0.5) 欧拉（Euler）通路/回路 1、基本概念： （1）定义 欧拉通路(欧拉迹)—通过图中每条边一次且仅一次，并且过每一顶点的通路。 欧拉回路(欧拉闭迹)—通过图中每条边一次且仅一次，并且过每一顶点的回路。 欧拉图—存在欧拉回路的图。欧拉图就是从一顶出发每条边恰通过一次又能回到出发顶点的那种图，即不重复的行遍所有的边再回到出发点。 通路和回路－称v i e1e2…e n v j为一条从v i到v j且长度为n的通路，其中长度是指通路中边的条数．称起点和终点相同的通路为一条回路。 简单图－不含平行边和自回路的图。 混合图－既有有向边，也有无向边的图 平凡图－仅有一个结点的图 完全图－有n个结点的且每对结点都有边相连的无向简单图，称为无向完全图；有n个结点的且每对结点之间都有两条方向相反的边相连的有向简单图为有向完全图。 （2）欧拉图的特征： 无向图 a）G有欧拉通路的充分必要条件为：G 连通，G中只有两个奇度顶点(它们分别是欧拉通路的两个端点)。 b）G有欧拉回路(G为欧拉图)：G连通，G中均为偶度顶点。 有向图 a）D有欧拉通路：D连通，除两个顶点外，其余顶点的入度均等于出度，这两个特殊的顶点中，一个顶点的入度比出度大1，另一个顶点的入度比出度小1。 b）D有欧拉回路(D为欧拉图)：D连通，D中所有顶点的入度等于出度。一个有向图是欧拉图，当且仅当该图所有顶点度数都是0。 2、弗罗莱（Fleury）算法思想－解决欧拉回路 Fleury算法： 任取v0∈V(G)，令P0=v0； 设Pi=v0e1v1e2…ei vi已经行遍，按下面方法从中选取ei+1： （a）ei+1与vi相关联； （b）除非无别的边可供行遍，否则ei+1不应该为Gi=G-{e1,e2, …, ei}中的桥（所谓桥是一条删除后使连通图不再连通的边）； （c）当（b）不能再进行时，算法停止。 可以证明,当算法停止时所得的简单回路Wm=v0e1v1e2….emvm(vm=v0)为G中的一条欧拉回路，复杂度为O(e*e)。 3、欧拉算法C语言描述 如下为算法的图示动态过程： 第8章 常微分方程数值解法 本章主要内容： 1．欧拉法、改进欧拉法. 2．龙格-库塔法。 3．单步法的收敛性与稳定性。 重点、难点 一、微分方程的数值解法 在工程技术或自然科学中，我们会遇到的许多微分方程的问题，而我们只能对其中具有较简单形式的微分方程才能够求出它们的精确解。对于大量的微分方程问题我们需要考虑求它们的满足一定精度要求的近似解的方法，称为微分方程的数值解法。本章我们主要 讨论常微分方程初值问题?????==00 )() ,(y x y y x f dx dy 的数值解法。 数值解法的基本思想是：在常微分方程初值问题解的存在区间[a,b]内，取n+1个节点a=x 0＜x 1＜…＜x N =b （其中差h n = x n –x n-1称为步长，一般取h 为常数，即等步长），在这些节点上把常微分方程的初值问题离散化为差分方程的相应问题，再求出这些点的上的差分方程值作为相应的微分方程的近似值（满足精度要求）。 二、欧拉法与改进欧拉法 欧拉法与改进欧拉法是用数值积分方法对微分方程进行离散化的一种方法。 将常微分方程),(y x f y ='变为() *+=?++1 1))(,()()(n x n x n n dt t y t f x y x y 1．欧拉法（欧拉折线法） 欧拉法是求解常微分方程初值问题的一种最简单的数值解法。 欧拉法的基本思想：用左矩阵公式计算（＊）式右端积分，则得欧拉法的计算公式为：N a b h N n y x hf y y n n n n -= -=+=+)1,...,1,0(),(1 欧拉法局部截断误差 11121 )(2 ++++≤≤''=n n n n n x x y h R ξξ或简记为O （h 2）。 常微分方程欧拉算法 Company Document number：WUUT-WUUY-WBBGB-BWYTT-1982GT 常微分方程欧拉算法 摘要：本文主要论述了常微分方程的欧拉算法的算法原理，误差分析，实例，程序，以及算法比较等内容。 关键词：常微分方程 显式欧拉法 隐式欧拉法 引言：微分方程初值问题模型是常见的一类数学模型。对于一些简单而典型的微分方程模型，譬如线性方程、某些特殊的一阶非线性方程等是可以设法求出其解析解的，并有理论上的结果可资利用。但在数学建模中碰到的常微分方程初值问题模型，通常很难，甚至根本无法求出其解析解，而只能求其近似解。因此，研究其数值方法，以便快速求得数值鳃有其重大意义。 一、欧拉算法原理 对于微分方程初值问题 的解在xy 平面上是一条曲线，称为该微分方程的积分曲线。积分曲线上一点(),x y 的切线斜率等于函数f 在点(),x y 的值，从初始点()000,P x y 出发，向该点的切线方向推进到下一个点()111,P x y ，然后依次做下去，得到后面的未知点。一般地，若知道(),n n n P x y 依上述方法推进到点()111,n n n P x y +++，则两点的坐标关系为： 即 这种方法就是欧拉（Euler ）方法（也叫显式欧拉法或向前欧拉法）。当初值0y 已知，则n y 可以逐步算出 对微分方程()=x y dy f dx ，从n x 到1n x +积分，那么有 现在用左矩形公式()(),n n hf x y x 代替()()1 ,n n x x f t y t dt +?，n y 代替()n y x ，1n y +代替() 1n y x +就得到了欧拉方法。如果用右矩形公式()()11,n n hf x y x ++去代替右端积分，则得到另外一 个公式，该方法就称为隐式欧拉法（或后退欧拉法），其公式为 欧拉公式与隐式欧拉公式的区别在于欧拉公式是关于1n y +的一个直接计算公式，然而隐式欧拉公式右端含有1n y +，所以它实际上是关于1n y +的一个函数方程。 二、实例 例 取h=，用Euler 方法解 利用MATLAB求解常微分方程数值解 目录 1. 内容简介 (1) 2. Euler Method(欧拉法)求解 (1) 2.1. 显式Euler法和隐式Euler法 (2) 2.2. 梯形公式和改进Euler法 (3) 2.3. Euler法实用性 (4) 3. Runge-Kutta Method(龙格库塔法)求解 (5) 3.1. Runge-Kutta基本原理 (5) 3.2. MATLAB中使用Runge-Kutta法的函数 (7) 4. 使用MATLAB求解常微分方程 (7) 4.1. 使用ode45函数求解非刚性常微分方程 (8) 4.2. 刚性常微分方程 (9) 5. 总结 (9) 参考文献 (11) 附录 (12) 1. 显式Euler法数值求解 (12) 2. 改进Euler法数值求解 (12) 3. 四阶四级Runge-Kutta法数值求解 (13) 4.使用ode45求解 (14) 1.内容简介 把《高等工程数学》看了一遍，增加对数学内容的了解，对其中数值解法比较感兴趣，这大概是因为在其它各方面的学习和研究中经常会遇到数值解法的问题。理解模型然后列出微分方程，却对着方程无从下手，无法得出精确结果实在是让人难受的一件事情。 实际问题中更多遇到的是利用数值法求解偏微分方程问题，但考虑到先从常微分方程下手更为简单有效率，所以本文只研究常微分方程的数值解法。把一个工程实际问题弄出精确结果远比弄清楚各种细枝末节更有意思，因此文章中不追求非常严格地证明，而是偏向如何利用工具实际求解出常微分方程的数值解，力求将课程上所学的知识真正地运用到实际方程的求解中去，在以后遇到微分方程的时候能够熟练运用MATLAB得到能够在工程上运用的结果。 文中求解过程中用到MATLAB进行数值求解，主要目的是弄清楚各个函数本质上是如何对常微分方程进行求解的，对各种方法进行MATLAB编程求解，并将求得的数值解与精确解对比，其中源程序在附录中。最后考察MATLAB中各个函数的适用范围，当遇到实际工程问题时能够正确地得到问题的数值解。 2.Euler Method(欧拉法)求解 Euler法求解常微分方程主要包括3种形式，即显式Euler法、隐式Euler法、梯形公式法，本节内容分别介绍这3种方法的具体内容，并在最后对3种方法精度进行对比，讨论Euler法的实用性。 本节考虑实际初值问题 使用解析法，对方程两边同乘以得到下式 i.常微分方程初值问题数值解法 常微分方程初值问题的真解可以看成是从给定初始点出发的一条连续曲线。差分法是常微分方程初值问题的主要数值解法，其目的是得到若干个离散点来逼近这条解曲线。有两个基本途径。一个是用离散点上的差商近似替代微商。另一个是先对微分方程积分得到积分方程，再利用离散点作数值积分。 i.1 常微分方程差分法 考虑常微分方程初值问题：求函数()u t 满足 (,), 0du f t u t T dt =<≤ （i.1a ） 0(0)u u = (i.1b) 其中(,)f t u 是定义在区域G : 0t T ≤≤, u <∞上的连续函数，0u 和T 是给定的常数。我们假设(,)f t u 对u 满足Lipschitz 条件，即存在常数L 使得 121212(,)(,), [0,]; ,(,)f t u f t u L u u t T u u -≤-?∈∈-∞∞ （i.2） 这一条件保证了（i.1）的解是适定的，即存在，唯一，而且连续依赖于初值0u 。 通常情况下，（i.1）的精确解不可能用简单的解析表达式给出，只能求近似解。本章讨论常微分方程最常用的近似数值解法-差分方法。先来讨论最简单的Euler 法。为此，首先将求解区域[0,]T 离散化为若干个离散点： 0110N N t t t t T -=<< <<= （i.3） 其中n t hn =，0h >称为步长。 在微积分课程中我们熟知，微商（即导数）是差商的极限。反过来，差商就是微商的近似。在0t t =处，在（i.1a ）中用向前差商 10()()u t u t h -代替微商du dt ，便得 10000()()(,())u t u t hf t u t ε=++ 如果忽略误差项0ε，再换个记号，用i u 代替()i u t 便得到 1000(,)u u hf t u -= 一般地，我们有 1Euler (,), 0,1, ,1n n n n u u hf t u n N +=+=-方法： （i.4） 从(i.1b) 给出的初始值0u 出发，由上式可以依次算出1,,N t t 上的差分解1,,N u u 。 1. Euler 公式 100(,)() i i i i y y hf x y y y x +=+??=? 实例： ，00(,),0,1,01f x y x y x y x =-==≤≤ 精确解为：1x y x e -=+- 程序代码： DIMENSION x(0:20),y(0:20),z(0:20),k(0:21) DOUBLE PRECISION x,y,z,k,h,x0,y0,z0,k0,n f(x,y)=x-y n=20 h=1/n x(0)=0 y(0)=0 DO i=0,n-1 y(i+1)=y(i)+f(x(i),y(i))*h x(i+1)=x(i)+h ENDDO k(0)=0 DO i=0,n z(i)=k(i)+exp(-k(i))-1 k(i+1)=k(i)+h END DO open(10,file='1.txt') WRITE(10,10) (x(i),y(i),z(i),i=0,20) WRITE(*,10) (x(i),y(i),z(i),i=0,20) 10 FORMAT(1x,f10.8,2x,f10.8,2x,f10.8/) END 输出结果： 0.00000000 0.00000000 0.00000000 0.05000000 0.00000000 0.00122942 0.10000000 0.00250000 0.00483742 0.15000000 0.00737500 0.01070798 0.20000000 0.01450625 0.01873075 0.25000000 0.02378094 0.02880078 ???=='00)(),(y x y y x f y ???=='0 0)(),(y x y y x f y Euler法解常微分方程 Euler法解常微分方程算法： Step 1 分别取积分上限、积分下限、步长 Step 2计算判断是否成立，成立转到Step 3，否则继续进行Step 4 Step 3 计算 Step 4 Euler法解常微分方程算程序： function euler2(fun,y0,A,h) %fun--y' %y0---初值 %A----x取值范围 %a----x左区间端点值 %b----x右区间端点值 %h----给定步长 x=min(A); b=max(A); y=y0; while x Step 3 （1）做显性Euler预测 （2）将带入 Step 4计算判断是否成立，成立返回Step 3，否则继续进行Step 5 Step 5 改进Euler法解常微分方程算程序： function gaijineuler2(fun,y0,A,h) %fun--y' %y0---初值 %A----x取值范围 %a----x左区间端点值 %b----x右区间端点值 %h----给定步长 a=min(A); b=max(A); x=a:h:b; y(1)=y0; for i=1:length(x)-1 w1=feval(fun,x(i),y(i)); y(i+1)=y(i)+h*w1; w2=feval(fun,x(i+1),y(i+1)); y(i+1)=y(i)+h*(w1+w2)/2; end x=x' y=y' 例：用改进Euler法计算下列初值问题（取步长h=0.25） 输入：fun=inline('-x*y^2') gaijineuler2(fun,2,[0 5],0.25) 得到： x = 0.2500 0.5000 0.7500 1.0000 1.2500 1.5000 1.7500 2.0000 2.2500 2.5000 2.7500 《离散数学》实验报告 （2015 / 2016 学年第一学期） 题目：图的随机生成及欧拉（回）路的确定 专业信息安全 学生姓名 班级学号 指导教师 指导单位计算机学院计算机科学与技术系 日期2015年12月29日 图的随机生成及欧拉（回）路的确定 一、实验内容和要求 内容：编程随机生成n个结点的无向图并能进行（半）欧拉图的判定，若是则给出欧拉（回）路。 要求：对给定n个结点，随机生成邻接矩阵以确定某无向简单图并进行欧拉图和半欧拉图的判定，若符合则给出至少一条欧拉回路或欧拉路。 二、实验目的 对给定n个结点，随机生成邻接矩阵以确定某无向简单图并进行欧拉图和半欧拉图的判定，若符合则给出至少一条欧拉回路或欧拉路。 三、实验任务 1、输入结点个数据生成随机图 2、进行（半）欧拉图的判定 3、若是则给出欧拉（回）路。 四、实验内容 #include 常微分方程初值问题数值解法 朱欲辉 （浙江海洋学院数理信息学院, 浙江舟山316004） [摘要]：在常微分方程的课程中讨论的都是对一些典型方程求解析解的方法.然而在生产实 际和科学研究中所遇到的问题往往很复杂, 在很多情况下都不可能给出解的解析表达式. 本篇文章详细介绍了常微分方程初值问题的一些数值方法, 导出了若干种数值方法, 如Euler法、改进的Euler法、Runge－Kutta法以及线性多步法中的Adams显隐式公式和预测校正 公式, 并且对其稳定性及收敛性作了理论分析. 最后给出了数值例子, 分别用不同的方法计算出近似解, 从得出的结果对比各种方法的优缺点. [关键词]：常微分方程;初值问题; 数值方法; 收敛性; 稳定性; 误差估计 Numerical Method for Initial-Value Problems Zhu Yuhui (School of Mathematics, Physics, and Information Science, Zhejiang Ocean University, Zhoushan, Zhejiang 316004) [Abstract]:In the course about ordinary differential equations, the methods for analytic solutions of some typical equations are often discussed. However, in scientific research, the problems are very complex and the analytic solutions about these problems can’t be e xpressed explicitly. In this paper, some numerical methods for the initial-value problems are introduced. these methods include Euler method, improved Euler method, Runge-Kutta method and some linear multistep method (e.g. Adams formula and predicted-corrected formula). The stability and convergence about the methods are presented. Some numerical examples are give to demonstrate the effectiveness and accuracy of theoretical analysis. [Keywords]:Ordinary differential equation; Initial-value problem; Numerical method; Convergence; Stability;Error estimate 常微分方程欧拉算法 Prepared on 22 November 2020 常微分方程欧拉算法 摘要：本文主要论述了常微分方程的欧拉算法的算法原理，误差分析，实例，程序，以及算法比较等内容。 关键词：常微分方程 显式欧拉法 隐式欧拉法 引言：微分方程初值问题模型是常见的一类数学模型。对于一些简单而典型的微分方程模型，譬如线性方程、某些特殊的一阶非线性方程等是可以设法求出其解析解的，并有理论上的结果可资利用。但在数学建模中碰到的常微分方程初值问题模型，通常很难，甚至根本无法求出其解析解，而只能求其近似解。因此，研究其数值方法，以便快速求得数值鳃有其重大意义。 一、欧拉算法原理 对于微分方程初值问题 的解在xy 平面上是一条曲线，称为该微分方程的积分曲线。积分曲线上一点(),x y 的切线斜率等于函数f 在点(),x y 的值，从初始点()000,P x y 出发，向该点的切线方向推进到下一个点()111,P x y ，然后依次做下去，得到后面的未知点。一般地，若知道(),n n n P x y 依上述方法推进到点()111,n n n P x y +++，则两点的坐标关系为： 即 这种方法就是欧拉（Euler ）方法（也叫显式欧拉法或向前欧拉法）。当初值0y 已知，则n y 可以逐步算出 对微分方程()=x y dy f dx ，从n x 到1n x +积分，那么有 现在用左矩形公式()(),n n hf x y x 代替()()1 ,n n x x f t y t dt +?，n y 代替()n y x ，1n y +代替() 1n y x +就得到了欧拉方法。如果用右矩形公式()()11,n n hf x y x ++去代替右端积分，则得到另外一 个公式，该方法就称为隐式欧拉法（或后退欧拉法），其公式为 欧拉公式与隐式欧拉公式的区别在于欧拉公式是关于1n y +的一个直接计算公式，然而隐式欧拉公式右端含有1n y +，所以它实际上是关于1n y +的一个函数方程。 二、实例 例 取h=，用Euler 方法解 第八章 常微分方程数值解法 考核知识点： 欧拉法，改进欧拉法，龙格-库塔法，单步法的收敛性与稳定性。 考核要求： 1． 解欧拉法，改进欧拉法的基本思想；熟练掌握用欧拉法，改进欧拉法、求微 分方程近似解的方法。 2． 了解龙格-库塔法的基本思想；掌握用龙格-库塔法求微分方程近似解的方 法。 3． 了解单步法的收敛性、稳定性与绝对稳定性。 例1 用欧拉法，预估——校正法求一阶微分方程初值问题 ? ??=-='1)0(y y x y ，在0=x （0，1）0.2近似解 解 （1）用1.0=h 欧拉法计算公式 n n n n n n x y y x y y 1.09.0)(1.01+=-+=+，1.0=n 计算得 9.01=y 82.01.01.09.09.02=?+?=y （2）用预估——校正法计算公式 1,0)(05.01.09.0)0(111)0(1=???-+-+=+=++++n y x y x y y x y y n n n n n n n n n 计算得 91.01=y ，83805.02=y 例2 已知一阶初值问题 ???=-='1 )0(5y y y 求使欧拉法绝对稳定的步长n 值。 解 由欧拉法公式 n n n n y h y h y y )51(51-=-=+ n n y h y ~)51(~1-=+ 相减得01)51()51(e h e h e n n n -==-=-Λ 当 151≤-h 时，4.00≤常微分方程作业欧拉法与改进欧拉法
求欧拉回路的Fleury算法
常微分方程的Euler解法
求欧拉回路的Fleury算法
常微分方程边值问题的数值解法
欧拉及改进的欧拉法求解常微分方程
欧拉路径(欧拉回路)算法
欧拉回路的求解matlab
求欧拉回路,Fleury算法的C语言实现[1]
第8章 常微分方程数值解法 本章主要内容： 1.欧拉法
常微分方程欧拉算法
MATLAB求解常微分方程数值解
常微分方程数值解法
fortran下欧拉法求解常微分方程(实例)
MATLABEuler法解常微分方程
图的随机生成及欧拉(回)路的确定
常微分方程初值问题数值解法
常微分方程欧拉算法
常微分方程数值解法