实验09 数值微积分与方程数值解(第6章)
常微分方程的数值解法实验报告
常微分方程的数值解法专业班级：信息软件 姓名：吴中原 学号：120108010002 一、实验目的1、熟悉各种初值问题的算法，编出算法程序；2、明确各种算法的精度与所选步长有密切关系；通过计算更加了解各种 算法的优越性。
二、实验题目1、根据初值问题数值算法，分别选择二个初值问题编程计算；2、试分别取不同步长，考察某节点j x处数值解的误差变化情况； 3、试用不同算法求解某初值问题，结果有何异常； 4、分析各个算法的优缺点。
三、实验原理与理论基础（一） 欧拉法算法设计对常微分方程初始问题(6-1)(6-2)用数值方法求解时，我们总是认为(6-1)、(6-2)的解存在且唯一。
欧拉法是解初值问题的最简单的数值方法。
从(6-2)式由于y (x 0) = y 0已给定，因而可以算出),()('000y x f x y =。
设x 1 = h 充分小，则近似地有：),()(')()(00001y x f x y hx y x y =≈-(6-3)记 ,n ,,i x y y i i 10 )(== 从而我们可以取),(0001y x hf y y ==作为)(1x y 的近似值。
利用1y 及f (x 1, y 1)又可以算出)(2x y 的近似值：),(1112y x hf y y +=一般地，在任意点()h n x n 11+=+处)(x y 的近似值由下式给出),(1n n n n y x hf y y +=+(6-4)这就是欧拉法的计算公式，h 称为步长。
⎪⎩⎪⎨⎧==)( ),(d d 00y x y y x f x y（二）四阶龙格-库塔法算法设计：欧拉公式可以改写为：()111,i i i i y y k k hf x y +=+⎧⎪⎨=⎪⎩，它每一步计算(),f x y 的值一次，截断误差为()2o h 。
改进的欧拉公式可以改写为：()()()11212112,,i i i i i i y y k k k hf x y k hf x h y k +⎧=++⎪⎪=⎨⎪=++⎪⎩，它每一步要计算(),f x y 的值两次，截断误差为()3o h 。
实验报告——常微分方程的数值解法
实验报告实验项目名称常微分方程的数值解法实验室数学实验室所属课程名称微分方程数值解实验类型上机实验实验日期2013年3月11日班级10信息与计算科学学号2010119421姓名叶达伟成绩实验概述：【实验目的及要求】运用不同的数值解法来求解具体问题，并通过具体实例来分析比较各种常微分方程的数值解法的精度，为以后求解一般的常微分方程起到借鉴意义。
【实验原理】各种常微分方程的数值解法的原理，包括Euler法，改进Euler法，梯形法，Runge-Kutta方法，线性多步方法等。
【实验环境】（使用的软硬件）Matlab软件实验内容：【实验方案设计】我们分别运用Euler法，改进Euler法，RK方法和Adams隐式方法对同一问题进行求解，将数值解和解析解画在同一图像中，比较数值解的精度大小，得出结论。
【实验过程】（实验步骤、记录、数据、分析）我们首先来回顾一下原题：对于给定初值问题：1. 求出其解析解并用Matlab画出其图形；2. 采用Euler法取步长为0.5和0.25数值求解（2.16）,并将结果画在同一幅图中，比较两者精度;3. 采用改进Euler法求解（2.16），步长取为0.5；4. 采用四级Runge-Kutta法求解（2.16），步长取为0.5；5. 采用Adams四阶隐格式计算（2.16），初值可由四级Runge-Kutta格式确定。
下面，我们分五个步骤来完成这个问题：步骤一，求出（2.16）式的解析解并用Matlab 画出其图形； ，用Matlab 做出函数在上的图像，见下图：00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015y=exp(1/3 t 3-1.2t)exact solution图一 初值问题的解析解的图像步骤二，采用Euler 法取步长为0.5和0.25数值求解(2.16),并将结果画在同一幅图中，比较两者精度;我们采用Euler 法取步长为0.5和0.25数值求解，并且将数值解与解析解在一个图中呈现，见下图：00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Numerical solution of Euler and exact solutionexact solution h=0.5h=0.25图二 Euler 方法的计算结果与解析解的比较从图像中不难看出，采用Euler 法取步长为0.5和0.25数值求解的误差不尽相同，也就是两种方法的计算精度不同，不妨将两者的绝对误差作图，可以使两种方法的精度更加直观化，见下图：00.51 1.52 2.53 3.54 4.550.511.522.533.5x 1015Absolute error of numerical solution and exact solutionh=0.5h=0.25图三 不同步长的Euler 法的计算结果与解析解的绝对误差的比较 从图像中我们不难看出，步长为0.25的Euler 法比步长为0.5的Euler 法的精度更高。
常微分方程数值解实验
多步法，Gear’s反向
数值积分，精度中等
若ode45失效时，
可尝试使用
ode23s
刚性
一步法，2阶Rosebrock算法，
低精度。
当精度较低时，
计算时间比ode15s短
odefx为显式常微分方程 中的 ，t为求解区间，要获得问题在其他指定点 上的解，则令t=[t0,t1,t2,…]（要求 单调），y0初始条件。
MATLAB 中有几个专门用于求解常微分方程的函数,它们的设计思想基于Runge-Kutta方法,基本设计思想为：从改进的欧拉方法比欧拉方法精度高的缘由着手,如果在区间[ x1, xi+1]多取几个点的斜率值,然后求取平均值,则可以构造出精度更高的计算方法。 这些函数主要包括：ode45、ode23、ode15s、ode113、ode23s、ode23t、ode23tb. 其中最常用的是函数ode45,该函数采用变步长四阶五阶Runge-Kutta法求数值解,并采用自适应变步长的求解方法。ode23采用二阶三阶Runge-Kutta法求数值解,与ode45类似,只是精度低一些。ode15s用来求刚性方程组。
43
4月22日
588
666
28
46
4月23日
693
782
35
55
4月24日
774
863
39
64
4月25日
877
954
42
73
4月26日
988
1093
48
76
4月27日
1114
1255
56
78
4月28日
1199
1275
59
78
4月29日
《数学软件》实验报告-数值微积分与方程数值求解
.附件二：实验项目列表附件三：实验报告（八）系：专业：年级：姓名：学号：实验课程：实验室号：_ 实验设备号：实验时间：指导教师签字：成绩：1. 实验项目名称：数值微积分与方程数值求解2. 实验目的和要求1.掌握利数据统计和分析的方法2.掌握数值插值与曲线拟合的方法及其应用3.掌握多项式的常用运算3. 实验使用的主要仪器设备和软件方正商祺N260微机；MATLAB7. 0或以上版本4. 实验的基本理论和方法（1）sym（x）：定义符号变量（2）det（X）：矩阵行列式的值（3）polyder（P）:多项式的导函数（4）[l,n]=quad(‘fnsme’,a,b,tol,trace):求定积分（5）直接解法：x=A\b（6）矩阵分解求法：[L,U]=lu(A);x=U\(L\b)（7）迭代解法：[x,n]=jacobi(A,b,[0,0,0,0]',1.0e-6)（8）[x,y]=line_solution(A,b)：线性方程组的通解（9）fzero（filename,x0,tol,trace）:单变量非线性方程求解（10）fsolve（filename,x0，option）：非线性方程组的求解（11）[x,fval]=fminbnd(filename,x1,x2,option):求（x1，x2）区间的极小值点x和最小值fval （12）[x,fval]=fminsearch(filename,x0，option )：基于单纯形算法求多元函数极小值点x和最小值fval（13）[t,y]=ode45(filename,tspan,y0):龙格-库塔法求微分方程的数值解（14）subplot （m ，n ，p ）：子图函数 （15）plot （x ，y ）：绘图函数 5. 实验内容与步骤（描述实验中应该做什么事情，如何做等，实验过程中记录发生的现象、中间结果、最终得到的结果，并进行分析说明）（包括：题目，写过程、答案）题目：1. 求函数在指定点的数值导数32162321)(232，，，==x xx x x x x x ffunction dsx=input('请输入x 的值：'); p=6*x^2 >> x=sym('x');>> f=det([x,x.^2,x.^3;1,2.*x,3.*x.^2;0,2,6.*x]) f = 2*x^3>> f=[2,0,0,0]; >> p=polyder(f) p =6 0 0 >> ds请输入x 的值：1 p = 6 >> ds请输入x 的值：2 p = 24 >> ds请输入x 的值：3 p =542. 用数值方法求定积分(1) ⎰++=π202211)2sin(4cos dtt t I 的近似值。
微分方程数值解法实验报告
微分方程数值解法实验报告2篇微分方程数值解法实验报告（一）在实际科学与工程问题中，我们经常会遇到微分方程的求解。
然而，许多微分方程往往没有解析解，这就需要我们利用数值方法来获得近似解。
本篇实验报告将介绍两种常见的微分方程数值解法：欧拉方法和改进的欧拉方法。
一、欧拉方法欧拉方法是最简单的微分方程数值解法之一。
其基本原理为离散化微分方程，将微分方程中的导数用差商代替。
设要求解的微分方程为dy/dx = f(x, y)，步长为h，则可用以下公式进行递推计算：y_{n+1} = y_n +hf(x_n, y_n)二、算法实现为了对欧拉方法进行数值实验，我们以一阶线性常微分方程为例：dy/dx = x - y, y(0) = 1步骤如下：（1）选择合适的步长h和求解区间[a, b]，这里我们取h=0.1，[a, b] = [0, 1]；（2）初始化y_0 = 1；（3）利用欧拉方法递推计算y_{n+1} = y_n + 0.1(x_n - y_n)；（4）重复步骤（3），直到x_n = 1。
三、实验结果与讨论我们通过上述步骤得到了在[0, 1]上的近似解y(x)。
下图展示了欧拉方法求解的结果。
从图中可以看出，欧拉方法得到的近似解与精确解有一定的偏差。
这是因为欧拉方法只是通过递推计算得到的近似解，并没有考虑到更高阶的误差。
所以在需要高精度解时，欧拉方法并不理想。
四、改进的欧拉方法针对欧拉方法的不足，我们可以考虑使用改进的欧拉方法（也称为改进的欧拉-柯西方法）。
它是通过利用前后两个步长欧拉方法得到的结果作为差商的中间项，从而提高了求解精度。
一阶线性常微分方程的改进欧拉方法可以表示为：y_{n+1} = y_n + h(f(x_n, y_n) + f(x_{n+1}, y_n + hf(x_n,y_n)))/2五、算法实现与结果展示将改进的欧拉方法应用于前述的一阶线性常微分方程，我们同样选择h=0.1，[a, b] = [0, 1]。
微分方程数值解实验报告
微分方程数值解实验报告实验目的：掌握微分方程数值解的基本方法，能够利用计算机编程求解微分方程。
实验原理：微分方程是自然科学与工程技术中常见的数学模型，它描述了变量之间的关系及其随时间、空间的变化规律。
解微分方程是研究和应用微分方程的基础，但有很多微分方程无法找到解析解，只能通过数值方法进行求解。
本实验采用欧拉方法和改进的欧拉方法求解微分方程的初值问题：$$\begin{cases}\frac{dy}{dt}=f(t,y)\\y(t_0)=y_0\end{cases}$$其中，$f(t,y)$是给定的函数，$y(t_0)=y_0$是已知的初值条件。
欧拉方法是最基本的数值解法，其步骤如下：1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$，其中$h$是步长3. 计算$y_{i+1}=y_i+hf(t_i,y_i)$4.重复步骤2、3直到达到终止条件改进的欧拉方法是对欧拉方法进行改进，通过利用函数$y(t)$在$t+\frac{1}{2}h$处的斜率来更准确地估计$y_{i+1}$，其步骤如下：1.给定$t_0$和$y_0$2.计算$t_{i+1}=t_i+h$，其中$h$是步长3. 计算$y_*=y_i+\frac{1}{2}hf(t_i,y_i)$4. 计算$y_{i+1}=y_i+hf(t_i+\frac{1}{2}h,y_*)$5.重复步骤2、3、4直到达到终止条件实验步骤：1.编写程序实现欧拉方法和改进的欧拉方法2.给定微分方程和初值条件3.设置步长和终止条件4.利用欧拉方法和改进的欧拉方法求解微分方程5.比较不同步长下的数值解与解析解的误差6.绘制误差-步长曲线，分析数值解的精度和收敛性实验结果：以一阶常微分方程$y'=3ty+t$为例，给定初值$y(0)=1$，取步长$h=0.1$进行数值求解。
利用欧拉方法求解微分方程得到的数值解如下：\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Euler}} & \text{误差} \\ \hline0.0&1.000&1.000&0.000\\0.1&1.035&1.030&0.005\\0.2&1.104&1.108&0.004\\0.3&1.212&1.217&0.005\\0.4&1.360&1.364&0.004\\0.5&1.554&1.559&0.005\\0.6&1.805&1.810&0.005\\0.7&2.131&2.136&0.005\\0.8&2.554&2.560&0.006\\0.9&3.102&3.107&0.006\\1.0&3.807&3.812&0.005\\\end{array}利用改进的欧拉方法求解微分方程得到的数值解如下：\begin{array}{cccc}t & y_{\text{exact}} & y_{\text{Improved Euler}} & \text{误差} \\\hline0.0&1.000&1.000&0.000\\0.1&1.035&1.035&0.000\\0.2&1.104&1.103&0.001\\0.3&1.212&1.211&0.001\\0.4&1.360&1.358&0.002\\0.5&1.554&1.552&0.002\\0.6&1.805&1.802&0.003\\0.7&2.131&2.126&0.005\\0.8&2.554&2.545&0.009\\0.9&3.102&3.086&0.015\\1.0&3.807&3.774&0.032\\\end{array}误差-步长曲线如下：实验结论：通过对比欧拉方法和改进的欧拉方法的数值解与解析解的误差，可以发现改进的欧拉方法具有更高的精度和收敛性。
第6章数值微积分与常微分方程求解研究报告
是两个等长的向量：X = (x1，x2，…，xn)，Y = (y1， y2，…，yn)，并且x1<x2<…<xn，积分区间是[x1，xn]。
6.2.3 多重定积分的数值求解实现
定积分的被积函数是一元函数，积分范围是一个区间；而重
6.2.2 定积分的数值求解实现 在MATLAB中可以使用quad或quadl来进行数值积分。 1．自适应辛普生法 MATLAB提供了基于自适应Simpson法的quad函数和自适应Lobatto法的
quadl函数来求定积分。函数的调用格式为
[I,n]=quad(@fname,a,b,tol,trace)
6.2 数 值 积 分
6.2.1 数值积分的原理
求解 定积 分 的数 值 方法 多种 多 样， 如 矩形 （ Rectangular） 法、 梯 形（Trape zia） 法 、辛 普生
(Simpson)• 法等都是经常采用的方法。它们的基本思想都是将整个积分区间 [a，b]分成 n 个子区间
[xi，xi + 1]，i = 1，2，…，n，其中 x1 = a，xn + 1 = b。这样求定积分问题就分解为下面的求和问题 ：
常微分方程初值问题的数值解法，首先要解决是建立求数值 解的递推公式。递推公式通常有两类，一类是计算yi +1时 只用到xi +1、xi和yi，即前一步的值，此类方法称为单步 法，其代表是龙格-库塔（Runge-Kutta）法。另一类是计 算yi+1时，需要前面k步的值，此类方法称为多步法，其 代表是亚当姆斯（Adams）法。这些方法都是基于把一个 连续的定解问题离散化为一个差分方程来求解，是一种步 进式的方法。
实验09 数值微积分与方程数值解(第6章)
实验09 数值微积分与方程数值求解（第6章 MATLAB 数值计算）一、实验目的二、实验内容1. 求函数在指定点的数值导数232()123,1,2,3026x x x f x x xx x==2. 用数值方法求定积分(1) 210I π=⎰的近似值。
程序及运行结果：《数学软件》课内实验王平(2) 2221I dx x π=+⎰程序及运行结果：3. 分别用3种不同的数值方法解线性方程组6525494133422139211x y z u x y z u x y z u x y u +-+=-⎧⎪-+-=⎪⎨++-=⎪⎪-+=⎩ 程序及运行结果：4. 求非齐次线性方程组的通解1234123412342736352249472x x x x x x x x x x x x +++=⎧⎪+++=⎨⎪+++=⎩5. 求代数方程的数值解(1) 3x +sin x -e x =0在x 0=1.5附近的根。
程序及运行结果（提示：要用教材中的函数程序line_solution ）：(2) 在给定的初值x 0=1，y 0=1，z 0=1下，求方程组的数值解。
23sin ln 70321050y x y z x z x y z ⎧++-=⎪+-+=⎨⎪++-=⎩6. 求函数在指定区间的极值(1) 3cos log ()xx x x xf x e ++=在(0,1)内的最小值。
(2) 33212112122(,)2410f x x x x x x x x =+-+在[0,0]附近的最小值点和最小值。
7. 求微分方程的数值解，并绘制解的曲线2250(0)0'(0)0xd y dyy dx dx y y ⎧-+=⎪⎪⎪=⎨⎪=⎪⎪⎩程序及运行结果（注意：参数中不能取0，用足够小的正数代替）：令y 2=y,y 1=y '，将二阶方程转化为一阶方程组：'112'211251(0)0,(0)0y y y x x y y y y ⎧=-⎪⎪=⎨⎪==⎪⎩8. 求微分方程组的数值解，并绘制解的曲线123213312123'''0.51(0)0,(0)1,(0)1y y y y y y y y y y y y =⎧⎪=-⎪⎨=-⎪⎪===⎩程序及运行结果：三、实验提示四、教程：第6章 MATLAB 数值计算（2/2）6.2 数值微积分 p155 6.2.1 数值微分1. 数值差分与差商对任意函数f(x)，假设h>0。
微分方程数值解实验报告
微分方程数值解实验报告微分方程数值解实验报告一、引言微分方程是数学中一类重要的方程，广泛应用于各个科学领域。
在实际问题中，往往难以得到微分方程的解析解，因此需要借助数值方法来求解。
本实验旨在通过数值解法，探索微分方程的数值解及其应用。
二、数值解法介绍常用的微分方程数值解法有欧拉法、改进欧拉法、四阶龙格-库塔法等。
在本实验中，我们将采用改进欧拉法进行数值解的求取。
改进欧拉法是一种一阶的显式迭代法，其基本思想是将微分方程的导数用差商来近似表示，并通过迭代逼近真实解。
具体迭代公式如下：\[y_{n+1} = y_n + h \cdot f(x_n + \frac{h}{2}, y_n + \frac{h}{2} \cdot f(x_n, y_n))\]其中，\(y_n\)表示第n步的近似解，\(h\)表示步长，\(f(x_n, y_n)\)表示微分方程的导数。
三、实验步骤1. 确定微分方程及初始条件在本实验中，我们选择经典的一阶常微分方程：\[y' = -2xy\]并给定初始条件\(y(0) = 1\)。
2. 设置步长和迭代次数为了得到较为准确的数值解，我们需要合理选择步长和迭代次数。
在本实验中，我们将步长设置为0.1，迭代次数为10。
3. 迭代计算数值解根据改进欧拉法的迭代公式，我们可以通过编写计算程序来求解微分方程的数值解。
具体计算过程如下：- 初始化：设定初始条件\(y_0 = 1\)，并给定步长\(h = 0.1\)。
- 迭代计算：使用改进欧拉法的迭代公式，通过循环计算得到\(y_1, y_2, ...,y_{10}\)。
- 输出结果：将计算得到的数值解输出，并进行可视化展示。
四、实验结果与分析通过以上步骤，我们得到了微分方程的数值解\(y_1, y_2, ..., y_{10}\)。
将这些数值解进行可视化展示，可以更直观地观察到解的变化趋势。
根据实验结果，我们可以发现随着迭代次数的增加，数值解逐渐逼近了真实解。
(完整版)MATLAB)课后实验答案[1]
![(完整版)MATLAB)课后实验答案[1]](https://img.taocdn.com/s3/m/f6c08c5902d276a201292ece.png)
1 + e2 (2) z = 1 ln( x + 1 + x 2 ) ，其中 x = ⎡⎢ 2⎣-0.45 ⎦2 2 ⎪t 2 - 2t + 1 2 ≤ t <3 ⎨实验一MATLAB 运算基础1. 先求下列表达式的值，然后显示 MATLAB 工作空间的使用情况并保存全部变量。
(1) z = 2sin 8501221 + 2i ⎤5 ⎥(3) z = e 0.3a - e -0.3asin(a + 0.3) + ln 0.3 + a ,a = -3.0, - 2.9, L , 2.9, 3.03⎧t 2 0 ≤ t < 1 (4) z = ⎪t 2 - 11 ≤ t <2 ，其中 t=0:0.5:2.5 4⎩解：M 文件:z1=2*sin(85*pi/180)/(1+exp(2))x=[2 1+2*i;-.45 5];z2=1/2*log(x+sqrt(1+x^2))a=-3.0:0.1:3.0;3=(exp(0.3.*a)-exp(-0.3.*a))./2.*sin(a+0.3)+log((0.3+a)./2)t=0:0.5:2.5;z4=(t>=0&t<1).*(t.^2)+(t>=1&t<2).*(t.^2-1)+(t>=2&t<3) .*(t.^2-2*t+1)4.完成下列操作：(1)求[100,999]之间能被21整除的数的个数。
(2)建立一个字符串向量，删除其中的大写字母。
解：(1)结果：m=100:999;n=find(mod(m,21)==0);length(n)ans=43(2).建立一个字符串向量例如：ch='ABC123d4e56Fg9';则要求结果是：ch='ABC123d4e56Fg9';k=find(ch>='A'&ch<='Z');ch(k)=[]ch=⎣O2⨯3⎥，其中E、R、O、S分别为单位矩阵、随机矩阵、零矩S⎦阵和对角阵，试通过数值计算验证A=⎢⎥。
《微分方程的数值解》课件
谱方法：将微分方程离散化为谱方程， 然后求解
边界元法：将微分方程离散化为边界 元方程，然后求解
有限元法：将微分方程离散化为有限 元方程，然后求解
网格法：将微分方程离散化为网格方 程，然后求解
数值解法的步骤
确定微分方程的初值 和边界条件
选择合适的数值解法， 如欧拉法、龙格-库塔 法等
实解
应用：广泛应 用于工程、物 理、化学等领
域
优缺点：优点 是计算速度快， 缺点是精度较
低
非线性方程的数值解法
牛顿法：通过迭 代求解非线性方 程
拟牛顿法：通过 迭代求解非线性 方程，比牛顿法 收敛更快
割线法：通过迭代 求解非线性方程， 适用于求解单变量 非线性方程
迭代法：通过迭 代求解非线性方 程，适用于求解 多维非线性方程
05 数值解法的实现
M AT L A B 编 程 实 现
MATLAB简介： MATLAB是一种高 级编程语言，广泛 应用于科学计算、 数据分析等领域
数值解法：包括欧 拉法、龙格-库塔 法、四阶龙格-库 塔法等
MATLAB实现：使 用MATLAB编写程 序，实现数值解法 的计算
示例代码：给出 MATLAB实现数值 解法的示例代码， 并解释其含义和作 用
设定时间步长和空间 步长
计算微分方程的解， 并进行误差分析
绘制解的图形，并进 行结果分析
对比不同数值解法的 优缺点，选择最优解 法
04 常用的数值解法
欧拉方法
基本思想：将微分 方程转化为差分方 程，然后求解差分 方程
优点：简单易行， 适用于初值问题
缺点：精度较低， 稳定性较差
改进方法：改进欧 拉方法，如改进欧 拉方法、龙格-库 塔方法等
微分方程数值解实验报告
微分方程数值解法课程设计报告班级：_______姓名： ___学号：__________成绩：2017年 6月 21 日摘要自然界与工程技术中的很多现象，可以归结为微分方程定解问题。
其中，常微分方程求解是微分方程的重要基础内容。
但是，对于许多的微分方程，往往很难得到甚至不存在精确的解析表达式，这时候，数值解提供了一个很好的解决思路。
，针对于此，本文对常微分方程数值解法进行了简单研究，主要讨论了一些常用的数值解法，如欧拉法、改进的欧拉法、Runge—Kutta方法、Adams法以及椭圆型方程、抛物型方程的有限差分方法等，通过具体的算例，结合MATLAB求解画图，初步给出了一般常微分方程数值解法的求解过程。
同时，通过对各种方法的误差分析，让大家对各种方法的特点和适用范围有一个直观的感受。
关键词：微分方程数值解、MATLAB目录摘要 (2)目录 (3)第一章常微分方程数值解法的基本思想与原理 (4)1.1 常微分方程数值解法的基本思路 (4)1.2用matlab编写源程序 (4)1.3 常微分方程数值解法应用举例及结果 (5)第二章常系数扩散方程的经典差分格式的基本思想与原理 (6)2.1 常系数扩散方程的经典差分格式的基本思路 (6)2.2 用matlab编写源程序 (7)2.3 常系数扩散方程的经典差分格式的应用举例及结果 (8)第三章椭圆型方程的五点差分格式的基本思想与原理 (10)3.1 椭圆型方程的五点差分格式的基本思路 (10)3.2 用matlab编写源程序 (10)3.3 椭圆型方程的五点差分格式的应用举例及结果 (12)第四章总结 (12)参考文献 (12)第一章常微分方程数值解法的基本思想与原理1.1常微分方程数值解法的基本思路常微分方程数值解法(numerical methods forordinary differential equations)计算数学的一个分支.是解常微分方程各类定解问题的数值方法.现有的解析方法只能用于求解一些特殊类型的定解问题，实用上许多很有价值的常微分方程的解不能用初等函数来表示，常常需要求其数值解.所谓数值解，是指在求解区间内一系列离散点处给出真解的近似值.这就促成了数值方法的产生与发展.1.2用matlab编写源程序龙格库塔法：M文件:function dx=Lorenz(t,x)%r=28,sigma=10,b=8/3dx=[-10*(x(1)-x(2));-x(1)*x(3)+28*x(1)-x(2);x(1)*x(2)-8*x(3)/3];运行程序：x0=[1,1,1];[t,y]=ode45('Lorenz',[0,100],x0);subplot(2,1,1) %两行一列的图第一个plot(t,y(:,3))xlabel('time');ylabel('z');%画z-t图像subplot(2,2,3) %两行两列的图第三个plot(y(:,1),y(:,2))xlabel('x');ylabel('y'); %画x-y图像subplot(2,2,4)plot3(y(:,1),y(:,2),y(:,3))xlabel('x');ylabel('y');zlabel('z');%画xyz图像欧拉法：h=0.010;a=16;b=4;c=49.52;x=5;y=10;z=10;Y=[];for i=1:800x1=x+h*a*(y-x);y1=y+h*(c*x-x*z-y);z1=z+h*(x*y-b*z);x=x1;y=y1;z=z1;Y(i,:)=[x y z];endplot3(Y(:,1),Y(:,2),Y(:,3));1.3常微分方程数值解法的应用举例及结果应用举例：⎪⎪⎪⎩⎪⎪⎪⎨⎧-=--=--=)()()()()()()()()())()(()(t bz t y t x dt t dz t z t x t y t rx dt t dy t y t x a dt t dx a=10,b=8/3,0<r<+∞,当1<r<24.74时，Lorenz 方程有两个稳定的不动点c()1(-r b ，)1(-r b ，r-1)和c '（-)1(-r b ，-)1(-r b ，r-1），一个稳定的不动点0=（0，0，0），当r>24.74时，c 和c '都变成不稳定的，此时存在混沌和奇怪吸引子。
第6章 数值微积分
第6章 数值微积分6.1 引言在科学研究和生产实际中，经常遇到求函数()f x 的导数（或微分）和求函数()f x 在闭区间[,]a b 上的定积分等问题. 由高等数学知道：若函数()f x 为初等函数，则可由基本求导公式和运算法则求出()f x 的导函数()f x '；若函数()f x 在区间[,]a b 上连续且原函数为()F x ，则由Newton-Leibniz 公式()d ()()ba f x x Fb F a =-⎰可求得()f x 在区间[,]a b 上的定积分.但是，在实际问题中遇到的求导数、求定积分的计算中，经常会有这样的情况：(1) 函数()f x 的原函数无法用初等函数给出，例如定积分121sin sin sin d ,d x xe x x x x⎰⎰等，从而无法用Newton-Leibniz 公式计算定积分。
(2) 只给出了函数()f x 的若干数据，函数表达式未知，因而无法用导数公式或定积分公式。
(3) 函数()f x 的导函数或原函数虽然能够求出，但形式过于复杂，不便使用.由此可见，利用求导公式及求导法则求导数或利用原函数求定积分有它的局限性。
因此，研究如何利用函数()f x 在有限个点处的函数值计算其导数的近似值，以及定积分的近似值的数值计算方法——数值微分(numerical differentiation)和数值积分(numerical integration)，既有理论意义又有实际应用价值。
这些方法也是微分方程和积分方程数值解法的基础。
6.2 数值微分6.2.1 插值型的数值微分公式已知数据(,)0,1,,i i x y i n = . 记 [,]a b 是包含插值节点01,,,n x x x 的区间,(),0,1,,i i y f x i n == ，则()f x 的n 次插值多项式为∑==ni i i n x l y x p 0)()(,其中()()()()n i i ni x l x x x x ωω='-, 0()()()n n x x x x x ω=-- ，误差为 (1)()()()()(), (,)(1)!n n n n f R x f x p x x a b n ξωξ+=-=∈+,故(1)(1)()()d ()()()()(1)!(1)!d n n n n n f x f x p x x f n n x ξωωξ++⎡⎤'''-=+⎢⎥++⎣⎦. (6.2.1) 因为ξ是x 的未知函数，所以无法对上式右端的第二项作出判断，因而对于任意给定的点x ，误差()nR x '是无法预估的. 如果只是求某个节点(0,1,,)i x i n = 上的导数值，这时有 (1)()()()()(1)!n i n i ni f f x p x x n ξω+'''-=+. (6.2.2) 以2n =为例，为简便起见，设01020,,,2x x x h x x h =+=+，已知001122(),(),()f x y f x y f x y ===.设满足插值条件200211222(),(),()p x y p x y p x y ===的2次插值多项式为2()p x ，即0201122012010210122021020112012222()()()()()()()()()()()()()()()()()()()()()(),22x x x x x x x x x x x x p x f x f x f x x x x x x x x x x x x x x x x x x x x x x x x x y y y h h h------=++------------=++-则0201122012222222()22x x x x x x x x x p x y y y h h h ------'=++-. (6.2.3) 于是得三点插值型数值求导公式1220012210202220121()(34)()231()()(),(,)(0,1,2)261()(43)()23i h f x y y y f h h f x y y f x x i h h f x y y y f h ξξξξ⎧⎪⎪⎪⎨⎪⎪⎪⎩=-+-+''''=-+-∈=''''=-++'''' (6.2.4)类似上述推导，在等距节点的情形，即0,0,1,,i x x ih i n =+= 时，可得如下常用的数值求导公式.(1) 两点公式(two-point formulas)0101101()()(),21()()().2hf x y y f h hf x y y f h ξξ⎧⎪⎪⎨⎪⎪⎩=--'''=-+''' (6.2.5)(2) 三点公式(three-point formulas)一阶求导公式：200122102220121()(34)(),231()()(),261()(43)().23h f x y y y f h hf x y y f h h f x y y y f h ξξξ⎧⎪⎪⎪⎨⎪⎪⎪⎩=-+-+''''=-+-''''=-++'''' (6.2.6) 二阶求导公式：2(4)20012122(4)10122(4)22012121()(2)()(),61()(2)(),121()(2)()().6h f x y y y hf f h h f x y y y f h h f x y y y hf f h ξξξξξ⎧⎪⎪⎪⎨⎪⎪⎪⎩=-+-+'''''=-+-''=-++-''''' (6.2.7) (3) 五点公式(five-point formulas)一阶求导公式：4(5)0012344(5)1012344(5)201344(5)301234401()(254836163)(),1251()(310186)(),12201()(88)(),12301()(618103)(),12201()(31612h f x y y y y y f h h f x y y y y y f h h f x y y y y f h h f x y y y y y f h f x y h ξξξξ=-+-+-+'=--+-+-'=-+-+'=-+-++-'=-'4(5)1234364825)().5h y y y y f ξ⎧⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎩+-++ (6.2.8) 二阶求导公式：(5)32001234(5)32101234(5)422012344(5)230123415()(351041145611)(),12611()(112064)(),121211()(163016)(),12901()(462011)(1220f x y y y y y h f h f x y y y y y h f h f x y y y y y h f h h f x y y y y y f h ξξξξ=-+-+-''=-++--''=-+-+-+''=-++-+-''(5)32401234),15()(115611410435)().126f x y y y y y h f h ξ⎧⎪⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎪⎩=-+-++''(6.2.9)注 上述数值求导公式不仅是计算导数近似值的有效方法，而且在微分方程数值解及样条插值的数据扩充方面都很有用处.明：五点公式比三点公式与两点公式准确，步长h越小结果越准确. 一般情况下，这个结论也是对的. 但由误差表示式可见，如果高阶导数无界，或者舍入误差超过截断误差时，这个结论就不对了。
matlab实验9答案详解
实验九数值微积分与方程数值求解一、实验目的1.掌握求数值导数和数值积分的方法2.掌握代数方程数值求解方法3.掌握常微分方程数值求解方法二、实验内容1.程序代码x=1;i=1;f=inline('det([x x^2 x^3;1 2*x 3*x^2;0 2 6*x])');while x<=3.01g(i)=f(x);i=i+1;x=x+0.01;endg;t=1:0.01:3.01;dx=diff(g)/0.01;f1=dx(1)f2=dx(101)f3=dx(length(g)-1)运行结果f1 =6.0602f2 =24.1202f3 =54.18022.程序代码f=inline('sqrt(cos(t.^2)+4*sin(2*t).^2+1)');I1=quad(f,0,2*pi)g=inline('log(1+x)./(1+x.^2)');I2=quad(g,0,2*pi)运行结果I1 =10.4285I2 =0.99973.程序代码A=[6 5 -2 5;9 -1 4 -1;3 4 2 -2;3 -9 0 2]; b=[-4 13 1 11]';x=A\by=inv(A)*b[L,U]=lu(A);z=U\(L\b)运行结果x =0.6667-1.00001.5000-0.0000y =0.6667-1.00001.50000.0000z =0.6667-1.00001.5000-0.00004.程序代码function [x,y]=line_solution(A,b)[m,n]=size(A);y=[ ];if norm(b)>0 %非齐次方程组if rank(A)==rank([A,b])if rank(A)==ndisp('有唯一解x');x=A\b;elsedisp('有无穷个解，特解x，基础解系y'); x=A\b;y=null(A,'r');endelsedisp('无解');x=[ ];endelse%齐次方程组disp('有零解x');x=zeros(n,1);if rank(A)<ndisp('有无穷个解，基础解系y');y=null(A,'r');endend命令format ratA=[2 7 3 1;3 5 2 2;9 4 1 7];b=[6 4 2]';[x,y]=line_solution(A,b)运行结果有无穷个解，特解x，基础解系yWarning: Rank deficient, rank = 2, tol = 8.6112e-015.> In line_solution at 11x =-2/1110/11y =1/11 -9/11-5/11 1/111 00 1 说明121/119/112/115/111/1110/11100010X k k --⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥=++⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦，其中12,k k 为任意常数。
微分方程数值解法-课件
只要f ( x, y)在a,bR1上连续, 且关于 y 满足
Lipschitz 条件，即存在与 x , y 无关的常数 L 使
|f ( x ,y 1 ) f ( x ,y 2 ) | L |y 1 y 2 |
对任意定义在 a , b 上的 y1x,y2x都成பைடு நூலகம்，
（2）用代数的方法求出解函数 y解y(在的x)唯点一的x性k近似值 y k y (x k ) 工k 程 师1 ,关2 , 注,n 解yk的光滑性 解*的振动y( 性xk )
解的周期性
数学界关注
yy解(x的) 稳定性 解的混沌性
……
所谓数值解法：
求函数 y(x) 在一系列节点 a = x0< x1<…< xn= b 处的近似值
称 f ( x, y) 在区域D上对 y 满足Lipschitz条件是指:
L0 s.t. f(x,y1)f(x,y2) Ly1y2 ,
x[a,b],y1,y2[y(x),y(x)]
4、 迭代格式的构造
(1) 构造思想：将连续的微分方程及初值条件离散为线性方程 组加以求解。由于离散化的出发点不同，产生出各种不同的数 值方法。基本方法有：有限差分法（数值微分）、有限体积法 （数值积分）、有限元法（函数插值）等等。
离变量法”等特殊方法求得初等函数形式的解，绝大
部分方程至今无法理论求解。
如
y ' sx i2 n )y , y '( 1 x sy i,ny ' e x 2 xy 等等
2、数值解的思想
如果找不到解函数
（1）将连续变量 x[a离,b散] 为 数学界还关注：
a x 0 x 1 x k x n 解 的b 存在性
微分方程数值解法实验报告
微分方程数值解法实验报告实验题目：数值解微分方程的实验研究引言：微分方程是描述自然现象、科学问题和工程问题中变量之间的关系的重要数学工具。
然而，大部分微分方程很难找到解析解，因此需要使用数值方法来近似求解。
本实验旨在研究数值解微分方程的方法和工具，并通过实验验证其有效性和准确性。
实验步骤：1.了解微分方程的基本概念和求解方法，包括欧拉法、改进的欧拉法和龙格-库塔法。
2. 配置实验环境，准备实验所需的工具和软件，如Python编程语言和相关数值计算库。
3.选择一种微分方程进行研究和求解，可以是一阶、二阶或更高阶的微分方程。
4.使用欧拉法、改进的欧拉法和龙格-库塔法分别求解选定的微分方程，并比较其结果的准确性和稳定性。
5.计算数值解与解析解之间的误差，并进行误差分析和讨论。
6.对比不同数值解法的性能，包括计算时间和计算精度。
7.结果展示和总结，根据实验结果对数值解方法进行评价和选取。
实验结果：以一阶线性常微分方程为例，我们选择经典的“衰减振荡”问题进行实验研究。
该问题的微分方程形式为：dy/dt = -λy其中，λ为正实数。
我们首先使用Python编程语言实现了欧拉法、改进的欧拉法和龙格-库塔法。
进一步，我们选择了λ=0.5和初始条件y(0)=1，使用这三种数值解法求解该微分方程，并比较结果的准确性。
通过对比数值解和解析解可以发现，在短时间内，欧拉法、改进的欧拉法和龙格-库塔法的结果与解析解非常接近。
但随着时间的增加，欧拉法的结果开始偏离解析解，而改进的欧拉法和龙格-库塔法仍然能够提供准确的近似解。
这是因为欧拉法采用线性逼近的方式，误差随着步长的增加而累积，而改进的欧拉法和龙格-库塔法采用更高阶的逼近方式，可以减小误差。
为了更直观地比较不同方法的性能，我们还计算了它们的计算时间。
实验结果显示，欧拉法计算时间最短，而龙格-库塔法计算时间最长。
这表明在计算时间要求较高的情况下，可以选择欧拉法作为数值解方法。
数学建模-- 常微分方程数值解及实验
注意:
1、在解n个未知函数的方程组时，x0和x均为n维向量， m-文件中的待解方程组应以x的分量形式写成. 2、使用Matlab软件求数值解时，高阶微分方程必须 等价地变换成一阶微分方程组.
例
d 2x dx 2 x 0 2 ( x 1) dt dt x (0 ) 2; x '(0 ) 0
实际应用时，与欧拉公式结合使用：
0 y i( 1) y i hf ( x i , y i ) h ( k 1) (k ) y i 1 y i [ f ( x i , y i ) f ( x i 1 , y i 1 )] k 0 ,1, 2 , 2
•欧拉法是一阶公式，改进的欧拉法是二阶公式。
•龙格-库塔法有二阶公式和四阶公式。 •线性多步法有四阶阿达姆斯外插公方程的数值解
[t，x]=solver（’f’,ts,x0,options）
自变 量值 函数 值
ode45 ode23 ode113 ode15s ode23s
1、建立m-文件vdp1.m如下： function dy=vdp (t,y) dy=zeros(2,1); dy(1)=y(2); dy(2)=(1-y(1)^2)*y(2)-y(1); 2、取t0=0，tf=20，输入命令： [T,Y]=ode15s('vdp1',[0 20],[2 0]); plot(T,Y(:,1),'-')
y ( x i 1 ) y ( x i )
xi 1
f ( t , y ( t )) dt
x i 1 x i 2
xi
[ f ( x i , y ( x i )) f ( x i 1 , y ( x i 1 ))]
常微分方程的数值解法实验报告
常微分方程的数值解法实验报告实验报告：常微分方程的数值解法摘要：常微分方程（ODE）是描述动力学系统中物理量随时间变化的数学方程，广泛应用于自然科学和工程领域。
然而，对于一些复杂的非线性ODE，很难找到解析解。
因此，我们需要数值解法来求解这些方程。
本实验报告将介绍四种常见的常微分方程数值解法，分别是欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法，并通过数值实验比较它们的精度和效率。
1.引言在实际问题中，许多物理量的变化规律可以由常微分方程描述。
然而，对于复杂的非线性ODE，很难找到解析解。
因此，为了解决这类问题，我们需要借助数值方法来求解。
2.方法本实验采用四种常见的常微分方程数值解法：欧拉法、改进的欧拉法、四阶龙格-库塔法和自适应步长的龙格-库塔法。
（1）欧拉法是最简单的数值解法，通过将微分方程转化为差分方程，使用离散的步长来近似微分方程。
（2）改进的欧拉法在欧拉法的基础上进行了改进，使用预估-校正的方法来提高精度。
（3）四阶龙格-库塔法是一种经典的数值解法，通过利用不同步长处的斜率来近似微分方程，具有较高的精度。
（4）自适应步长的龙格-库塔法是在四阶龙格-库塔法的基础上改进而来的，根据步长的大小自适应地选择不同的步长，同时保证精度和效率。
3.实验设计为了比较这四种数值解法的精度和效率，我们设计了两个实验。
实验一是求解一阶常微分方程：dy/dx = -2x，初始条件y(0) = 1，解析解为y = 1 - x^2、实验二是求解二阶常微分方程：d^2y/dx^2 + y = 0，初始条件y(0) = 0，dy/dx(0) = 1，解析解为y = sin(x)。
4.结果与分析实验一中，比较四种数值解法在不同步长下的近似解和解析解，计算其误差。
实验结果表明，四阶龙格-库塔法和自适应步长的龙格-库塔法具有较高的精度，而欧拉法和改进的欧拉法的精度较低。
实验二中，我们比较四种数值解法在不同步长下的近似解和解析解，并计算其误差。
微分方程数值解试验报告
中国矿业大学（北京）理学院微分方程数值解实验报告实验名称：常微分方程数值解法学号： XXX专业年级：信息与计算科学11级姓名： XXX实验时间： 2013年3月成绩：一、实验目的，内容 二、相关背景知识介绍 三、代码四、数值结果与分析 六、计算结果的分析一、实验目的，内容掌握一阶常微分方程的初值问题常用数值解法编写程序，利用Euler 法、改进Euler 法、Adams 外插内插法及Runge-Kutta 法解决相关问题。
二、相关背景知识介绍研究一阶常微分方程的初值问题的数值解1、Euler 法和改进Euler 法 （1）Euler 法用差商近似导数问题转化为（2）改进Euler 法]2/),(),([111++++=+n n n n n u t f u t f h y y n00(,)()dy f x y a x bdx y x y ⎧=⎪≤≤⎨⎪=⎩dxdy h y y n n ≈-+1,...)3,2,1,0()(),(001=⎩⎨⎧=+=+n x y y y x hf y y n n n n（3）线性多步法（Adams 外插内插法）j n k j j j n k j j f h u+=+=∑∑=00βα （4）Runge-Kutta 法用f 在一些点的值表示ψ（tn,un,h ）,使单步局部截断误差的阶和Taylor 展开法相等。
将初值问题写成积分形式：di i u i f t u h t u h t t ⎰++=+))(,()()(三、代码（Matlab ）题目：用Euler 法、改进Euler 法、Adams 外插内插法及Runge-Kutta 法解决下面一阶常微分方程初值问题。
代码:方程：function y=f(t,u)y=2*t*u/(1+t^2);endEuler 法:function y=euler(n,u0,x0)t=1;u=zeros(n+1,1);h=t/n;u(1)=u0;x(1)=x0;for i=1:nx(i+1)=x0+h*i;u(i+1)=u(i)+f(x(i),u(i))*h;endPlot(x,u,'r*')end改进Euler法:function y=imeuler(n,u0,x0)err=10^-6;t=1;u=zeros(n+1,1);x=zeros(n+1,1);h=t/n;u(1)=u0;x(1)=x0;for i=1:nx(i+1)=x0+h*i;u(i+1)=u(i);ut=30;while abs(u(i+1)-ut)>errut=u(i+1);u(i+1)=u(i)+(f(x(i),u(i))+f(x(i+1),u(i+1)))*0.5*h;endendplot(x,u)endAdams内插法：k=3function y=iadams(n,u0,x0)t=1;u=zeros(n+1,1);h=t/n;u(1)=u0;u(2)=u0;u(3)=u0;x(1)=x0;x(2)=x0+h;x(3)=x0+2*h;for i=3:nx(i+1)=x0+h*i;u(i+1)=u(i)+h*(9*f(x(i+1),u(i+1))+19*f(x(i),u(i))-5*f(x(i-1),u(i-1))+f(x(i-2),u(i-2)))/24;endplot(x,u)endRunge-Kutta法：（本文采用heun法）function y=heun(n,u0,x0)t=1;u=zeros(n+1,1);h=t/n;u(1)=u0;x(1)=x0;for i=1:nx(i+1)=x0+h*i;k1=f(x(i),u(i));k2=f(x(i)+0.3333*h,u(i)+0.33333*h*k1);k3=f(x(i)+0.666663*h,u(i)+0.666666*h*k2);u(i+1)=u(i)+0.25*h*(k1+3*k3);endplot(x,u)end四、数值结果与分析对于所给方程，各个方法均给出了良好的解，其实heun方法最好10步就出了的结果。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
实验09 数值微积分与方程数值求解（第6章 MATLAB 数值计算）一、实验目的1. 掌握求数值导数和数值积分的方法。
2. 掌握代数方程数值求解的方法。
3. 掌握常微分方程数值求解的方法。
二、实验内容1. 求函数在指定点的数值导数232()123,1,2,3026x x x f x x x x x==程序及运行结果：2. 用数值方法求定积分(1) 22210cos 4sin(2)1I t t dt π=++⎰的近似值。
程序及运行结果：《数学软件》课内实验王平(2) 2221I dx xπ=+⎰程序及运行结果：3. 分别用3种不同的数值方法解线性方程组6525494133422139211x y z u x y z u x y z u x y u +-+=-⎧⎪-+-=⎪⎨++-=⎪⎪-+=⎩ 程序及运行结果：4. 求非齐次线性方程组的通解1234123412342736352249472x x x x x x x x x x x x +++=⎧⎪+++=⎨⎪+++=⎩5. 求代数方程的数值解(1) 3x +sin x -e x =0在x 0=1.5附近的根。
程序及运行结果（提示：要用教材中的函数程序line_solution ）：(2) 在给定的初值x 0=1，y 0=1，z 0=1下，求方程组的数值解。
23sin ln 70321050y x y z x z x y z ⎧++-=⎪+-+=⎨⎪++-=⎩6. 求函数在指定区间的极值(1) 3cos log ()xx x x xf x e ++=在(0,1)内的最小值。
(2) 33212112122(,)2410f x x x x x x x x =+-+在[0,0]附近的最小值点和最小值。
7. 求微分方程的数值解，并绘制解的曲线2250(0)0'(0)0xd y dyy dx dx y y ⎧-+=⎪⎪⎪=⎨⎪=⎪⎪⎩程序及运行结果（注意：参数中不能取0，用足够小的正数代替）：令y 2=y,y 1=y '，将二阶方程转化为一阶方程组：'112'211251(0)0,(0)0y y y x x y y y y ⎧=-⎪⎪=⎨⎪==⎪⎩8. 求微分方程组的数值解，并绘制解的曲线123213312123'''0.51(0)0,(0)1,(0)1y y y y y y y y y y y y =⎧⎪=-⎪⎨=-⎪⎪===⎩程序及运行结果：三、实验提示四、教程：第6章 MATLAB 数值计算（2/2）6.2 数值微积分 p155 6.2.1 数值微分1. 数值差分与差商对任意函数f(x)，假设h>0。
➢ 向前差分：()()()f x f x h f x ∆=+- ➢ 向后差分：()()()f x f x f x h ∇=--➢ 中心差分：()(/2)(/2)f x f x h f x h δ=+-- 当步长h 充分小时，有➢ 向前差商：()'()f x f x h ∆≈➢ 向后差商：()'()f x f x h ∇≈➢ 中心差商：()'()f x f x hδ≈2. 数值微分的实现MATLAB 没有提供求数值导数的函数，只有计算向前差分的函数diff 。
➢ DX=diff(X)：计算向量X 的向前差分，DX(i)=X(i+1)-X(i)，i=1,2,…,n -1。
➢ DX=diff(X,n)：计算X 的n 阶向前差分。
例如，diff(X,2)=diff(diff(X))。
➢ DX=diff(A,n,dim)：计算矩阵A 的n 阶差分，dim=1时(缺省状态)，按列；dim=2，按行。
例6.18 （向前差分）求1~3阶差分p156设x 由[0,2π]间均匀分布的10个点组成，求sin x 的1~3阶差分。
例6.19 （数值导数）3种方法求导p157=+()52f x x用不同的方法求f(x)的数值导数，并在同一个坐标系中做出f '(x)的图像。
6.2.2 数值积分 p1571. 数值积分基本原理 求解定积分的数值方法： ➢ 梯形法➢ 辛普生(Simpson)•法➢ 牛顿－柯特斯(Newton-Cotes)法 基本思想：将整个积分区间[a ,b ]分成n 个子区间 [x i, x i +1]，i =1,2,…,n ，其中x 1=a ，x n +1=b这样求定积分问题就分解为求和问题。
2. 数值积分的实现(1) 被积函数是一个解析式调用格式： quad(fname,a,b,tol,trace)quadl(fname,a,b,tol,trace)➢ fname 是被积函数名。
➢ a 和b 分别是定积分的下限和上限。
➢ tol 用来控制积分精度，默认tol=10-6。
➢ trace 控制是否展现积分过程，若取非0则展现积分过程，取0则不展现，默认trace=0。
例6.20（解析式函数数值积分）两种方法p158用两种不同方法求 12x I edx -=⎰调用格式：trapz(X,Y) ➢ Y=f(X)➢ X=[x 1,x 2,...,x n ], x i <x i+1(i=1,2,...,n-1)例6.21（表格函数数值积分）p159用trapz 函数计算 120x I e dx -=⎰(3) 二重积分数值求解(,)dbcaI f x y dxdy =⎰⎰调用格式：I=dblquad(f,xmin,xmax,ymin,ymax,tol,trace)例6.22（二重数值积分）p160212/2212sin()xI e x y dxdy ---=+⎰⎰表 数值微积分函数和命令 p155~1606.3 离散傅立叶变换 p160 6.3.1 离散傅立叶变换算法简述在某时间片等距抽取N 个抽样时间t m 处的样本值f(t m )，记f(m)，m=1,2,...,N 。
➢ 称向量F(k)(k=1,2,...,N)为f(m)的一个离散傅里叶变换，公式：2(1)(1)/1()(),1,2,,Nj m k N m F k f m e k N π---===∑又称F(k)为f(m)的离散频谱。
➢ 由F(k)逆求f(m)的过程，称离散傅立叶逆变换，公式：2(1)(1)/1()(),1,2,,Nj m k N k f m F k e m N π--===∑6.3.2 离散傅立叶变换的实现 p161一维离散傅立叶变换函数：(1) fft(X)：返回向量X 的离散傅立叶变换。
➢ 设X 的元素个数为N 。
➢ 若N 为2的幂次，则为以2为基数的快速傅立叶变换。
➢ 否则为运算速度很慢的非2幂次的算法。
➢ 对于矩阵X ，它应用于矩阵的每一列。
(2) fft(X,N)：计算N 点离散傅立叶变换。
➢ 它限定向量的长度为N 。
➢若X的长度小于N，则不足部分补零。
➢若大于N，则删去超出N的那些元素。
➢对于矩阵X，它应用于矩阵的每一列。
(3)fft(X,[ ],dim)或fft(X,N,dim)：前者与fft(X)基本相同，后者与fft(X,N)基本相同。
➢dim=1时，作用于X的每一列。
➢dim=2时，作用于X的每一行。
注：当已知给出的样本数N0不是2的幂次时，可取一个N使它大于N0且是2的幂次，然后利用函数格式fft(X,N)或fft(X,N,dim)便可进行快速傅立叶变换。
一维离散傅立叶逆变换函数，调用格式：(1)ifft(F)：返回F的一维离散傅立叶逆变换。
(2)ifft(F,N)：为N点逆变换。
(3)ifft(F,[ ],dim)或ifft(F,N,dim)：则由N或dim确定逆变换的点数或操作方向。
例6.23 （快速傅里叶变换）p161给定数学函数x(t)=12sin(2π×10t+π/4)+5cos(2π×40t)取N=128，试对t从0~1秒采样，用fft作快速傅立叶变换，绘制相应的振幅-频率图。
在0~1秒时间范围内采样128点，确定采样周期和采样频率。
由于离散傅立叶变换时的下标应是从0到N-1，故在实际应用时下标应该前移1。
ix=real(ifft(X));%求逆变换，结果只取实部plot(t,x,t,ix,':')%逆变换结果和原函数的曲线norm(x-ix)%逆变换结果和原函数之间的距离ans =2.2403e-014表 离散傅立叶变换函数和命令 p160~163函数/命令说 明fft 一维离散傅立叶变换函数 ifft一维离散傅立叶逆变换函数6.4 线性方程组求解 p163 6.4.1 直接解法1. 利用左除运算符的直接解法对于线性方程组Ax =b ，可利用左除运算符“\”求解：x=A\b➢ A 为方阵，若是奇异的，则给出警告信息。
➢ b 为列向量，得一个解。
➢ b 为矩阵，得多个解。
例6.24（左除）求解线性方程组p16312341242341234251357926640x x x x x x x x x x x x x x +-+=⎧⎪-+=-⎪⎨+-=⎪⎪+--=⎩2. 利用矩阵的分解求解线性方程组矩阵分解指按一定的原理用某种算法将一个矩阵分解成若干个矩阵的乘积。
常见的矩阵分解： ➢ LU 分解 ➢ QR 分解➢ Cholesky 分解 ➢ Schur 分解➢ Hessenberg 分解 ➢ 奇异分解 (1) LU 分解将一个矩阵表示为一个下三角阵和一个上三角阵的乘积形式。
已经证明，只要方阵A 是非奇异的，LU 分解是可行的。
调用格式：➢ [L,U]=lu(A)：产生一个上三角阵U 和一个变换形式的下三角阵L (行交换)，满足A=LU 。
➢ [L,U,P]=lu(A)：产生一个上三角阵U 和一个下三角阵L 及一个置换矩阵P ，满足PA=LU 。
线性方程组Ax=b 的解：x=U\(L\b) 或 x=U\(L\P*b) 注意：使用第一种格式时，矩阵L 往往不是一个下三角矩阵，但可通过行交换成一个下三角阵。
例 LU 分解 p164111543211A -⎡⎤⎢⎥=-⎢⎥⎢⎥⎣⎦例6.25（LU 分解）求解线性方程组p16612341242341234251357926640x x x x x x x x x x x x x x +-+=⎧⎪-+=-⎪⎨+-=⎪⎪+--=⎩A 为方阵，调用格式：➢ [Q,R]=qr(A)：产生一个正交矩阵Q 和一个上三角矩阵R ，满足A=QR 。
➢ [Q,R,E]=qr(A)：产生一个正交矩阵Q 、一个上三角矩阵R 及一个置换矩阵E ，满足AE=QR 。
线性方程组Ax=b 的解：x=R\(Q\b) 或 x=E*(R\(Q\b))例 QR 分解 p1661115432710A -⎡⎤⎢⎥=-⎢⎥⎢⎥⎣⎦ans =1.0000 -1.0000 1.00005.0000 -4.0000 3.00002.0000 7.0000 10.0000>>例6.26（QR分解）求解线性方程组p16712341242341234251357926640x x x xx x xx x xx x x x+-+=⎧⎪-+=-⎪⎨+-=⎪⎪+--=⎩(3) Cholesky分解A是对称正定矩阵，调用格式：➢R=chol(A)：产生一个上三角阵R，使R'R=A。