MATLAB与数值分析课程总结

MATLAB与数值分析课程总结

姓名：董建伟 学号：2015020904027

一：MATLAB部分 1. 处理矩阵－容易

矩阵的创建

（1） 直接创建 注意 a中括号里可以用空格或者逗号将矩阵元素分开 b矩阵元素可以是任何MATLAB表达式，如实数复数等

c可以调用赋值过的任何变量，变量名不要重复，否则会被覆盖

（2） 用MATLAB函数创建矩阵 如：a空阵［］ b rand/randn——随机矩阵 c eye——单位矩阵 d zeros——0矩阵 e ones——1矩阵 f magic——产生n阶幻方矩阵 等

向量的生成

（1） 用冒号生成向量

（2） 使用linspace和logspace分别生成线性等分向量和对

数等分向量

矩阵的标识和引用 （1） 向量标识

（2） “0 1”逻辑向量或矩阵标识 （3） 全下标，单下标，逻辑矩阵方式引用 字符串数组

（1）字符串按行向量进行储存

（2）所有字符串用单引号括起来

（3）直接进行创建 矩阵运算

（1） 注意与数组点乘，除与直接乘除的区别，数组为乘方对应元素

的幂

（2） 左右除时斜杠底部靠近谁谁是分母

（3） 其他运算如，inv矩阵求逆，det行列式的值, eig特征值，diag

对角矩阵

2.绘图－轻松

plot－绘制二维曲线

（1） plot(x)绘制以x为纵坐标的二维曲线

plot(x,y) 绘制以x为横坐标，y为纵坐标的二维曲线 x,y为向量或矩阵

（2） plot（x1,y1,x2,y2，。。。 。。。）绘制多条曲线，不同字母代替不

同颜色:b蓝色，y黄色，r红色，g绿色

（3） hold on后面的plot图像叠加在一起

hold off解除hold on命令，plot将先冲去窗口已有图形

（4） 在hold后面加上figure，可以绘制多幅图形 （5） subplot在同一窗口画多个子图 三维图形的绘制

(1) plot3(x,y,z,’s’) s是指定线型，色彩，数据点形的字

符串

(2) ［X,Y］＝meshgrid(x,y)生成平面网格点 (3) mesh(x,y,z，c)生成三维网格点，c为颜色矩阵 (4) 三维表面处理mesh命令对网格着色，surf对网格片着色 (5) contour绘制二维等高线

(6) axis(［x1,xu,y1,yu］)定义x,y的显示范围

3．编程－简洁

（1）变量命名时可以由字母，数字，下划线，但是不得包含空格和标点

（2）最常用的数据类型只有双精度型和字符型，其他数据类型只在特殊条件下使用

（3）为得到高效代码，尽量提高代码的向量化程度，避免使用循环结构

（4）为得到最快的运行速度，在循环指令前尽量对数组进行预定义 （5）流程控制语句

如a break命令可以使包含break的最内层的for或while语句强制终止，跳出循环结构，执行end后面的语句。Break一般与if结构结合使用。

b continue命令用于结束本次for或while循环，只结束本次循环而进行下次循环。

c pause 命令用来使程序运行暂停，等待用户按任意键继续。但是pause(n)表示暂停n秒 等

4．常见字符含义

（1）disp和fprintf表示函数输出 （2）eps机器的精度 （3）inf零做了除数的输出 （4）NAN零被零除的输出 （5）int积分 （6）diff微分 （7）taglor泰勒展开

（8）syms声明符号变量 等

二：数值分析部分

基本概念

1. 误差与有效数字

误差、误差限、相对误差、相对误差限和有效数字的定义及相互关系；

误差的来源和误差的基本特性； 误差的计算（估计）的基本方法。 2. 算法的适定性问题

数值分析中的病态和不稳定性问题；

病态问题和不稳定算法的实例分析。 3. 数值计算的几个注意问题

数值计算的基本概念

? 误差概念和分析

误差的定义：

设x是精确值，p是近似值，则定义两者之差是绝对误差: ?a?x?p

由于精确值一般是未知的,因而Δ不能求出来,但可以根据测量误差或计算情况估计它的上限

|x-p|??

?称为绝对误差限。相对误差定义为绝对误差与精确值之比

?r??a x?a?r???x

? 误差的来源：

称为相对误差限

舍入误差

将无限位字长的精确数处理成有限位字长近似数的处理方法称为舍入方法。带来舍人误差。

截断误差

