algoc++之求解根号2

前言

目的：学习开源项目algorithm-cpp，记录学习过程，有问题欢迎各位批评指正!

作者：shouxieai

项目地址：https://github.com/shouxieai/algorithm-cpp

本次任务：利用梯度下降法和牛顿法求解x\scriptsize \sqrt{x}x

PS：由于CSDN公式渲染效果太差，有点影响体验。

1.问题描述

使用梯度下降法和牛顿法实现求解根号x
x=?\sqrt{x}=? x=?

2.梯度下降法

问题思考方式：

第一步：转化问题，将x\scriptsize \sqrt{x}x转化为L(t)=(t2−x)2\scriptsize L(t)=(t^2-x)^2L(t)=(t2−x)2，当L(t)=0\scriptsize L(t)=0L(t)=0时，t\scriptsize tt就是计算得出的结果

第二步：寻找合适的解，找t\scriptsize tt使得L(t)=(t2−x)2=0\scriptsize L(t)=(t^2-x)^2 = 0L(t)=(t2−x)2=0

第三步：找到的t\scriptsize tt，就是得出的结果

转化问题的目的是，求解L(t)\scriptsize L(t)L(t)的极小值，也就是导数为0的点

补充知识点

L(t)\scriptsize L(t)L(t)为偶函数，关于y轴对称，只需要考虑[0,+∞)\scriptsize [0,+∞)[0,+∞)
在[0,+∞)\scriptsize [0,+∞)[0,+∞)上L(t)\scriptsize L(t)L(t)为凹函数，参考自什么是凹函数和凸函数?
如果凹函数在区间上的局部极小值存在则一定是该区间最小值。参考自函数凹凸性和极值
L(t)=0\scriptsize L(t)=0L(t)=0时取最小值同时也是极小值

如何找t\scriptsize tt，使得L(t)=0\scriptsize L(t)=0L(t)=0，注意L(t)\scriptsize L(t)L(t)取0是L\scriptsize LL函数的极小值

核心思想:

第一步：随机给定一个初始值t0\scriptsize t_0t0

第二步：在这个位置观察，t\scriptsize tt应该如何调整能够使得函数值变小，即L(t1)<L(t0)\scriptsize L(t_1)<L(t_0)L(t1)<L(t0)

第三步：不停根据第二步做调整，直到L(tn)\scriptsize L(t_n)L(tn)足够接近0为止(或调整量足够小)

如何调整使函数值变小

具体如何调整呢？在t的领域内对比L(t−Δt)\scriptsize L(t-\Delta t)L(t−Δt)与L(t+Δt)\scriptsize L(t+\Delta t)L(t+Δt)

若L(t−Δt)<L(t)\scriptsize L(t-\Delta t)<L(t)L(t−Δt)<L(t)，则t1=t−Δt\scriptsize t1=t-\Delta tt1=t−Δt

若L(t+Δt)<L(t)\scriptsize L(t+\Delta t)<L(t)L(t+Δt)<L(t)，则t1=t+Δt\scriptsize t1=t+\Delta tt1=t+Δt

由于考虑的是t\scriptsize tt的领域内，因此Δt\scriptsize \Delta tΔt是大于0的无穷小数

可以改写为

若L(t)−L(t−Δt)>0\scriptsize L(t)-L(t-\Delta t)>0L(t)−L(t−Δt)>0，则t1=t−Δt\scriptsize t_1=t-\Delta tt1=t−Δt

若L(t+Δt)−L(t)<0\scriptsize L(t+\Delta t)-L(t)<0L(t+Δt)−L(t)<0，则t1=t+Δt\scriptsize t_1=t+\Delta tt1=t+Δt

合成后如下：
t1=t−L(t+Δt)−L(t)∣L(t+Δt)−L(t)∣Δtt_{1}=t-\frac{L(t+\Delta t)-L(t)}{|L(t+\Delta t)-L(t)|} \Delta t t1=t−∣L(t+Δt)−L(t)∣L(t+Δt)−L(t)Δt

t领域内

对式子进行变换
t1=t−L(t+Δt)−L(t)∣L(t+Δt)−L(t)∣Δt(1)t_{1}=t-\frac{L(t+\Delta t)-L(t)}{|L(t+\Delta t)-L(t)|} \Delta t \tag1 t1=t−∣L(t+Δt)−L(t)∣L(t+Δt)−L(t)Δt(1)

