非线性方程最终版
非线性方程求根
几何含义：求曲线 y = (x) 与直线 y = x 的交点
18:28:03
Numerical Analysis
9
y p0 x0 y
18:28:03
y=x
p1
y= (x)
✓
y p0
y=x
✓
p1 y= (x)
x x1 x2 x*
x0
(x) ( y) L x y
则(x) 在 [a,b] 上存在唯一的不动点 x*
证明：P216
18:28:03
Numerical Analysis
12
收敛性分析
不动点迭代的收敛性
定理：设 (x) C[a,b] 且满足
(1) 对任意的 x[a,b] 有 (x)[a,b]
(2) 存在常数 0<L<1，使得任意的 x, y[a,b] 有
for k = 1 : n x = g(x); fprintf('k=%2d, x=%.7f\n',k,x); if abs(x-xt)<tol, break, end
end xt = fzero(f,[3,4]);
fprintf('True solution: x = %.7f\n', xt)
% Steffenson 加速
性质：若
lim
k
xk
x *，则不动点迭代收敛，且 x*
是 f(x)=0 的解；否则迭代法发散。
18:28:03
Numerical Analysis
11
解的存在唯一性
解的存在唯一性
定理：设 (x) C[a,b] 且满足
非线性方程组的解法
非线性方程组的解法
基本思路：
分段线性化方法，将荷载划分成很多小步，逐 步施加
具体操作方法：
显式求解法（增量法） 隐式求解法（迭代法）
清华大学研究生课程——《钢筋混凝土有限元》
清华大学研究生课程——《钢筋混凝土有限元》
显式求解法
将荷载分成若干小步，逐步施加 认为在每个小步中，结构是线性的，同 一荷载步的刚度矩阵相同 不同荷载步的刚度矩阵可以不同 用一系列的折线去近似曲线
某次迭代位移改变量 同级荷载节点总位移 < 误差容限
能量收敛标准
一般以某次迭代的应变能增量为分析对 象，以同级荷载作用下总应变能为参考 标准
某次迭代应变能改变量 同级荷载总应变能 < 误差容限
一般使用无穷范数
清华大学研究生课程——《钢筋混凝土有限元》
清华大学研究生课程——《钢筋混凝土有限元》
注意事项
修正的欧拉折线法(Mid-point Method)
P K2 K1 Pn Pn-1
单元刚度矩阵
已知应力，应变，应变增量
[K n−1 ]
根据当前应力应变求切线刚度矩阵 求中点应力 {σ n '} = {σ n−1 } + [K n−1 ] {dε n }
2
根据中点应力和应变 {σ n '}， {ε n−1 } + n 求此时 2 的切线刚度矩阵 [K n−1 ']
判断对象
力收敛标准 位移收敛标准 能量收敛标准
范数
无穷范数 一范数 二范数
V
V
V
∞
= max Vi
= ∑ Vi
=
判断标准
相对误差 绝对误差
1
范数
无穷范数 一范数 二范数
非线性方程求解
非线性方程求解在数学中，非线性方程是一种函数关系，其表达式不能通过一次函数处理得到。
与线性方程不同，非线性方程的解决方案往往更具挑战性，因为它涉及到更复杂的计算过程。
尤其在实际应用中，非线性方程的求解是一个非常重要的问题。
本文将讨论几种常用的非线性方程求解方法。
二分法二分法，也称为折半法，是一种基本的求解非线性方程的方法之一。
它的核心思想是将区间一分为二并判断方程在哪一半具有根。
不断这样做直到最终解得精度足够高为止。
下面是利用二分法求解非线性方程的流程：1. 设定精度值和区间范围2. 取区间的中点并计算函数值3. 如果函数值为0或函数值在给定精度范围内，返回中点值作为精确解4. 如果函数值不为0，则判断函数值的正负性并缩小区间范围5. 重复步骤2-4直到满足给定精度为止当然，这种方法并不总是能够找到方程的解。
在方程存在多个解或者区间范围不合适的情况下，二分法可能会导致求解失败。
但它是一种很好的起点，同时也是更复杂的求解方法中的一个重要组成部分。
牛顿迭代法牛顿迭代法是一种更复杂的求解非线性方程的方法。
它利用泰勒级数和牛顿迭代公式，通过不断迭代来逼近根的位置。
下面是利用牛顿迭代法求解非线性方程的流程：1. 先取一个近似值并计算函数值2. 求出函数的导数3. 利用牛顿迭代公式，计算下一个近似根4. 检查下一个近似根的精度是否满足条件，如果满足，返回当前近似根5. 如果精度不满足，则将新的近似根带入公式，重复步骤2-5当然，牛顿迭代法的收敛性并不总是保证的。
如果迭代过程太过温和，它可能无法收敛到精确解。
如果迭代过程过于暴力，则会出现发散现象，使得求解变得不可能。
其他方法此外，还有一些其他的求解非线性方程的方法，例如黄金分割法、逆二次插值法、牛顿切线法等等。
其中每一种方法都有其优缺点，不同的情况下，不同的方法都可能比其他方法更加适合。
结论总体来说，求解非线性方程的方法非常复杂。
无论是哪种方法，都需要一定的数学基础和计算机知识。
常见非线性方程
杜芬方程来源：
无阻尼单摆
x sin x 0
有阻尼单摆
x 2x sin x 0
sin x x
杜芬方程
1 3 x 3 ！
x x x x 3 F cos t
化的说法就是所谓的蝴蝶效应。
1979年12月，洛伦兹（Lorenz）在华盛顿的美国科学促进会的再 一次讲演中提出：一只蝴蝶在巴西扇动翅膀，有可能会在美国的德克 萨斯引起一场龙卷风。他的演讲和结论给人们留下了极其深刻的印象。 从此以后，所谓“蝴蝶效应”之说就不胫而走，名声远扬了。
他说，一只南美洲亚马逊河 流域热带雨林中的蝴蝶，偶尔扇 动几下翅膀，可能在两周后在美 国德克萨斯引起一场龙卷风。其 原因在于：蝴蝶翅膀的运动，导 致其身边的空气系统发生变化， 并引起微弱气流的产生，而微弱 气流的产生又会引起它四周空气 或其他系统产生相应的变化，由此引起连锁反映，最终导致 其他系统的极大变化。洛伦兹把这种现象戏称做“蝴蝶效应” ， 意思即一件表面上看来毫无关系、非常微小的事情，可能带 来巨大的改变。
O
l 2 m x mgl sin x
有阻力时， 阻力与摆的速度 lx 成正比：
l
x
l 2 m x mgl sin x lx
取 2m ，得： g x 2 x sin x 0 l g 2 n l 2 x 2 x n sin x 0
2 n 1
m
x 2 x sin x 0
讨论非线性项
sin x x
1 3 1 5 1 7 x x x 3！ 5！ 7！
线性振动方程
第2章非线性方程求根
解：设最多需要迭代n次。
∵要求精确到小数点后3位，
∴误差限≤ 12×10-3，
∴由定理2.1得：
n=
lg(1
0)
lg( lg 2
1 2
103 )
1
=
0
( lg 2 lg 2
3)
1
=
0.301 0.301
3
1= 10 .97
1
=10，即最多需要迭代10次。
用对分法求f(x)=0在某区间的单根，最多迭代次数与函数f(x)曲线形状无 关。一般情况下，对分法的最多迭代次数比其他的变步长逐步搜索法要 少，因此对分法是用得最多的变步长逐步搜索法。
if(f(b)==0) x=b;
else for(begin=a,end=a+h;begin<b;begin=end,end+=h)
{
if(end>b)end=b;
if(f(begin)==0)
{x=begin;break; }
if(f(begin)*f(end)<0)
{x=(begin+end)/2;break; }}
Y
N
x=begin;
break;
f(begin)*f(end)<0
Y
N
(end－begin)/2<=ε
Y
N
x=(begin+end)/2;
h/=hnumber;
break;
end=begin;
输出方程f(x)=0的根x。
8
2.3 根的搜索
变步长逐步搜索法对应的程序
#include <stdio.h> double f(double x);
非线性方程
§4 牛顿法 /* Newton - Raphson Method */
原理：将非线性方程线性化 —— Taylor 展开 /* Taylor’s expansion */ 取 x0 ≈ x*，将 f (x)在 x0 做一阶Taylor展开: f ′′(ξ ) f ( x) = f ( x0 ) + f ′( x0 )( x − x0 ) + ( x − x0 )2 ，ξ 在 x0 和 x 之间。
an x n + an −1 x n −1 + + a1 x + a0 = 0 n ≠ 0) (a
为 n 次代数方程，当 n > 1 时，方程显然是非线性的。
ZUCC
浙江大学城市学院
根的计算步骤 判定根的存在性。即方程有没有根？如果有，有几个根？
思 可以利用闭区间上连续函数的零点定理解决。 路
ZUCC 浙江大学城市学院
§2 Bisection Method
Algorithm: Bisection
Find a solution to f (x) = 0 given the continuous function f on the interval [ a, b ], where f (a) and f (b) have opposite signs. Input: endpoints a, b; tolerance TOL; maximum number of iterations Nmax. Output: approximate solution x or message of failure. Step 1 Set i = 1; FA = f (a) ; Step 2 While ( i ≤ Nmax) do steps 3-6 Step 3 Set x = ( a + b ) / 2 ; /* compute xi */ FX = f (x) ; Step 4 If ( FX == 0 ) or ( b − a ) / 2 < TOL then Output (x); STOP; /* successful */ Step 5 Set i ++; Step 6 If FA·FX > 0 then set a = x ; FA = FX ; Else set b = x ; /* update ai, bi */ Step 7 Output (The method failed after Nmax iterations); /* unsuccessful */ STOP.
非线性方程组的解法
非线性方程组的解法
非线性方程组的解法包括：
（1）近似法。
近似法是根据所给非线性方程组，使用一定的数值方法，建立非线性方程组结果的拟合曲线，以此求解非线性方程组的常用方法，目前有贝塔、拉格朗日近似法和微分近似法等。
（2）多元分割法。
多元分割法根据非线性方程组的参数和变量空间，
将整个运算范围分割成多余小区间，利用各区间中只含有一个未知变
量的简单方程组，将非线性方程组转换成多个一元方程组，再用一次法、弦截法和二分法等算法求解，最终得出整个非线性方程组的解。
（3）迭代映射法。
迭代映射法是通过给定一个初始值，然后利用迭代，反复运算，最终达到收敛点的一种方法，主要包括牛顿法、收敛法、
弦截法、松弛法和隐函数法等。
（4）最小二乘法。
最小二乘法是将非线性方程组表示为残差函数，然
后求解残差函数最小值，获得未知变量的最优解，常用于数值分析中。
（5）特征法。
特征法是采用将非线性方程组表示为线性方程组特征值
和它们关于某一特征量的关系式，利用梯度下降法，最小化残差函数，求解非线性方程组的方法。
以上是非线性方程组的解法的简单综述，它们在一定程度上增加了解决非线性方程组的效率，但并非所有情况都能使用以上求解方法。
正确使用相应的求解方法就可以有效的求解非线性方程组，以便更好的解决实际问题。
一类高阶非线性差分方程最终正解的存在性
数 Ⅳ， 使得
由式 （ 又有 ４）
ｃ … 吉 一 薹
当 Ⅳ 一 ≤ ＜Ⅳ时 ， 二
×
茎 ～ 。 吉 ㈤
定 的 集 ｛∈ ： ≤， 义 子 ＝ ≤ ｌ
≥ ） 则 是 的 个 界 凸 集定 算 ０。 力 Ｅ 一 有 闭 子 ，义 子
： Ｅ如下
高 阶非线 性差 分方 程正 解 的存在 性 问题 。
考 虑方程 △ （ ＋ｐ ｘ一）＋ｇ 一） ：０ √ （） １
定理 １ 假设 条 件 （ 成 立 ， 有 Ｈ） 且
∑ｓ１ ＜∞ ｍｑ －
则 方程 （ ） 在有界 的最终正 解 。 １存
（） ３
证 令 ｍ ） ≤ ≤ ｐ， 明 ＝ ａ ： ｘ 二 一）
（） ３ ｆ∈ Ｃ Ｒ， （ Ｒ）， 满足 ≠ ０时 ｘ（ ）＞０， ｆｘ 且
存 在 李普 希兹 常数 Ｌ ＞０， 使
ｌ ）一 （ ） ≤Ｌ 一Ｙ ， ， － 厂 （ －Ｙ ｌ 厂 Ｉ Ｖ Ｙ∈Ｒ （ ） Ｉ ２
实序 列 ｛ ： 成 为差 分 方 程 （ ） Ｙ｝ 一 １ 的一 个解 ， 若 当 ｎ≥０时 ， 此序 列 的各 项满 足方 程 （ ） １ 。称 为它
关键词
高 阶差分方程 Ｏ７ ． ； １５ ７
存在性
最终正解 Ａ
中图法 分类 号
文献标志码
由 于计 算 机 科 学 、 物 数 学 、 代 物 理 等 自然 生 现 科学 与边 缘科 学 的迅 速 发 展 ， 时滞 差 分 方程 稳 定 对 性 理 论 的 研 究 近 年 来 十 分 活 跃 。一 些 作 者 对 低
） ｎⅣ ， ≥
（ ）， Ⅳ Ⅳ 一 ≤ ｎ＜Ⅳ （） ７
数学中的非线性方程求解
数学中的非线性方程求解非线性方程是指未知量与其函数之间不满足线性关系的方程。
解决非线性方程的问题一直是数学领域的研究重点之一，因为非线性方程在自然科学、工程技术以及金融经济等领域中具有广泛的应用。
在本文中，我们将探讨几种常见的非线性方程求解方法。
一、二分法二分法也称为区间二分法，是求解非线性方程最基本的方法之一。
该方法利用非线性方程连续性的特点，将方程的解所在的区间不断二分并缩小区间范围，最终找到非线性方程的解。
考虑一个一元非线性方程f(x)=0，其中f(x)在区间[a, b]上连续且f(a)与f(b)异号。
根据区间中值定理可知，存在一点c属于(a, b)，使得f(c)=0。
我们可以按以下步骤进行二分法的求解：步骤1：选择区间[a, b]，计算函数值f(a)与f(b)。
步骤2：如果f(a)与f(b)异号，则继续进行下一步。
否则，结束计算，方程无解。
步骤3：计算区间中点c=(a+b)/2，并计算f(c)。
步骤4：如果f(c)接近于0或满足终止条件，则c为方程解。
否则，根据f(a)与f(c)的符号确定新的区间[a, c]或[c, b]。
步骤5：重复步骤3和步骤4，直至满足终止条件。
二、牛顿法牛顿法是一种迭代逼近的方法，通过使用函数的一阶和二阶导数来逼近非线性方程的解。
该方法基于泰勒级数展开，通过不断迭代逼近函数零点的位置。
考虑一个一元非线性方程f(x)=0，我们可以按以下步骤进行牛顿法的求解：步骤1：选择一个初始近似值x0。
步骤2：计算函数f(x)的一阶导数f'(x)和二阶导数f''(x)。
步骤3：使用初始近似值x0和函数导数来进行迭代计算，得到新的近似值x1。
迭代公式为x1 = x0 - f(x0)/f'(x0)。
步骤4：计算函数f(x1)的值。
步骤5：如果f(x1)接近于0或满足终止条件，则x1为方程解。
否则，将x1作为新的近似值，重复步骤3和步骤4。
步骤6：重复步骤3至步骤5，直至满足终止条件。
第七章非线性方程.ppt
-2.5 -3
s in(x )2 ex p (-01 x )-05 abs (x )
-6
-4
-2
0
2
4
6
作业：
❖ 习题2
C 2 C .
2k
xk xk
C C
x0 x0
C C
记 q x0 C , x0 C
q2k xk C 2 C 1 q2k .
对任意x0 0,总有 q 1，故当k 时，xk C.
Newton迭代法的缺陷 1.被零除错误
方程: f(x)=x3 –3x + 2 = 0
在重根x*=1附近,f’(x)近 似为零
f(x) 0: f(x0 ) f '(x0 ) xx0 0
解出x作为近似根x1：
x1 x0 f( x0 )
f
'
(
x0
)
(
f
'
(
x0
)
0)
依次产生迭代格式，称 Newton 法：
x x x f x'
f ( ) ( ), n 0, 1, 2,
n1
n
n
n
Newton 法的几何解释
当
x0 在取后（在真根
（1）x x 1 1(x) 因 1.5 1.5 1 1(x) 2 1 2
且
1' (x) 2
1 1 1 x 1 2 2.5 3.162
(2)
x 1 1 x
2(x)
因
1.5
1
1 2
2
(x)
1
1 1.5
2
且
2' (x)
1 x2
1 1.52
非线性方程
方程根的数值计算步骤
判断根的存在 确定根的分布范围 根的精确化
计 根的存在定理(零点定理)：
f(x)为[a,b]上的连续函数，若 f(a)· f(b)<0，则 方 [a,b]中至少有一个实根。如果f(x)在[a,b]上还是单 法 调递增或递减的，则f(x)=0仅有一个实根。
算
计 算
根的分布范围:
k 1 k
方 法
取初值
x0 0
3 0
x1 2 x 1 1
x2 2x 1 3 ⇒
3 1
同样的方程
不同的迭代格式 有不同的结果
什么形式的迭代 法能够收敛呢? 如何构造迭代函 数呢?
x3 2 x 1 55
3 2
由此可见，这种迭代格式是发散的
若从任何可取的初值出发都能保证收敛，则称它 算 为大范围收敛。如若为了保证收敛性必须选取初值充 方 分接近于所要求的根，则称它为局部收敛。
h=(b-a)/n 其中n是正整数，在[a,b]内取定节点： xi=x0＋ih
(i=0,1,2，……，n)
计算f(xi)的值,依据函数值异号及实根的个数确 定隔根区间,通过调整步长，总可找到所有隔根 区间。
计 对于m次代数方程 f (x) = xm+am-1xm-1+ …+a1x+a0=0其根的 模的上下界有如下结论： 算 (1)若μ= max { |am-1| , ……, |a1| , |a0| }，则方程根的模小于μ+1 方 1 1 …… 法 (2)若ν= max {1， |am-1| , , |a1| },则方程根的模大于 1 a1 xk ε1
或
f ( x ) ε2
计 误差 分析： 算 方 法 第1步产生的
非线性方程
2 、 迭 代 法
基本思想：
使用某个固定公式，反复校正根的近似值， 使用某个固定公式 ， 反复校正根的近似值 ， 从而得到一个近似根的序列 从而得到一个近似根的序列 { xk }，使得该序列的 ， 极限就是方程的根。 极限就是方程的根。 就是方程的根
具体办法：
1、 f ( x ) = 0 、
⇔
x = g( x )
( k = 0,1L)
2、取初值 x0 ，使得 x k +1 = g ( x k ) 、
3、如 g(x) 连续，且序列 { xk } 收敛到 x*，则有 、 连续， 收敛到 ，
x∗ = lim x k +1 = lim g ( x k ) = g (lim x k ) = g ( x ∗ )
精确值 快速割线法三次迭代所求解 一般割线法三次迭代所求解 割线法与切线法的收 敛条件相似。 敛条件相似。
0
x4 x3
x2
x1
x0
x
编程的基本步骤： 割线法 编程的基本步骤：
Step 1：give x0, x1, eps1, eps2, function f(x) ：
注：循环计算时，常设置一最大循环次数。 循环计算时，常设置一最大循环次数。
离 散 化 逼 近 推 递
f ( x 0 ) + f ' ( x 0 )( x − x0 ) = 0
⇒
f ( x0) x = x0 − f '( x0)
第十章非线性方程及非线性方程组解法
例5 设
f1 ( x, y ) x 2 y 2 5 0 f 2 ( x, y ) ( x 1) y (3x 1) 0 用牛顿法求( x0 , y0 ) (1 ,1)附近的解
解：先计算出偏微商矩阵
f1 x f 2 x f1 2y y 2 x f 2 y 3 x 1 y
第十章 非线性方程及非线性方程组解法
在许多实际问题中常常会遇到求解非线性方程
或非线性方程组的问题。
例如，求n次代数方程
an x n an1 x n1 a1 x a0 0
的根，或求超越方程。
x
e
sin(
X
2
)0
的根，这些都可以表示为求f(x)=0的根，或称为 求函数f(x)的零点。
x 2 x1
x1 x0
f ( x0 )
f ( x0 )
再把 x1 作为节点，利用泰勒展式，令f ( x1 ) 0，
f ( x1 )
f ( x1 )
如此下去，得到牛顿法的迭代序列
x n 1 x n f ( xn ) f ( x n )
牛顿法的几何意义 切线方程：y f ( x0 ) f ( x0 )(x x.0 )
(10-24）
只要系数矩阵的行列式
f1 ( x0 , y0 ) x J f 2 ( x0 , y0 ) x f1 ( x0 , y0 ) y 0 f 2 ( x0 , y0 ) y
则方程组（10—24）的解可以写成
x1 x0 1 J0 1 y1 y0 J0 f1 ( x0 , y0 ) y f 2 ( x0 , y0 ) y f1 ( x0 , y0 ) f 2 ( x0 . y0 ) f1 ( x0 , y0 ) f 2 ( x0 , y0 ) f1 ( x0 . y0 ) x f 2 ( x0 , y0 ) x
非线性光学习题解答终极版
非线性光学习题解答[李春蕾(2011111772)]第一章 晶体光学简介 电光效应1.解答：由于矢量运算不受坐标系的影响，只是表示形式不同而已，不妨在直角坐标系下建立方程，设x x y y z z k k e k e k e =++ ，x y z r xe ye ze =++，x y z e e e x y z∂∂∂∇=++∂∂∂ 则exp[i()]x y z fk x k y k z =++ 于是i i i i x y z x x y y z zf f f fe e e fk e fk e fk e fk x y z∂∂∂∇=++=++=∂∂∂，问题得证。
对于平面波，设0000exp[i()]()exp[i()]x x y y z z EE t k r E e E e E e t k r ωω=−⋅=++−⋅其中，0exp[i()]x x E E t k r ω=−⋅，0exp[i()]y y E E t k r ω=−⋅ ，0exp[i()]z z E E t k r ω=−⋅()()()(i i )(i i )(i i )i x y zy y x xz z x y z x y z x y z z y y z x x z z x y y x E e e e E E E E E E e e e x y z y z z x x y E E E e k E k E e k E k E e k E k E k E∇×∂∂∂∂∂∂∂∂∂==−+−+−∂∂∂∂∂∂∂∂∂=−++−++−+=×同理，i H k H ∇×=×.2.证明：在选定主轴坐标系的情况下，物质方程可以写成0i i i D E εε=，1,2,3i =同时，将晶体光学第一基本方程写成分量形式，20[()]ii i D n E k k E ε=−⋅，1,2,3i =联立两式，整理得到02()11ii i k k E D nεε⋅=− 对于对应同一个k 的两个电位移矢量D ′ 、D ′′，建立它们的标量积2222312022222211223322210221()()111111111111()()()()()()()()()()()1D D k k k k E k E n n n n n n k n n k E k E n n εεεεεεεεε′′′⋅ ′′′=⋅⋅++−−−−−−′′′′′′′′′′′′′′′=⋅⋅′′′− 22222331222222221223311111111111()()()()()()k k k k k n n n n n n εεεεε −+−+−−−−−−− ′′′′′′′′′由23122011i i ik n n ==−∑，得到大括号中的第一、三、五项之和为零，第二、四、六项之和为零，所以0D D ′′′⋅=即对应同一个k 的两个电位移矢量D ′ 、D ′′相互垂直.3.解答：22011i i k nε=−∑是方程20[()]i ii i E n E k k E εε=−⋅的本征值方程，设其本征值为m n ，相应的本征解为()m E ，则可以得到，()2()()0[()]m m m mEn E k k E ε⋅=−⋅ε晶体中可以有两个本征解，设另一个为()n E ，用其点乘上式得到()()()()()201[()]n m m n m mEE k k E E E n ε⋅−⋅=⋅⋅ε 交换指标m 和n 后可以得到()()()()()201[()]m n n m n nEE k k E E E n ε⋅−⋅=⋅⋅ε 上两式相减，考虑到介电常数张量ε为对称张量，则可以得到()()2201110m n n m E E n n ε −⋅⋅=ε 如果m n n n ≠，则有()()0m n EE ⋅⋅=ε 如果m n n n =，显然方程成立。
第二章非线性方程(精)
95第二章 非线性方程求根代数方程求根问题是一个古老的数学问题，早在16世纪时找到了三次、四次方程的求根公式。
但直到19世纪才证明5≥n 次的一般代数方程式不能用代数公式求解。
因此需要研究用数值方法求得满足一定精度的代数方程式的近似解。
在工程和科学技术中许多问题常常归结为求解非线性方程式问题，例如在控制系统的设计领域，研究人口增长率等。
本章将介绍几种求解非线性方程近似解的数值方法。
最后一节简单介绍求解非线性方程组的方法。
§1 二 分 法设有非线性方程0)(=x f (1.1)其中，)(x f 为],[b a 上连续函数且设0)()(<b f a f （不妨设方程（1.1）于],[b a 内仅有一个实根）。
求方程（1.1）实根*x 的二分法过程，就是将含根区间],[b a 逐步分半，检查函数值符号的变化，以便确定含根的充分小区间。
二分法叙述如下：记b b a a ==11,。
第1步分半计算)1(=k ：将],[11b a 分半，计算中点2/)(111b a x +=及)(1x f ，如果 0)()(11<x f a f则根一定在区间],[],[2211b a x a ≡内，否则根一定在区间],[],[2211b a b x ≡内（若,0)(1=x f 则*=x x 1）。
于是得到长度缩小一半的含根区间],[22b a ，即0)()(22<b f a f ，且)(211122a b a b -=- 第k 步分半计算：重复上述过程，设已完成第1步，...，第1-k 步分半计算得到含根区间96 ],[],[2211b a b a ⊃],[...k k b a ⊃⊃ 且满足：（1）0)()(<k k b f a f ， 即],[k k b a x ∈*； （2）)(211a b a b k k k -=--；现进行第k 步分半计算： （3） 计算2/)(k k k b a x += 且有 )(212a b a b x x k k k k -=-≤-*（1.2） （4） 确定新的含根区间],[11++k k b a ，即如果0)()(<k k b f a f ，则根一定在],[],[11k k k k x a b a =++内，否则根一定在区间],[],[11k k k k b x b a =++内，且有)(2111a b a b k k k -=-++总之，由上述二分法得到一序列}{k x ，由（1.2），则有 *∞→=x x k k lim可用二分法求方程0)(=x f 实根*x 的近似值到任意指定的精度。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
非线性方程求解化11 叶亦欣 2011011786问题1：两个购房贷款的问题（1） 小张夫妇以按揭方式贷款买了一套价值20万元的房子，首付5万，每月还款1000元，15年还清，问贷款利率是多少？（2） 某人欲贷款50万元购房，有两家银行可供选择：第一家每月还4500元，15年还清，第二家每年还45000元，20年还清。
从利率方面看，哪家更优惠？（假设年利率=月利率*12）【数学模型】可以看到，两个问题都是每个月还款数量相同（（2）中第二家是按年，但道理相同），这称作等额本息还款法。
计算方法是每月利息按月初剩余欠款本金乘以月利率，再加上一定数额的本金作为每月还款总额，因此每个月所还本金和利息的数额随时间改变（但总额不变），这样到设定期限使本金恰好还清。
设贷款本金数额为A ，银行月利率为a ，一共还n 个月，每月还款额为x 。
则，第1个月后剩余欠款本金为x a A −+)1(；第2个月后剩余欠款本金为)1)1((2)^1()1)()1((++−+=−+−+a x a A x a x a A ； 第3个月后剩余欠款本金为)1)1(2)^1((3)^1()1)()1)()1(((++++−+=−+−+−+a a x a A x a x a x a A ；……第k 个月后剩余欠款本金为i a x k a A k i ^)1()^1(1∑−=+−+；而到第n 个月后，剩余欠款本金恰好为零，即0^)1()^1(1=+−+∑−=i a x n a A n i ，整理得0)^1)((=++−x n a x Aa ，其中a 为此单变量非线性方程的未知量，将各值代入即可求解单位时间的利率 【解答】（1） 小张夫妇贷款的月利率计算：将小张夫妇贷款的具体数额等值代入等额本息还款法的数学模型，得到如下方程：01000180)^1)(1000150000(=++−a a对于这个等额本息还款模型，显然a=0是方程的一个根，而我们所需要的是[0,1]上的另一个根，因此不能以0作为求解的初值，所以我们要先将方程左边的函数图象画出来以便有一个基本认识，然后来确定解的大致范围。
由常识可知，月利率一般不会超过1%，因此可以只看[0,0.001]区间上的图象[1]：发现我们需要的零点在(0，0.003]区间，但由于范围太大，无法从图中确定我们所需要的左面的闭区间节点，因此缩小范围：发现我们需要的根大致落在[0.002,0.0025]区间内，之后用本章所学的fzero命令求解方程的根[2]。
得此种情况下月利率约为0.208%（3）某人比较两家贷款银行的利率的问题：与（1）求解过程基本相同，将两种情况下的实际值分别代入模型，先画出比较大范围内的方程左边函数的图像以确定解的大致范围，再解出较精确的根：第一家银行：左端函数图像：发现我们所需要的根大致落在[0.005,0.007]区间，再用fzero命令求解，得这家银行月利率约为0.585%；第二家银行：左端函数图象（这家采用年利率，因此图象显示的区间相应应该扩大，变成[0,0.1]）：发现我们所需要的根大致落在[0.06,0.07]区间，再用fzero命令求解，得这家银行年利率约为6.395%；由题中假设，认为年利率=月利率*12，计算发现0.585*12=7.02>6.395，因此认为第二家银行更优惠。
【讨论】1．年利率=月利率*12的假设是否合理？假如第一家银行实际上是用年利率计算的，那么每年还款额为54000元，代入模型方程求解得其年利率为6.74%，仍然大于第二家，但是不如按照年利率=月利率*12算的多。
考虑到按常理说，按月支付小额欠款比按年支付较大额度对贷款者压力更小，折合成年利率时应该有所增加，而且也并没有影响比较结果，认为这一假设基本合理。
【命令代码】[1]fplot('(150000*a-1000)*(1+a)^180+1000',[0,0.01]);[2]function [y] = rate( a )y=(150000*a-1000)*(1+a)^180+15000；end %非线性方程函数x=fzero(@rate,[0.0018,0.0022]) %求解问题2由汽缸控制关闭的门，关闭状态的示意图如下左图，H为固定的可转动的点，与H相距b处有一门销，通过活塞与圆柱形汽缸相连，活塞半径r，汽缸长l0，汽缸内气体压强为p0。
当用力F推门，使门打开一个角度α（如下右图），活塞下降的距离为c，门销与H的水平距离b保持不变，于是汽缸内的气体被压缩，对活塞的压强增加。
在绝热条件下，气体的压强p和体积V满足过程方程pVγ=q,其中γ是绝热系数，q是常数。
利用开门的力矩和作用在活塞上的力矩相平衡的关系，求在一定的力F作用下，门打开的角度α。
a=0.8m,b=0.25m,r-0.04m,l0=0.5m,p0=10^4N/m^2,γ=1.4,F=25N.【数学模型】对于H点，F产生的力矩为F∗cosα×a，气体压强通过活塞产生的力矩为p∗πr2∗cosα×(b)，使两者相等则有p=F∗a∗cosαb∗πr,其中 cosα=b√22.另一方面，由绝热状态方程有p0V0γ=pVγ,其中V0=l0∗πr2,V=(l0−c)∗πr2，将p 代入整理得（1−c0）γ−p0∗πr2∗�b2+c2=0此方程为单变量非线性方程，未知量为c。
进而由cosα=b√22可求出cosα值。
【求解】将题中所给数值代入方程得下式：(1−2c)1.4−0.8π∗�0.0625+c2=0,我们知道0一定不是根，而且c的有效取值一定在（0，0.5]区间，因此先将左端函数在这一区间上的图像画出来以确定迭代初值[1]：发现根在0.1附近，大约在[0.05，0.15]上。
接下来用fzero命令求解，得到当F=25N时，c约为0.1155m，再通过cosα=b√22，求得cosα=0.9077,α=0.4329(rad)=24.81o[2].【讨论】1.推力F=0时，显然c=0，α=0，实际上是不满足以上所列方程的。
由常识，门的另一侧应有阻挡物以平衡气体压强带来的推力。
2.题图中画的F作用点位置与H的距离也没有改变，然而从常理推断应该是门的总长度不变才对，何况若距离不变，方程就是一个简单的幂函数方程，没有编程计算的必要了。
2.进一步研究F与α的关系，即，在题中所给条件下花费多少力气能使门打开一定的角度：整理上面用来求解的方程并将c与α关系代入得：F=5∗π∗1cosα∗（1−0.5∗tanα）1.4画出右端函数图像为可见当角度较小时，增大角度是比较容易的，而当α超过0.8rad左右时，每增大一点角度，所施加的力都要增加很多，这也是符合常识的。
【命令代码】[1]fplot('(1-2*c)^1.4-0.8*pi*sqrt(0.0625+c^2)',[0,0.5]);[2]function [ y ] = push( c )y=(1-2*c)^1.4-0.8*pi*sqrt(0.0625+c^2);end %非线性方程函数x=fzero(@push,[0.05,0.15]); %求解cosalfa=0.25/sqrt(0.0625+x^2); %根据c算cosαalfa=acos(ans) %得到αgrade=alfa*180/pi %α的角度表示法问题3用迭代公式x k+1=a∗x k∗exp (−b∗x k)计算序列{x k}分析其收敛性，其中a分别取5，11，15；b(>0)任意，初值x0=1。
观察是否有混沌现象出现，并找出前几个分岔点，观察分岔点的极限趋势是否符合Feigenbaum常数提示的规律.【分析与预测】非线性差分方程参数不同时，迭代序列会产生不同结果。
有可能收敛，有可能存在多个收敛的子列，甚至可能完全不收敛（混沌）。
本题中改变a的取值来研究这种现象。
先取几个分立的a值来观察它们的收敛情况，再在一定范围内改变一系列a值将它们的收敛序列以图像形式呈现出来，来观察收敛形式改变的节点位置和混沌现象。
先从理论上计算分岔点：收敛到唯一值的情况：对原迭代公式，令x k+1=x k=x，得到方程x=a∗x∗exp (−b∗x)，它的根x=ln a b就是这个非线性差分方程收敛到唯一值的平衡点，其稳定的条件是�f′(ln a b)�<1，代入得到1<a<e2（说明收敛情况与b无关）。
因此题中取a=5时，原序列只有一个极限。
收敛到两个值的情况：考查x k+2与x k的关系，并在等式中令x k+2=x k=x，可以列出方程如下x=a2xe−bx(1+ae−bx)这个方程我就不会解了，收敛到两个值以上的情况通过实验来得到结果。
【收敛性】如题意，分别令a=5，11，15，b=a，进行足够多次迭代，观察迭代序列的x值[1]：a=5时，数列很快地收敛到一个值，约为0.3219：a=11时，数列有两个收敛终点，分别为0.0747和0.3613（因为有两个终点，为了显示得更明显，增加了迭代次数）：a=15时，将迭代次数增加到120，仍然没有收敛到有限个值的迹象：【混沌现象】鉴于a=5时仍然可以收敛到唯一值，而a=15时已经出现了不能收敛的现象，将待观察a的取值范围定在[4,20]，出于计算量和方便观察两个角度，步长定为0.01，b暂定为5，参照书上的例子选择打印每个收敛序列中的第100到200个点[2]，得到的图形如下：明显出现了混沌现象。
第一个分岔点与计算相符，为e2≈7.4，第二个分岔点从图中读出约为12.5，第三个约为14.2，第四个约为14.6，后面就没办法看出2n的分岔现象了，下面来检验Feigenbaum 规律：b2−b1b3−b2=3b3−b2b4−b3=4.25后面的无法验证了，但是4.25与4.67已经比较接近了。
【讨论】书中提到后面会出现非2n收敛的窗口，而在观察混沌现象时，我确实发现在混沌区又出现了两条间断的白线，表明当a取这两种值时分别收敛到6个值和5个值，因为序数越大趋势越明显，为了看得更仔细，我打印了每个收敛序列里的第800到1000个点如下：确实有这种情况。
【命令代码】[1]a=5;n=40; %指定迭代次数x(1)=1;for k=1:(n-1)x(k+1)=a*x(k)*exp(-5*x(k)); %循环迭代end[2]function [ y ] = iter01( x,a,b )y=a*x*exp(-b*x);end %迭代公式函数function chaos( iter_fun,r1,r2,n )x0=1;kr=0;for rr=r1(1):r1(3):r1(2);kr=kr+1;y(kr,1)=feval(iter_fun,x0,rr,r2);for i=2:n(2);y(kr,i)=feval(iter_fun,y(kr,i-1),rr,r2);end %每个a值迭代一次，结果保存在y的一行里endplot([r1(1):r1(3):r1(2)],y(:,n(1)+1:n(2)),'k.');%对应每个a值，将它对应的迭代序列里第n（1）到n（2）的值都打印出来endchaos(@iter01,[4,20,0.01],11,[100,200]); %打印混沌图形。