近代物理计算机模拟第03章.pdf

刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 第三章第三章 方程求根方程求根 求方程 0 xf 的一个或多个根是应用数学最常见的问题之一 大多情形下求根问题没有解析解 于是 人们试图给出求根问题的数值方法 求根问题的数值方法主要是迭代法 第一节第一节 有根区间的判定有根区间的判定 方程求根问题的数学描述是 0 xf 求x 这里 是一个函数 既没有对函 数作任何规定 也没有对求根的范围作任何限制 表面上看 这似乎是一个不太难解 决的问题 实际上 这是一个非常棘手的问题 我们必须对函数和求根的范围作一些 限制 以便能够给出一些有用的数值计算方法 xf xf xf 我们作的限制是 1 函数是连续函数 2 只求函数的实根 xf xf 假设函数在区间上连续 如果 xf ba0 bfaf 根据介值定理可知 方程 在区间内至少有一个实根 这时称区间是方程的有根区间 0 xf ba ba0 xf 在使用数值计算方法求根之前 确定所有的有根区间 特别是有一个单根或者一个重 根的区间 是非常有用的 如果我们能够确定 函数在 xf a 和 b内没有实根 那么 我们可以采用逐 步搜索的办法来确定有根区间 先将区间等分 取一个正整数 令 ban n ab h khaxk 这样 nk 1 0 就将区间等分为 等个区间 然后 在每一个小区间 上判断其是否是有根区间 这就是所谓的逐步搜索 ba 10 xx 21 xx 1nn xx n 1 kk xx 例 3 1 1例 3 1 1 设 tanf xxx 4 34 6x 求方程0 xf的有根区间 解解 分别对 4 和求函数4 3x 4 454 6 f x的值 4 32 01415212 4 41 30367622 4 50 13733205 4 64 26017490 xf x 由于函数 f x在区间 的两个端点的函数值异号 所以区间 是方程4 4 4 5 4 4 4 5 0 xf 的有根区间 48 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 例 3 1 2例 3 1 2 找出方程的区间长度不大于1的有根区间 432 44310f xxxxx 解解 2 45f 1 3f 0 1f 1 9f 2 105f 而当和时 所以 区间和 0是方程 2x 0f x 1 0 1 0 xf的有根区间 如果取的非常大 使得很小 而区间被判定是有根区间 那么 这个区间 的两个端点和中点都可以作为根的一个近似 但是 靠这种方法求根 计算量非常大 是 不合适的 nh 1 kk xx 用逐步搜索的办法来确定有根区间 有时会丢根 例如 n取的过小 使得很大 使 偶数个有单根区间包含在一个区间内 从而 该区间端点的函数值同号 该区间被误认为 不是有根区间 再比如 如果一个根是偶次重根 那么 用逐步搜索的办法是很难找到这 个根所在的区间 除非小区间的端点恰好是这个根 h 第二节第二节 试位法试位法 上一节我们讲过 用逐步搜索的办法找有根区间是个好方法 但是 用逐步搜索的办法 求根不是个好方法 试位法可以看作是它的一种改进 这里 我们介绍几种常用的试位法 二分法 黄金分割法和弦线法 一 二分法 一 二分法 设区间是方程的有根区间 即 ba0 xf0 bfaf 令 取区间 aa 0 bb 0 的中点 2 00 0 ba x 求该点的函数值 如果我们很幸运 恰巧 那么 0 xf0 0 xf 0 x 即是我们所求 这样的可能性太小了 于是 我们来看一下函数值与函数值是 否同号 如果同号 区间为有根区间 令 0 xf 0 af 00 bx 01 xa 01 bb 否则 区间为有 根区间 令 无论是哪一种情况 新的有根区间是原有根区间的一 半 将这个过程继续下去 得到了一个有根区间序列 满足 00 xa 01 aa 01 xb 11 ba kk ba 1100kk bababa 其中每一个区间都是前一个区间的一半 这就是二分法名称的由来 由于区间的长度 kk ba k kk ab ab 2 当 k时趋于零 所以这个二分过程无限继续 的结果使有根区间最后收缩为一点 这个点即是我们所求 每次二分后 取有根区间的中点 kk ba 2 kk k ba x 作为根的近似值 那么 在二分过 程中 我们得到了一个近似根的序列 该序列收敛 且收敛到根 0 x 1 x k x 在实际计算中 我们只计算有限步 对于给定的误差限 要使近似根与精确根的 k x x 49 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 误差满足 k xx 只要 1 2k ab 因为这时 1 22 k kk k abab xx 由于 b a 为已知 这样可以确定二分的次数 例 3 2 1例 3 2 1 用二分法求方程在区间内的一个实根 要求准确 到小数点后的第位 01 3 xxxf 5 1 0 1 2 解解 这里 0 1 a5 1 b01 bf 根据题的要求 005 0 2 1 k ab 即6 k时 005 0 k xx 我们将计算结果列表如下 3242 1 3203 1 3281 1 3438 1 3125 1 3750 1 2500 1 3281 1 3438 1 3750 1 5000 1 3203 1 3125 1 2500 1 0000 1 6 5 4 3 2 1 0 kkkk xfxbak 我们求得方程在区间内的一个实根的近似值是 01 3 xxxf 5 1 0 1 3242 1 二分法的优点是直观 简单 收敛性总能保证 但是 二分法收敛较慢 二 黄金分割法 二 黄金分割法 二分法的二分过程是每次将有根区间二等分 黄金分割法的分割过程也是每次分割有 根区间 与二分法不同的是 它不是等分有根区间 而是将有根区间分成两部分 一部分 是原有根区间的倍 而另一部分是原有根区间的618 0618 01 倍 设区间是方程的有根区间 即 ba0 xf0 bfaf 令 取分割 点为 求该点的函数值 检查函数值与函数值 是否同号 如果同号 区间为有根区间 令 aa 0 bb 0 000 618 0 1 618 0 bax 0 xf 0 xf 0 af 00 bx 01 xa 01 bb 否则 区间为 有根区间 令 无论是哪一种情况 新的有根区间是原有根区间的 倍或者倍 将这个过程继续下去 我们也得到了一个有根区间序列 满足 00 xa 01 aa 01 xb 11 ba 618 0618 01 kk ba 1100kk bababa 区间的长度当 kk ba kk ab k时也趋于零 所以这个分割过程无限继续的结果使 有根区间最后也收缩为一点 每次分割后 取有根区间的分割点 kk ba kkk bax 618 0 1 618 0 作为根的近似值 那么 在分割过程中 我们也得到了一个近似根的序列 该序列收 敛 且收敛到根 0 x 1 x k x 在实际计算中 我们只计算有限步 对于给定的误差限 要使近似根与精确根的 k x x 50 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 误差满足 k xx 只要 1kk xx 因为这时 1 kkk xxxx 例 3 2 2例 3 2 2 用黄金分割法求方程在区间内的一个实根 要求 准确到小数点后的第位 01 3 xxxf 5 1 0 1 2 解解 这里 0 1 a5 1 b01 bf 根据题的要求 判断 005 0 1 kk xx是否成立 如果成立 则停止计算 我们将计算结果列表如下 0040 0 0065 0 0172 0 0451 0 0730 0 1180 0 3222 1 3262 1 3197 1 3369 1 3820 1 3090 1 1910 1 3262 1 3369 1 3820 1 5000 1 3197 1 3090 1 1910 1 0000 1 6 5 4 3 2 1 0 1 kkkkkk xxxfxbak 我们求得方程在区间内的一个实根的近似值是 01 3 xxxf 5 1 0 1 3222 1 尽管在我们的例子中 黄金分割法和二分法的收敛速度相同 但大多数情形下 黄金 分割法比二分法的收敛速度快 不过仅仅快一点而已 在实际应用中 就是这一点点的快 有时会产生很大的技术革命 三 弦线法 三 弦线法 无论是二分法还是黄金分割法 每次之前都计算出了有根区间端点的函数值 但在计 算分割点时没有充分利用区间端点的函数值 弦线法利用了区间端点的函数值 期望得到 收敛速度较快的求根方法 设区间是方程的有根区间 即 ba0 xf0 bfaf 令 过点 和作一条直线 点到的线段称为曲线的 一条弦线 由此我们把这种方法称为弦线法 取分割点为弦线与 aa 0 bb 0 00 afa 00 bfb 00 afa 00 bfb xfy x轴的交点 0 00 00 00 af afbf ab ax 求分割点的函数值 检查函数值与函数值是否同号 如果同号 区间 为有根区间 令 否则 区间为有根区间 令 0 xf 0 xf 0 af 00 bx 01 xa 01 bb 00 xa 01 aa 01 xb 无论是哪一种情况 新的有根区间都比原有根区间小 将这个过程继续下去 我们 得到了一个近似根的序列 11 ba 0 x 1 x k x 在实际计算中 我们与黄金分割法相同 只计算有限步 对于给定的误差限 要使近 似根与精确根的误差满足 k x x k xx 只要 1kk xx 弦线法似乎比二分法和黄金分割法都好 事实上 它开始时收敛很快 但越来收敛越 慢 51 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 四 试位法 四 试位法 设区间是方程的有根区间 即 ba0 xf0 bfaf 令aa 0 取分割 点为 这里 bb 0 000 1 bwwax 10 w 称为试位因子 求分割点的函数值 检 查函数值与函数值是否同号 如果同号 区间为有根区间 令 w 0 xf 0 xf 0 af 00 bx 01 xa 否则 区间为有根区间 令 01 bb 00 xa 01 aa 01 xb 无论是哪一种情况 新的有 根区间是原有根区间的倍或者 11 baww 1倍 将这个过程继续下去 我们也得到了一个 有根区间序列 满足 kk ba 1100kk bababa 每次分割后 取有根区间的分割点 kk ba kkk bwwax 1 作为根的近似值 那么 在分割过程中 我们也得到了一个近似根的序列 0 x 1 x k x 在实际计算中 我们只计算有限步 对于给定的误差限 要使近似根与精确根的 误差满足 k x x k xx 只要 1kk xx 这就是所谓的试位法 在试位法中 取5 0 w 即是二分法 取618 0 w 即是黄金分割法 弦线法的试位 因子每次是不同的 与有根区间端点的函数值有关 afbf bf w 第三节第三节 迭代法迭代法 上一节讲的试位法是迭代法的一种 我们这里说的迭代法是另一类 不包括试位法 对于方程 构造其等价方程0 xf xx 给定一个方程根的猜测值 由迭代公 式 0 x 1kk xx 得到一个近似根的序列 这就是一般的迭代法 1 0 k k x 如果序列有极限 我们说这个迭代法收敛 否则 我们就说这个迭代法 发散 我们期望构造出的迭代法是收敛的 k x k k xx lim 例 3 3 1例 3 3 1 求方程在附近的根 01 3 xxxf5 1 解解 我们选择两种迭代法进行比较 取 3 1 1 xx 计算结果如下 1 3 2 xx 8 1903 396 12 3750 2 5000 1 3247 1 3247 1 3247 1 3248 1 3249 1 3259 1 3309 1 3572 1 5000 1 8 7 6 5 4 3 2 1 0 2 1 kk xxk 52 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 计算结果表明 第一种方法是收敛的 而第二种方法是发散的 那么 如何来选择迭代 法使它能够收敛就显得很重要了 假设函数 x 在上可导 我们有如下的结论 ba 如果函数 x 满足 1 对任意 bax bax 2 存在正数1pp阶收敛 当1 p 并且1 xCx 0 根据中学知识 当 时 0 a0 b abba2 而等号只有在时才能成立 由假设 我们有 ba C x C xx 2 1 0 01 假设 Cxk 我们有 C x C xx k kk 2 1 1 由数学归纳法 Cxk 0 k 另外 56 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 0 2 2 1 2 1 x 2 1 1 k kk x C xx 就收敛 并收敛到C 可见 这个迭代公式对初值的要求很松 由于 2 1 2 1 2 1 Cx x C x C xCx k kk kk 所以 CxCx Cx k k k 2 1 2 1 2 1 k 这说明 这个迭代公式是2阶收敛的 记对 k xC的相对误差为e k rk Cx x C 则有 22 1 1 1 e e 22 k rkkrk k k CxC xCxx xCxC 假设Cx 0 并且的相对误差限是 即 0 x1 0 0 0 1 r e x 那么 我们有 22 10 0 1 e e 0 10 5 10 22 rr C xx x 2 222 21 1 1 e e 0 5 10 1 25 10 22 rr C xx x 5 252 32 2 1 e e 1 25 10 7 8125 10 22 rr C xx x 11 211 2 43 3 1 e e 7 8125 10 3 051757813 10 22 rr C xx x 21 从这个例子看出 用Newton公式求一个正数的平方根是非常有效的 事实上 在实际 应用中就是这样计算一个正数的平方根的 我们再来看一个例子 例 3 4 3例 3 4 3 对于0 C 用Newton法求方程0 1 C x xf的根 解解 这实际上是要推导求一个非零数倒数的计算公式 由于 57 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 C x xf 1 则 2 1 x xf 于是迭代函数为 2 1 1 2 Cxx x C x x xf xf xx 这样 求一个非零数倒数的迭代公式为 2 1kkk Cxxx 我们来考察一下这个迭代公式 由于 2 1 1 1 2 1 C xC C Cxx C x kkkk 所以 C C x C x k k 2 1 1 1 因此 只要这个迭代公式收敛 它的收敛阶是2 下面 我们确定初值的取值范围 0 x 为了确定初值的取值范围 令 0 x kk Cxr 1 则 22 11 1 2 11 kkkkkk rCxCxCxCxr 于是就有 k rrk 2 0 如果要这个迭代公式收敛 必须要 0 k r 也就是要 1 0 r 即 11 0 Cx 58 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 所以 初值满足 0 x C x 2 0 0 m 其中 xgxxxf m 0 xg 于是 1 1 xgxxxmgxxxgxxxgxxmxf mmm 这时 迭代函数为 xgxxxmg xgxx x xf xf xx 而迭代函数的导数为 1 xhxxxhx 其中 xgxxxmg xg xh 于是 0 1 1 1 m xhx 根据收敛阶的判断规则 用Newton法求方程0 xf的重根时 仅仅是线性收敛 第五节第五节 割线法与抛物线法割线法与抛物线法 上节介绍的Newton法 在求时 不仅要给出函数值 同时要给出导数值 当函数比较复杂时 给出它的导数值往往很困难 甚至有一些函数根本就没 有导数 这时 我们希望能多用一些函数值来回避用导数值 1 k x k xf k xf xf 一 割线法 一 割线法 我们知道 Newton法从几何意义上讲就是切线法 而切线是割线的极限 所以 很容 59 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 易想到用割线来代替切线 设和是方程根的两个近似值 过点和的割线方 程为 k x 1 k x0 xf kk xfx 11 kk xfx 1 1 1 1 k kk kk k xx xx xfxf xfy 此割线方程的根为 1 1 1 12kk kk k kk xx xfxf xf xx 这就是割线法的迭代公式 这个迭代公式要求给出两个初值和 它从形式上看 很像是弦线法 但是 他们 有很大的不同 割线法每次只管迭代计算 依次求出 而弦线法则要判断函数 值的正负 从而选择下次迭代谁参加运算 0 x 1 x 3 x 4 x 例 3 5 1例 3 5 1 用割线法求方程在附近的根 01 3 xxxf5 1 解解 可以取1 4和1 5做初始值 也可以取1 5和1 6做初始值 这里 先取1 4和1 5做初 始值 令 0 1 4x 1 1 5x 1 1 1 12kk kk k kk xx xfxf xf xx 1 计算 得 0k 012345 1 40001 50001 33521 32621 32471 3247 k k x 再取1 5和1 6做初始值 令 0 1 5x 1 1 6x 1 1 1 12kk kk k kk xx xfxf xf xx 1 计算得 0k 0123456 1 50001 60001 35911 33201 32491 32471 3247 k k x 从这个例子可以看出 割线法比Newton法要慢 割线法与Newton法一样 也是局部 收敛的 它的收敛阶为618 1 2 15 p 我们可以从另一个角度来认识割线法 设和是方程根的两个近似值 用点和作一次多 项式 其函数形式就是过这两点的割线方程 用这个一次多项式的根来近似方程 k x 1 k x0 xf kk xfx 11 kk xfx 0 xf根 这样 同样可以推导出割线法 把这种思想用到三个点的情形 就可推出抛物线法 60 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 二 抛物线法 二 抛物线法 设 和是方程 k x 1 k x 2 k x0 xf根的三个近似值 过 和 三点作二次多项式 得到 kk xfx 11 kk xfx 22 kk xfx 1 21 1 kk kk kk f xf x pxf xS xxxx xx k 其中 211 211 2 kkk kkkk kk k f xf xf xf x xxxx S xx 用这个二次多项式的根来近似方程0 xf根 2 32 2 2 2 4 k kk k f x xx Sf x 这里 21 21 21 kk kk kk f xf x S xx xx 设 2k x 更接近根 取正负号使 3k x 较接近 2k x 所以 根号前取与 相同的符号 这样就推导出一种新的方法 因为二次多项式的曲线代表一条抛物线 所以这样推导 出的方法叫做抛物线法 抛物线法也是局部收敛的 它的收敛阶为1 840p 第六节第六节 多项式方程求根多项式方程求根 在本节中 我们将函数限定为多项式 即为多项式 求方程的根 也就是要求多项式方程的根 xf xf0 xf 一 多项式的求值 一 多项式的求值 设给定多项式 n nx axaaxf 10 千万不能将系数的值和自变量的值直接代入上式来求此多项式的值 这样做 一方面 计算工作量很大 另一方面可能会产生很大的误差 甚至不能进行计算 例如 通常计算 机做数值计算时 只能算以内的数 超过这个数叫做溢出 计算机会停止运算并指出 溢出错误 而当计算多项式最高项时取10 38 10 nx取5000就会产生溢出 下面我们推导一个有效的计算多项式值的方法 用除多项式 记商为 余数显然等于 即有 0 xx xf xg 0 xf 00 xfxgxxxf 61 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法计算机数值方法 第三章第三章 方程求根方程求根 设 1 110 n n xbxbbxg 代入上式 比较两边同次幂的系数 得 1 1 00 0 00 1 1 ni xb xb xfa ba ba iii nn 从而有 0 2 000 01 0 1 1 ni xba xb xf ab ab iii nn 这种计算多项式函数值的有效算法称为秦九韶法 62 4 例 3 6 1例 3 6 1 已知 23 1f xxxx x 用秦九韶法计算的值 0 5 f 解解 根据题意 0 1a 1 1a 2 1a 3 1a 4 1a 0 0 5x 由秦九韶法得 34 1ba 2330 1 0 50 5bab x 1220 1 0 250 75bab x 0110 1 0 3750 625bab x 0000 1 0 31250 6875f xab x 二 多项式的求导 二 多项式的求导 用同样的方法可以推出多项式导数值的计算公式 因为 00 xfxgxxxf 所以 0 xgxxxgxf 00 xgxf 用除多项式 记商为 余数显然等于 即 于是有 0 xx xg xh 0 xg 0 xf 00 xfxhxxxg 设 2 210 n n xcxccxh 代入上式 比较两边同次幂的系数 得 1 2 00 0 00 1 21 ni xc xc xfb cb cb iii nn 从而有 刘刘 悦 近代物理计算机模拟悦 近代物理计算机模拟 第一部分第一部分 计算机数值方法