t1=t−L((t+Δt)−L(t)Δt∣L(t+Δt)−L(t)Δt∣Δt=t−L(t+Δt)−L(t)ΔtΔt∣L(t+Δt)−L(t)Δt∣(2)t_{1}=t-\frac{\frac{L((t+\Delta t)-L(t)}{\Delta t}}{\left|\frac{L(t+\Delta t)-L(t)}{\Delta t}\right|} \Delta t=t-\frac{L(t+\Delta t)-L(t)}{\Delta t} \frac{\Delta t}{\left|\frac{L(t+\Delta t)-L(t)}{\Delta t}\right|} \tag2 t1=t−∣∣∣ΔtL(t+Δt)−L(t)∣∣∣ΔtL((t+Δt)−L(t)Δt=t−ΔtL(t+Δt)−L(t)∣∣∣ΔtL(t+Δt)−L(t)∣∣∣Δt(2)

t1=t−αL(t+Δt)−L(t)Δt(α→0+)(3)t_{1}=t-\alpha \frac{L(t+\Delta t)-L(t)}{\Delta t} \quad\left(\alpha \rightarrow 0^{+}\right) \tag3 t1=t−αΔtL(t+Δt)−L(t)(α→0+)(3)
得出迭代式，t1\scriptsize t_1t1只需要取t−αL′(t)\scriptsize t-\alpha L^{\prime}(t)t−αL′(t)，就能够保证L(t1)≤L(t0)\scriptsize L(t_1) \leq L(t_0)L(t1)≤L(t0)，并使得L\scriptsize LL逐渐下降，直到逼近0为止，注意这里的α\scriptsize \alphaα必须取趋近于0的无穷小时，该行为才能恒成立，α\scriptsize \alphaα也被视作步长

然而真实世界是

t1=t−αL(t+Δt)−L(t)Δt(α→0+)t_{1}=t-\alpha \frac{L(t+\Delta t)-L(t)}{\Delta t} \quad\left(\alpha \rightarrow 0^{+}\right) t1=t−αΔtL(t+Δt)−L(t)(α→0+)

1.当α\scriptsize \alphaα取无穷小时，虽然一定保证下降，但效率太慢
2.日常设计的很多函数，可以运行使用相对大一些的步长，比如α=0.01\scriptsize \alpha=0.01α=0.01。原因在于虽然步长大了可能跳过合适位置，使得L(t1)>L(t0)\scriptsize L(t_1)>L(t_0)L(t1)>L(t0)，但是在下一个时刻，依旧可能跳回来使得L(t2)<L(t1)\scriptsize L(t_2)<L(t_1)L(t2)<L(t1)
3.大的步长不能保证一定收敛，但是大部分时候是可以很好的工作
4.因此步长α\scriptsize \alphaα，我们常称之为学习率，通常会给一个相对小的数字，但不会太小

代码如下

float x = 2;
float t = x / 2;
float L = (t * t - x) * (t * t - x);
float a = 0.01;while(L > 1e-5){float delta = 2 * (t * t - x) * 2 * t;t = t - a * delta;L = (t * t - x) * (t * t - x);
}
cout << "result: " << t << endl;

第一步：初始化t=x2t=\frac{x}{2}t=2x，t\scriptsize tt也可以初始化为随机值，不过加入先验给一个更好的初始化，使得求解更迅速

第二步：计算L(t)=(t2−x)2\scriptsize L(t)=(t^2-x)^2L(t)=(t2−x)2

第三步：若L(t)>1e−5\scriptsize L(t)>1e-5L(t)>1e−5，则继续迭代，否则停止迭代(注：delta足够小时也可以停止迭代)

第四步：计算导数dL(t)dt=2(t2−x)⋅2t\scriptsize \frac{d L(t)}{d t}=2\left(t^{2}-x\right) \cdot 2 tdtdL(t)=2(t2−x)⋅2t

第五步：更新t+=t−α⋅dL(t)dt\scriptsize t^{+}=t-\alpha \cdot \frac{d L(t)}{d t}t+=t−α⋅dtdL(t)，这里α\scriptsize \alphaα取0.01

第六步：继续执行第二步

从泰勒展开来理解梯度下降法

对于f(x)\scriptsize f(x)f(x)的一阶泰勒展开，表示为：
f(x)=f(x0)+f′(x0)(x−x0)f(x)=f\left(x_{0}\right)+f^{\prime}\left(x_{0}\right)\left(x-x_{0}\right) f(x)=f(x0)+f′(x0)(x−x0)
对应的L\scriptsize LL在t0\scriptsize t_0t0处的一阶泰勒展开近似表示为(t1\scriptsize t_1t1与t0\scriptsize t_0t0都是领域)：
L(t1)=L(t0)+L′(t0)(t1−t0)L\left(t_{1}\right)=L\left(t_{0}\right)+L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right) L(t1)=L(t0)+L′(t0)(t1−t0)
为了使得L(t1)<L(t0)\scriptsize L(t_1)<L(t_0)L(t1)<L(t0)
L(t0)+L′(t0)(t1−t0)<L(t0)L′(t0)(t1−t0)<0L\left(t_{0}\right)+L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right)<L\left(t_{0}\right) \\ L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right)<0 L(t0)+L′(t0)(t1−t0)<L(t0)L′(t0)(t1−t0)<0