用数值法求解数学模型时，往往用简单代替复杂,或者用有限过程代替无限过程所引起的误差。

? 有效数字

对于

a=a0 a1 ? am . am+1 ? am+n(a0≠0)

的近似数， 若|Δ|≤0.5x10-n，

则称a为具有m+n+1位有效数字的有效数，其中每一位数字都叫做a的有效数字。有效数和可靠数的最末位数字称为可疑数字

有效数位的多少直接影响到近似值的绝对误差与相对误差的大小。 推论1 对于给出的有效数，其绝对误差限不大于其最末数字的半个单位。

x??0.a1a2?an?10m?1??x?x??10m?n2

推论2 对于给出的一个有效数，其相对误差限可估计如下：

x??0.a1a2?an?10m5?r(x)??10?na1? 数值计算的算法问题

“良态”问题和“病态”问题

在适定的情况下，若对于原始数据很小的变化δX，对应的参数误差δy也很小，则称该数学问题是良态问题；若δy很大，则称为病态问题。

病态问题中解对于数据的变化率都很大，因此数据微小变化必将导致参数模型精确解的很大变化。

数学问题的性态完全取决于该数学问题本身的属性，在采用数值方法求解之前就存在，与数值方法无关。

稳定算法和不稳定算法

如果用数值方法计算时,误差在计算过程中不扩散的算法称为稳定算法。否则称为不稳定算法。

? 数值计算应注意的问题

避免相近二数相减； 避免小分母； 避免大数吃小数； 选用稳定的算法。

绝对误差的运算：

?(x1?x2)??(x1)??(x2)

?(x1x2)?x1?(x2)?x2?(x1)

x??(x?)?x??(x?)x??()? ?x?x??(f(x))?f?(x)?(x)

线性方程组求解的数值方法：

高斯消元法

（1）性质1：若A的所有顺序主子式均不为0，则高斯消元无需换行即可进行到底，得到唯一解。

（2）性质2：只要 A 非奇异，即 A 存在，则可通过逐次消元及行交换，将方程组化为三角形方程组，求出唯一解

（3）列主元消去法：如果对角线元素为0，则需要交换元素。列主元消去法值选择该列 k~n 元素中最大值。

（4）全主元消去法在第k步消去前， 在系数矩阵右下角的n-k+1阶主子阵中，选绝对值最大的元素作为主元素

?1

矩阵分解

矩阵LU分解的一般计算公式； 利用LU分解的线性方程组求解方法； Cholesky分解；

Matlab的Cholesky分解函数。

向量范数与矩阵范数

向量范数及其性质； 矩阵函数及其性质； 常用范数形式。

线性方程组的迭代法求解

迭代求解的思路； Jacobi迭代法； 高斯_赛德尔迭代法； 迭代法的收敛性。

方程组的病态问题与误差分析

线性方程组解的误差分析； 条件数和方程组的病态性。

消元法：

问题：

消去法是按照系数矩阵的主对角线上的元素（主元）进行消元。从而可能出现：

（1）某个主元为零，导致消元过程无法进行。 （2）当某个主元的绝对值很小时，计算结果误差很大。 全主元消去法

每一步选绝对值最大的元素为主元素。

|aikjk|?max|aij|?0;k?i,j?n列主元消去法

省去换列的步骤，每次仅选一列中最大的元。

|aik,k|?max|aik|?0k?i?n矩阵三角分解法

三角矩阵的性质：

1. 上（下）三角方阵的行列式的值等于对角线元素的乘积； 2. 上（下）三角方阵的转置为下（上）三角矩阵；

3. 上（下）三角方阵的逆矩阵为上（下）三角矩阵，且对角元是原三角矩阵对角元的倒数；

4. 两个上（下）三角方阵的乘积也是上（下）三角矩阵，且对角元是原三角矩阵对角元的乘积。

