『笔记』求解方程的根

『笔记』求解方程的根

面试时一个经常被问到的小编程题就是求方程的根。可以参考：HJ107 求解立方根.

计算一个浮点数的立方根，不使用库函数。
保留一位小数。

解决这个问题主要有两种方法：二分法和使用高斯牛顿

二分法：不必多说，可以设置自变量的上界和下界，以二分的形式不断求mid，根据mid的f(mid)值判断下一步进入的区间，直到收敛至满足threshold.

// half search
int main() {
    double val;
    cin >> val;
    double thr = 1e-4;
    double left = -3;
    double right = 3;
    while (true) {
        double mid = (left + right) / 2;
        if (abs(mid * mid * mid - val) <= thr) {
            cout << fixed << setprecision(1) << mid << endl;
            return 0;
        }
        if (mid * mid * mid > val) {
            right = mid;
        }
        else {
            left = mid;
        }
    }

高斯牛顿法：再次回顾，我们知道牛顿法即，为了找到\(\min_xf(x)\)，\(f(x)\)是non-linear and non-convex且without constraints. 使用泰勒二次展开对局部函数进行估计，然后发现可以转换成unconstrained的quadratic问题，即可通过-2a分之b的简单方式得到\(x_{k+1}=x_k-\hat{H}_k^{-1}\nabla f(x_k)\)，其就像一个根据已知具体函数的信息智能化根据二阶导（Hessian）调制了步长的steepest descent. 如果\(f(x)\)是一个least square problem，即，经常以一个cost function的平方出现的话，类似\(f(x)=e^2(x)\)，求Hessian时，取到它求它不求加它不求它求时，第二项能够由于\(e(x)\)值较小而省略，Hessian成为\(\hat{H} (x) = 2\nabla e(x_k) e(x_k)\)，最终得到很好用的形式\(x_{k+1}=x_k-(\nabla e(x_k) \nabla^Te(x_k))^{-1} \nabla e(x_k) e(x_k)\)，我们发现不需要显式地求Hessian了。

那么对于我们的题的情况，首先1. 求根即让\(f(x)=x^3-t=0\)，于是这个立方值和target的差距就像一个cost function一样，于是构造的问题就是\(\min_x (x^3-t)^2\)了，\(e(x)=x^3-t\)，可以直接套了。2. 这些只是scalar而不是vector了，那么形式就更简单了：\(x_{k+1}=x_k- \frac{e(x_k)}{\nabla e(x_k)}\)

个人小tip：牛顿法是f(x)的gradient除以Hessian，高斯牛顿是e(x)的值自己除以e(x)的gradient. 为了简单地记忆（无数学意义），f(x)变成了特定问题下的e(x)，以及两个值都升了一个“档次”

#include <iostream>
#include <iomanip>

using namespace std;

// Gaussian-Newton algorithm
int main() {
    double val;
    cin >> val;
    double thr = 1e-4;
    double x = 3;
    while (abs(x * x * x - val) > thr) {
        x = x - (x * x * x - val) / (3 * x * x);
    }
    cout << fixed << setprecision(1) << x << endl;
}