由于t0\scriptsize t_0t0与t1\scriptsize t_1t1是领域，则∣t1−t0∣=Δt\scriptsize \left|t_{1}-t_{0}\right|=\Delta t∣t1−t0∣=Δt，又L′(t0)(t1−t0)<0\scriptsize L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right)<0L′(t0)(t1−t0)<0，则必有L′(t0)\scriptsize L^{\prime}\left(t_{0}\right)L′(t0)与t1−t0\scriptsize t_1-t_0t1−t0符号相反

令
t1−t0=−L′(t0)∣L′(t0)∣Δtt_{1}-t_{0}=-\frac{L^{\prime}\left(t_{0}\right)}{\left|L^{\prime}\left(t_{0}\right)\right|} \Delta t t1−t0=−∣L′(t0)∣L′(t0)Δt
则有
t1=t0−L′(t0)∣L′(t0)∣Δtt1=t0−αL′(t0)(α→0+)t_{1}=t_{0}-\frac{L^{\prime}\left(t_{0}\right)}{\left|L^{\prime}\left(t_{0}\right)\right|} \Delta t \\ t_{1}=t_{0}-\alpha L^{\prime}\left(t_{0}\right) \quad\left(\alpha \rightarrow 0^{+}\right) t1=t0−∣L′(t0)∣L′(t0)Δtt1=t0−αL′(t0)(α→0+)
总结：

1.梯度下降法是通过观察局部，决定如何调整的算法。如果函数具有多个极值，则可能陷入局部极值，此时初始点的选择直接影响收敛结果
2.大的步长在一定程度上可能跨过局部极值，但也可能造成震荡导致不收敛
3.步长的选择，需要根据函数的特性来找到合适取值，若导数特别大时，则步长取小，导数小时，步长取大。否则很容易造成收敛问题
4.存在一类算法，可以在一定范围内搜索一个合适步长，使得每一次迭代更加稳定

3.牛顿法

问题思考方式：

第一步：转化问题，将x\scriptsize \sqrt{x}x转化为L(t)=(t2−x)\scriptsize L(t)=(t^2-x)L(t)=(t2−x)，当L(t)=0\scriptsize L(t)=0L(t)=0时，t\scriptsize tt就是计算得出的结果

第二步：寻找合适的解，找t\scriptsize tt使得L(t)=(t2−x)=0\scriptsize L(t)=(t^2-x) = 0L(t)=(t2−x)=0

第三步：找到的t\scriptsize tt，就是得出的结果

具体实现：

考虑L(t)\scriptsize L(t)L(t)在t0\scriptsize t_0t0处的切线与x轴交点作为t1\scriptsize t_1t1，不断逼近零点
如下图所示，考虑以o2(t0,L(t0))\scriptsize o2(t_0,L(t_0))o2(t0,L(t0))为原点，则切线可表示为k=L′(t0),b=0\scriptsize k=L^{\prime}\left(t_{0}\right),b=0k=L′(t0),b=0，而与x轴交点可表示为：

k=d(o2,t0)d(t1,t0)=L(t0)t0−t1=L′(t0)k=\frac{d\left(o 2, t_{0}\right)}{d\left(t_{1}, t_{0}\right)}=\frac{L\left(t_{0}\right)}{t_{0}-t_{1}}=L^{\prime}\left(t_{0}\right) k=d(t1,t0)d(o2,t0)=t0−t1L(t0)=L′(t0)