由 Gauss 消去法加上列主元 ( 或全主元）有 LU 分解： ?LUA?1?l211 ?L?? l31l321??????ln1ln2??u11???? ? U? ??????ln(n?1)1???uu1222

??2n???????unn????1nuu 计算公式：

u?a , j?1, ?,nl?au , i?2, ?, n1j1ji1i111对 k?2, 3, ?, n 计算 ukj?akj??lkmumj j?k , ? , nm?1k-1k-1 lik?(aik??limumk )/ukk i?k?1 , ? , nm?1

éêêêêêê?aaaa11213141aaaa1112223242aaaa13233343aaaa121214243444ùéúê1úêl21ú=êúêl31úêúêl??4101001ll3242l43é0ùêu11ú0úê0úê0úê0ê1úê0ú??uu122200uuuuuu0u13233314211414443444ùúúúúúú?ùúúúúúú?éêê=êêêê ?uululu+ululu+lululu+lu21112131113141121232411142222222ulu+ulu+lu+ulu+lu+lu132113233113322341134223433333ulu+ulu+lu+ulu+lu+lu+u24311432243441144224433444·?·?·?·?·?·?·

按颜色顺序依次计算

不是所有矩阵都可分解为A=LU，更一般形式PA=LU

MATLAB中lu为LU分解函数，调用方式：[L,U,P] = lu(A)

Cholesky分解：

定理：

设矩阵A对称正定，则存在非奇异下三角阵L使得A?LLT 。若限定L对角元为正，则分解唯一。

Matlab中的Cholesky分解函数：

chol（）

? 正定阵的判定：

– A 亦对称正定，且 aii > 0 – A 的顺序主子阵 Ak 亦对称正定 – A 的特征值 ?i > 0 – det ( Ak ) > 0

?1

向量和矩阵的范数

为了研究线性方程组近似解的误差估计和迭代法的收敛性，引进向量（矩阵）的范数的概念。 向量范数 定义：

x,y?R满足下列条Rn空间的向量范数 || · || 对任意

件：

??n????（正定性) (1)||x||?0;||x||?0?x?0??(2)||?x||?|?|?||x||(齐次性)

??C????(3)||x?y||?||x||?||y||常用范数：

n（三角不等式)

x1??i?1xi1/2?2?x2???xi??i?1?n

定义：

xp?p????xi??i?1?n1/px矩阵范数

??maxxii

R件：

n?mn?mA,B?R空间的向量范数 || · || 对任意 满足下列条

(1)||A||?0;||A||?0?A?0

(2)||?A||?|?|?||A||

(3)||A?B||?||A||?||B||

(4) || AB || ? || A || · || B ||

常用矩阵范数： Frobenius范数：

由向量范数 || · ||p导出关于矩阵 A ? Rn?n的p范数:

||A||F???|ai?1j?1nnij|2?||Ax||p?||A||p?max?max||Ax||p????x?0||x||p?1||x||p特别有：

||A||??max?|aij|1?i?nnn（行和范数） （列和范数） （谱范数 ）

j?1||A||1?max?|aij|1?j?ni?1T||A||2??max(AA) 谱半径： 矩阵A的谱半径记为

? (A) =maxi?i，其中?i为A的特征根。 定理：

对任意算子范数 || · || 有 ?(A)?||A||定理：

若A对称，则有

A2??(A)

定理：

若矩阵B对某个算子范数满足 ||B|| < 1，则必有(1)(I?B)可逆。

（2）?I?B??1?11?||B||

解线性方程组的迭代法 研究内容：

? 如何建立迭代公式？ ? 收敛速度？

? 向量序列的收敛条件？ ? 误差估计？

（1）雅可比(Jacobi)迭代法

（

2

）

高

斯

—

塞

德

尔

迭

代

法

迭代法的收敛性

迭代矩阵M谱半径小于1.

（3）超松弛迭代 / SOR）

SOR是Gauss-Seidel 迭代法的一种加速法。 (k?1)假设:x(k)已知，xG为Gauss-Seidel 迭代法结果，(k?1)定义:?x?xG-x(k)，则,SOR迭代法的表达式为： x(k?1)?x(k)???x(k?1)(k?1)(k?1)x(k?1)?x(k)???x?x(k)??xG??x(k)? xG?(??1)(xG?x(k))其中，?称为松弛因子。当??1时称为低松弛；??1是Gauss-Seidel迭代；??1时称为超松弛法。

线性方程组迭代收敛性

1.2. 3. 4.

推论2： 松弛法收敛的必要条件是0???2。证明：设松弛法的迭代矩阵M有特征值?1,?2,...,?n。因为 det(M)??1?2...?n?[?(M)]-1n由定理，松弛法收敛必有 det(M)?1又因为 det(M)?(D-?L)(1-?)D??U1 (D-?L)?a11a22...ann-1 (1-?)D??U?(1-?)na11a22...ann?det(M)?(1-?)n?1?0???2。

5.

设有线性方程组Ax?b,下列结论成立：1.若A为严格对角占优阵或不可约弱对角占优阵，则Jacobi迭代法和Gauss-Seidel迭代法均收敛。2.若A为严格对角占优阵，0???1,则松弛法收敛。3.若A为对称正定阵，则松弛法收敛的充要条件为0???2。

