1. 概述
介绍割线法求解含单变量非线性方程的过程。
2. 原理步骤
2.1 原理
求解非线性方程
f
(
x
)
=
0
f\left(x\right) =0
f(x)=0时,将函数在点处作Taylor展开:
f
(
x
)
=
f
(
x
k
)
+
f
′
(
x
k
)
(
x
−
x
k
)
+
f
′
′
(
ξ
)
2
!
(
x
−
x
k
)
2
\begin{equation} f(x)=f\left(x_{k}\right)+f^{\prime}\left(x_{k}\right)\left(x-x_{k}\right)+\frac{f^{\prime \prime}(\xi)}{2 !}\left(x-x_{k}\right)^{2} \end{equation}
f(x)=f(xk)+f′(xk)(x−xk)+2!f′′(ξ)(x−xk)2
取前两项近似代替
f
(
x
)
f\left(x\right)
f(x),则非线性方程
f
(
x
)
=
0
f\left(x\right)=0
f(x)=0用线性方程式(2)近似。
f
(
x
k
)
+
f
′
(
x
k
)
(
x
−
x
k
)
=
0
\begin{equation} f\left(x_{k}\right)+f^{\prime}\left(x_{k}\right)\left(x-x_{k}\right)=0 \end{equation}
f(xk)+f′(xk)(x−xk)=0
式(2)变换得到:
x
=
x
k
−
f
(
x
k
)
f
′
(
x
k
)
x=x_{k}-\frac{f\left(x_{k}\right)}{f^{\prime}\left(x_{k}\right)}
x=xk−f′(xk)f(xk)
得到迭代公式:
x
k
+
1
=
x
k
−
f
(
x
k
)
f
′
(
x
k
)
\begin{equation} x_{k+1}=x_{k}-\frac{f\left(x_{k}\right)}{f^{\prime}\left(x_{k}\right)} \end{equation}
xk+1=xk−f′(xk)f(xk)
公式(3)称Newton迭代公式,要求
f
′
(
x
k
)
≠
0
f^{\prime}\left(x_{k}\right) \neq 0
f′(xk)=0,实际中
f
(
x
k
)
f\left(x_{k}\right)
f(xk)一般比较复杂或者不存在解析函数,无法求出
f
′
(
x
k
)
f^{\prime}\left(x_{k}\right)
f′(xk)的表达式,改用
f
(
x
k
)
−
f
(
x
k
−
1
)
x
k
−
x
k
−
1
\frac{f\left(x_{k}\right)-f\left(x_{k-1}\right)}{x_{k}-x_{k-1}}
xk−xk−1f(xk)−f(xk−1)代替,得到迭代公式:
x
k
+
1
=
x
k
−
f
(
x
k
)
f
(
x
k
)
−
f
(
x
k
−
1
)
(
x
k
−
x
k
−
1
)
,
k
≥
2
\begin{equation} x_{k+1}=x_{k}-\frac{f\left(x_{k}\right)}{f\left(x_{k}\right)-f\left(x_{k-1}\right)}\left(x_{k}-x_{k-1}\right), k \geq 2 \end{equation}
xk+1=xk−f(xk)−f(xk−1)f(xk)(xk−xk−1),k≥2
公式(4)称割线法,在计算前需要提供两个初值
x
0
x_{0}
x0和
x
1
x_{1}
x1。
2.2 几何意义
回到问题的起点:如何找到满足
f
(
x
)
=
0
f\left(x\right)=0
f(x)=0的
x
x
x?
一般人会这么做:第一步随便试了个数,看看返回多少
x
1
→
f
(
x
1
)
≠
0
x_{1}\rightarrow f\left(x_{1}\right)\neq0
x1→f(x1)=0
再试一把
x
2
→
f
(
x
2
)
≠
0
x_{2}\rightarrow f\left(x_{2}\right)\neq0
x2→f(x2)=0
当然这并不是白费力气,你会发现当
x
x
x变化量为1时
f
(
x
)
f\left(x\right)
f(x)变化了多少,体育老师也会说的斜率
f
(
x
2
)
−
f
(
x
1
)
x
2
−
x
1
\frac{f\left(x_{2}\right)-f\left(x_{1}\right)}{x_{2}-x_{1}}
x2−x1f(x2)−f(x1)
倒过来变成
x
2
−
x
1
f
(
x
2
)
−
f
(
x
1
)
\begin{equation} \frac{x_{2}-x_{1}}{f\left(x_{2}\right)-f\left(x_{1}\right)} \end{equation}
f(x2)−f(x1)x2−x1
式(5)表示
f
(
x
)
f\left(x\right)
f(x)变化量为1时
x
x
x的变化量,那么
f
(
x
)
f\left(x\right)
f(x)从
f
(
x
2
)
→
0
f\left(x_{2}\right)\rightarrow0
f(x2)→0,
x
x
x变化量不就是
x
2
−
x
1
f
(
x
2
)
−
f
(
x
1
)
⋅
(
0
−
f
(
x
2
)
)
\frac{x_{2}-x_{1}}{f\left(x_{2}\right)-f\left(x_{1}\right)}\cdot\left (0-f\left(x_{2}\right)\right)
f(x2)−f(x1)x2−x1⋅(0−f(x2))
所以得到了第三个点:
x
3
=
x
2
+
x
2
−
x
1
f
(
x
2
)
−
f
(
x
1
)
⋅
(
0
−
f
(
x
2
)
)
=
x
2
−
f
(
x
2
)
f
(
x
2
)
−
f
(
x
1
)
(
x
2
−
x
1
)
\begin{equation} x_{3}=x_{2}+\frac{x_{2}-x_{1}}{f\left(x_{2}\right)-f\left(x_{1}\right)}\cdot\left(0-f\left(x_{2}\right)\right)=x_{2}-\frac{f\left(x_{2}\right)}{f\left(x_{2}\right)-f\left(x_{1}\right)}\left(x_{2}-x_{1}\right) \end{equation}
x3=x2+f(x2)−f(x1)x2−x1⋅(0−f(x2))=x2−f(x2)−f(x1)f(x2)(x2−x1)
这就得到了割线法的递归公式,用数学老师的话讲,割线法就是用线性方程
y
=
f
(
x
k
)
+
f
(
x
k
)
−
f
(
x
k
−
1
)
x
k
−
x
k
−
1
(
x
−
x
k
)
y=f\left(x_{k}\right)+\frac{f\left(x_{k}\right)-f\left(x_{k-1}\right)}{x_{k}-x_{k-1}}\left(x-x_{k}\right)
y=f(xk)+xk−xk−1f(xk)−f(xk−1)(x−xk)
的零点逐步近似曲线方程
f
(
x
)
=
0
f\left(x\right)=0
f(x)=0的零点。
2.3 步骤
步骤如下:
3. 结论
介绍了割线法的原理和步骤。
4. 详细代码
测试:函数
f
(
x
)
=
(
x
−
2
)
2
−
4
f(x)=(x-2)^{2}-4
f(x)=(x−2)2−4
二分法根据隔根区间搜索,步数多但结果较可控;割线法根据初值不同会取到意想不到的结果,比如x0为0.1,dx0分别为1和4时结果为0和4。
double MyFunction(double x) {
return (x - 2)*(x - 2) - 4;
}
int main()
{
cout << "SecantSearch:" << endl;
CSecantSearch SecantSearch(100);
try {
SecantSearch.Search(0.1, 4, 20, 0.01, MyFunction);
auto secantRecord = SecantSearch.getRecord();
cout << "i=" << secantRecord.index << endl;
cout << "x=" << secantRecord.data[secantRecord.index].x << endl;
cout << "f=" << secantRecord.data[secantRecord.index].fx << endl;
}
catch (invalid_argument &e) {
cout << e.what() << endl;
//......
}
catch (out_of_range &e) {
cout << e.what() << endl;
//......
}
return 0;
}
CBinarySearch.h
/*
* Auther: sanfan66
*/
#pragma once
class CSecantSearch
{
public:
struct SData {
double x, fx, deltaX;
};
struct SRecord {
SData *data;
int maxNum;
int index;
};
CSecantSearch(int mMaxNum) {
Record.maxNum = mMaxNum;
Record.data = new SData[Record.maxNum]{};
}
~CSecantSearch() {
delete[] Record.data;
Record.data = nullptr;
}
double Search(const double x0, const double deltaX0, const double maxDeltaX, const double ObjFunEps, double(*ObjFunction)(double));
double RangLimit(const double minValue, const double maxValue, double *value);
SRecord getRecord() {
return Record;
}
private:
SRecord Record;
};
CBinarySearch.cpp
#include "CSecantSearch.h"
#include <iostream>//cout
using namespace std;//cout
double CSecantSearch::Search(const double x0, const double deltaX0, const double maxDeltaX, const double ObjFunEps, double(*ObjFunction)(double)) throw(exception)
{
Record.data[0] = { x0, (*ObjFunction)(x0), deltaX0 };
if (abs(Record.data[0].fx) < ObjFunEps) {
Record.index = 0;
return 0;
}
Record.index = -1;//判断是否有解
for (int i = 1; i < Record.maxNum; i++) {
Record.data[i].x = Record.data[i - 1].x + Record.data[i - 1].deltaX;
Record.data[i].fx = (*ObjFunction)(Record.data[i].x);
Record.data[i].deltaX = -Record.data[i].fx / (Record.data[i].fx - Record.data[i - 1].fx)* (Record.data[i].x - Record.data[i - 1].x);
RangLimit(-maxDeltaX, maxDeltaX, &(Record.data[i].deltaX));
if (abs(Record.data[i].fx) < ObjFunEps) {
Record.index = i;
break;
}
}
if (Record.index < 0) {
cout << "Record.maxNum=" << Record.maxNum << endl;
throw out_of_range(":割线法迭代次数超过设定值");
}
return 0;
}
double CSecantSearch::RangLimit(const double minValue, const double maxValue, double *value)
{
if (*value > maxValue) {
*value = maxValue;
}
else if (*value < minValue) {
*value = minValue;
}
return 0;
}