t0−t1=L(t0)L′(t0)t1=t0−L(t0)L′(t0)t_{0}-t_{1}=\frac{L\left(t_{0}\right)}{L^{\prime}\left(t_{0}\right)} \\ t_{1}=t_{0}-\frac{L\left(t_{0}\right)}{L^{\prime}\left(t_{0}\right)} t0−t1=L′(t0)L(t0)t1=t0−L′(t0)L(t0)

代码如下

float x = 2;
float t = x/2;
float L = t * t - x;while(abs(L) > 1e-5){float dL = 2 * t;t = t - L / dL;L = t * t - x;
}
cout << "simple_nt1: " << t << endl;

第一步：随机给定一个初始值t=x2t=\frac{x}{2}t=2x

第二步：计算L(t)\scriptsize L(t)L(t)和L′(t)\scriptsize L^{\prime}(t)L′(t)

第三步：更新t+=t−L(t)L′(t)\scriptsize t^{+}=t-\frac{L(t)}{L^{\prime}(t)}t+=t−L′(t)L(t)

第四步：不停根据规则2做调整，直到L(t)\scriptsize L(t)L(t)足够接近0为止(或调整量足够小)

从泰勒展开来理解牛顿法

对于f(x)\scriptsize f(x)f(x)的一阶泰勒展开，表示为：
f(x)=f(x0)+f′(x0)(x−x0)f(x)=f\left(x_{0}\right)+f^{\prime}\left(x_{0}\right)\left(x-x_{0}\right) f(x)=f(x0)+f′(x0)(x−x0)
对应的L\scriptsize LL在t0\scriptsize t_0t0处的一阶泰勒展开近似表示为(t1\scriptsize t_1t1与t0\scriptsize t_0t0都是领域)：
L(t1)=L(t0)+L′(t0)(t1−t0)L\left(t_{1}\right)=L\left(t_{0}\right)+L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right) L(t1)=L(t0)+L′(t0)(t1−t0)
令L(t1)=0\scriptsize L(t_1)=0L(t1)=0
L(t0)+L′(t0)(t1−t0)=0t1−t0=−L(t0)L′(t0)t1=t0−L(t0)L′(t0)L\left(t_{0}\right)+L^{\prime}\left(t_{0}\right)\left(t_{1}-t_{0}\right)=0 \\ t_{1}-t_{0}=-\frac{L\left(t_{0}\right)}{L^{\prime}\left(t_{0}\right)} \\ t_{1}=t_{0}-\frac{L\left(t_{0}\right)}{L^{\prime}\left(t_{0}\right)} L(t0)+L′(t0)(t1−t0)=0t1−t0=−L′(t0)L(t0)t1=t0−L′(t0)L(t0)
若函数二阶可导，则可考虑导函数零点时方程的根

根据t+=t−L(t)L′(t)\small t^{+}=t-\frac{L(t)}{L^{\prime}(t)}t+=t−L′(t)L(t)，可令L=f′(t)\small L=f^{\prime}(t)L=f′(t)，得出t+=t−f′(t)f′′(t)\small t^{+}=t-\frac{f^{\prime}(t)}{f^{\prime\prime}(t)}t+=t−f′′(t)f′(t)。则函数为L(t)=(t2−x)2\small L(t)=(t^2-x)^2L(t)=(t2−x)2时

代码如下

float x = 2;
float t = x / 2;
float L = (t * t - x) * (t * t - x);while(L > 1e-5){float d1L = 2 * (t * t - x) * 2 * t;float d2L = 4 * t * 2 * t + 4 * (t * t - x);t = t - d1L / d2L;L = (t * t - x) * (t * t - x);
}
cout << "simple_nt2: " << t << endl;

第一步：随机给定一个初始值t=x2\small t=\frac{x}{2}t=2x

第二步：计算L(t)\small L(t)L(t)和L′(t)\small L^{\prime}(t)L′(t)

第三步：更新t+=t−L(t)L′(t)\small t^{+}=t-\frac{L(t)}{L^{\prime}(t)}t+=t−L′(t)L(t)

第四步：不停根据规则2做调整，直到L(t)\scriptsize L(t)L(t)足够接近0为止(或调整量足够小)

结语

学习开源项目algorithm-cpp，仅供自己参考。

参考

https://github.com/shouxieai/algorithm-cpp
什么是凹函数和凸函数?
函数凹凸性和极值