迭代法的误差估计：

if x(k?1)?Mx(k)?g,M?1,?x??x,then(k)* x If(k)(k)?x?*Mk1?Mx(1)?x(0)

x(k?1)?Mx(k)?g,M?1,*?x??x,then x(k)?x?*M1?Mx(k)?x(k?1)

解线性方程组的误差分析

问题：

设A为非奇异矩阵，非齐次方程组Ax?b的准确解为x,当A和b有一个小的扰动?A,?b时，方程组有准确解x??x,即(A??A)(x??x)?b??b我们需要研究?x与?A,?b之间的关系。

? b有扰动，A无扰动

||?x||||?b||?(||A||||A?1||) ||x||||b||

? A有扰动，b有扰动

||?x||?||(I?A?1?A)?1||||A?1||||?A||||x|| ||?x||?(I?A?1?A)?1A?1||x||||A?1||||A||||A?1||||?A||?A?1?||A?1?A||?A?A?||A||?||A?1||||A||?A||A||1?||A?1||||A||||A||

? A有扰动，b有扰动

类似推导，可得||?x||??||x||?A ||A||?b(?)?A||b||||A||1?||A?1||||A||||A||||A?1||||A||?A 定义：

条件数: cond(A) = ||A–-1 || ||A||

条件数与所取的矩阵范数有关。常用的条件数有：（）1cond(A)??A?A?1（2）cond(A)2?A2?

A?12?max(ATA)??min(ATA)当条件数很大时,方程组 Ax = b是病态问题;当条件数较小时,方程组 Ax = b是良态问题。

函数的数值逼近：

1. 代数多项式插值问题

插值基函数和插值多项式； 插值的误差分析；

多项式插值的Runge现象。

2. 分段低次插值

分段线性插值；

Hermite插值和分段Hermite插值。

3. 三次样条插值

样条插值的定义； 三次样条函数的计算；

4. 曲线拟合的最小二乘法

曲线拟合的最小二乘法法； 多项式拟合方法；

Matlab中的多项式拟合函数；

5. 最佳平方逼近

权内积；

正交多项式的最佳平方逼近。

插值问题：

函数解析式未知，或计算复杂，用函数p(x)去近似代替它，使得

p(xi)= yi (i=0,1,2,?,n)

多项式的插值问题 构造n次多项式

Pn(x)= a20 + a1x + a2x+?+ anxn 使满足 Pn(xi)= yi (i=0,1,2,?,n)， 拉格朗日插值

插值多项式的存在唯一性： 结论

通过n+1个节点的n阶插值多项式唯一存在。 拉格朗日插值多项式的一般形式：

l(x)?(x?x0)...(x?xi?1)(x?xi?1)...(x?xn)i(xi?x0)...(xi?xi?1)(xi?xi?1)...(xi?xn)插值公式：

Pn(x)?y0l0(x)?y1l1(x)?...?ynln(x)n??y

klk(x)k?0插值的误差分析

：

分段低次插值

分段线性插值

Ih(x)?x?xk?1x?xkfk?fk?1xk?xk?1xk?1?xkn(xk?x?xk?1)若用插值基函数表示，则：Ih(x)＝?fjlj(x)j?0x?[a,b]其中：lj(x)满足lj(xk)??jk?x?xj?1??xj?xj?1??x?xj?1lj(x)???xj?xj?1?0???xj?1?x?xjxj?x?xj?1x?[xj?1,xj?1]

收敛性：

f(x)?Ih(x)?lk(x)f(x)?fk?lk?1(x)f(x)?fk?1?(lk(x)?lk?1(x))?(hk)??(h)'''其中，?(h)?maxf(x)?f(x)'''x?x?h当f(x)?C[a,b],则lim?(h)?0h?0

因此，只要f(x)?C[a,b]，limIh(x)?f(x)h?0埃尔米特插值

H2n?1(x)?a0?a1x?...?a2n?1x2n?1

mi表示yi的一阶导数

求插值基函数?i(x),?i(x),(i?0,1,...,n)满足：（）1?i(x),?i(x)是2n?1次多项式。(2)(3)于是：H2n?1(x)??[yi?i(x)?mi?i(x)]i?0n?i(xk)??i,k,且?i'(xk)?0,(i,k?0,1,...,n)?'i(xk)??i,k,且?i(xk)?0,(i,k?0,1,...,n)

利用拉格朗日插值基函数

li(x)?(x?x0)...(x?xi?1)(x?xi?1)...(x?xn)

(xi?x0)...(xi?xi?1)(xi?xi?1)...(xi?xn)得到

?n1??i(x)?1?2(x?xi)??k?0xi?xk?k?i??i(x)?(x?xi)li2(x)??l2(x)i? ??H2n?1(x)??[yi?i(x)?mi?i(x)]i?0n

Hermite插值的余项

f(2n+2)(x)2R(x)=f(x)-H2n+1(x)=wn+1(x)(2n+2)!Hermite插值也存在Runge现象

例题

步骤3、写出Hermite多项式插值结果：H5(x)?y0H0(x)?y1H1(x)?y2H2(x)?y'0H'0(x)?y'1H'1(x)?y'2H'2(x)

步骤3、写出Hermite多项式插值结果：H5(x)?y0H0(x)?y1H1(x)?y2H2(x)?y'0H'0(x)?y'1H'1(x)?y'2H'2(x)

Hermite插值结果，保留了f’(x)的特征 多项式插值结果，没有保留了f’(x)的特征

样条插值

三次样条函数是由分段三次曲线连接而成，且二阶导数连续的函数。

样条插值的构建： 由Sj(x)是[xj?1,xj]上的三次多项式，则 S''j(x)在[xj?1,xj]是线性函数（即每段是直线段）。S(x)?Mj?1''jxj-xhj?Mjx-xj?1hjhj?xj-xj-1?a(xj-x)hj?b(x-xj?1)hj

其中：S''j(xj)?MjSj(x)?Mj?1(xj-x)36hjj?0,1,...,n?Mj6hj(x-xj?1)3该公式有四个参数需要确定：Mj?1，Mj，a，b由Sj(xj)?f(xj),Sj-1(xj)?f(xj),可求出a、b，即:?Mj?1hj2?a?f(xj?1)-?6 ?2Mjhj?b?f(x)-j?6?通过上述公式，a，b可由Mj ，Mj-1计算， 上式中Mj可利用一阶导数在节点的连续性求解， 即由：S'j(xj-0)?S'j-1(xj?0)

Hilbert矩阵，非奇异，但问题病态。

数值积分

1. 插值型求积公式

求积公式的代数精度； 求积公式的误差分析； 复合求积公式； 高斯求积公式；

MATLAB中的数值积分函数。

2. 积分方程的数值求解

积分方程的数值求解的思路分析； 积分方程的数值求解方法介绍。

n次代数精度

?af(x)dx??Akf(xk)

k?0bn

对于任意不超过n次的代数多项式都准确成立, 而对某一个n+1次代数多项式不成立,。由于龙格现象，阶数越高，不稳定性越大，积分误差可能增加。

梯形公式 b b?a f(x)dx?[f(a)?f(b)]辛普森公式

?a2 b-aa?b S=[f(a)?4f()?f(b)]62：

应用高阶型插值求积公式计算积分会出现数值不稳定，

而低阶公式（如梯形、辛普生公式）又因积分区间步长过大使得离散误差大.

解决办法：缩小步长，即把积分区间分成若干子区间，

a公式，再将结果加起来，这种公式称为复合求积公式。

复化梯形公式

n?1h???f(a)?f(b)?2f(a?kh)?Tn2??

k?1??截断误差： RN(f)?复化辛普森公式

b?a?hM?， M2＝maxf??(x)

a?x?b??n?1n?1hSn?[f(a)?4?f(xk?1/2)?2?f(xk)?f(b)]

6k?0k?0截断误差： ?RN[f]??

高斯求积公式

定义

b?a2880h4M4f(?)(x) M??maxa?x?b，

求积公式 ?f(x)dx??Akf(xk) 含有2n?2个ak?0bn待定参数xk,Ak(k?0,1,?,n),适当选择这些参 数使其具有2n+1次代数精度。这类求积公式称为高斯公式。xk(k?0,1,?,n)是高斯点。

二点高斯公式具有三阶代数精度

n点高斯积分

若构造的n+1个节点的插值求积公式，则可将 f (x) = 1, x, x, …, x

2

2n+1

代入求积公式可求解，

得到的公式具有2n+1 次代数精度。这样的节点称为Gauss 点，公式称为Gauss 型求积公式。 注意: 高斯积分是不等距划分插值型求积公式

高斯点的确定方法

Matlab 积分函数

函数名 quad quad8 功能 采用Simpson计算积分。精度高，较常用 采用8样条Newton-Cotes公式计算积分。精度高，最常用 trapz 采用梯形法计算积分。精度差，速度快

积分方程的数值求解

Fredholm积分方程：by(t)???k(t,s)y(s)ds?f(t)a

用数值积分代替积分

由?k(x,s)y(s)ds??Ak(x,x)y(x)

kkkak?1bny(x)???Akk(x,xk)y(xk)?f(x)

k?1ny(xj)???Akk(xj,xk)y(xk)?f(xk)k?1n

??(I??K)y?f

常微分方程初值问题

1. 求解方法

欧拉方法；龙格—库塔方法 2. 稳定性与收敛性分析

（前向）欧拉公式：

?y(xk??)?yk???yk?hf(xk,yk)(k??,?,?,...,n??) ??xk?x??kh后向欧拉公式：

ì?y(x)?yk+1=yk+hf(xk+1,yk+1) ík+1(k=0,1,2,...,n-1)??xk=x0+kh

局部截断误差是O(h2).前向欧拉方法具有一阶精度；

同样可以证明，后向欧拉方法具有一阶精度。欧拉方法误差随着迭代次数的增加而增加。 梯形公式:

目标：求解一阶方程组：y'?f(x,y),y(x(0))?0根据牛顿-莱布尼兹公式，有： y(x(n?1))-y(x(n))??x(n?1)x(n)y'(x)dx??x(n)?hx(n)f(x,y)dx利用梯形积分公式，可得：x(n)?hh ?y'(x)dx?(f(x(n),y(n))?f(x(n?1),y(n?1)))x(n)2则：h y(n?1)?y(n)?f?x(n),y(n)??f?x(n?1),y(n?1)?2??

梯形公式几何意义：是欧拉折线法与后向欧拉法的算术平均。

精度高于前向/后向欧拉算法。

改进欧拉公式：

?预报值 ???校正值?yk?1?yk?hf(xk,yk)yk?1 h?yk?[f(xk,yk)?f(xk?1,yk?1)]2hf(xk,yk))] 2

或表示成： yk?1?yk?[f(xk,yk)?f(xk?1,yk? 平均形式：

??yp?yk?hf(xk,yk)? ??yc?yk?hf(xk??,yp)

??yk????(yp?yc)???h2 局部截断误差是O(h3). 龙格—库塔方法

一般地，RK方法设近似公式为p??yn?1?yn?h?ciKii?1???K1?f(xn,yn)?i?1?Ki?f(xn?aih,yn?h?bijKj)(i?2,3?,p)?j?1?其中ai，bij,ci都是参数，确定它们的原则是使近似公式

在(xn,yn)处的Taylor展开式与y(x)在xn处的Taylor展开式的前面项尽可能多地重合。根据泰勒定理，可表示为：h2hn(n)hn?1(n?1)y(x(n?1))?y(x(n))?hy'(x(n))?y''(x(n))?...?y(x(n))?y(?)2n!(n?1)!又由于：y'(x)?f(x,y)y''(x)?f'(x,y)...y(k)(x)?f(k-1)(x,y)则n阶泰勒法公式为：h2hn(n-1)y(x(n?1))?y(x(n))?hf(x(n),y(n))?f'(x(n),y(n))?...?f(x(n),y(n))2n!

前向欧拉方法为一阶泰勒方法，n阶泰勒方法具有n

阶精度。

例：利用二阶泰勒方法计算微分方程初值问题： y'?y-t2?1解：f(t,y)?y-t2?1,则有： f'(t,y)?y'-2t?y-t2?1-2t f''(t,y)?y'-2t-2?y-t2-2t-1利用二阶泰勒方法，得到迭代公式为：y(n?1)h2h322?y(n)?h(y(n)?t?1)?(y(n)?tn?1?2tn)?(y(n)?tn?1?2tn)26

2n目标：求解一阶方程组：y'?f(x,y),y(x(0))?0根据牛顿-莱布尼兹公式，有： y(x(n?1))-y(x(n))??x(n?1)x(n)y'(x)dx??x(n)?hx(n)y'(x)dx利用simpson积分公式，可得：x(n)?hhh ?y'(x)dx?[y'(x(n))?4y'(x(n)?)?y'(x(n?1))]x(n)62h利用f(x,y)得到y'(x(n))， 4y'(x(n)?)，y'(x(n?1))的近似公式：2

4阶龙格库塔方法,其局部误差具有四阶精度。y(x(n))-y(n)￡Ch4

本文来源：https://www.bwwdw.com/article/yafg.html