牛顿迭代法实现sqrt函数——学习笔记

牛顿迭代法又称为牛顿-拉弗森方法，它是一种在实数域和复数域上近似求解方程的方法，简而言之就是一种快速寻找近似零点的方法。牛顿法具体我不多介绍，可以直接看这篇文章，非常详细。

由于很多的高阶方程没法直接求解，就是利用一个公式来不断的近似逼近零点：

sqrt函数对应于方程组 的解，由于我们没法直接求得的值，所以借助于牛顿迭代法来计算。可设：

带入牛顿迭代法公式可得：

OK，公式推出来了，用代码来实现一下：

double mySqrt(double a) {
	double x = a;

	while (x*x > a) {
		x = (x + a / x) / 2;
	}

	return x;
}

 

实现了平方根以后，那么我们想开k阶根怎么办？同样利用牛顿迭代法就可以推出一个k阶根公式：

代码如下：

double Sqrt(double a,int k) {
	if (k == 1) {
		return a;
	}

	double x = a;

	while ( pow(x, k) > a) {
		x = (x*(k - 1) + a / pow(x, k - 1)) / k;
	}

	return x;
}

 

PS

这里提一个非常神奇的算法，为什么称之为神奇的算法呢，因为这个算法采用了一堆看不懂的计算方法，而且效率还比单纯用牛顿法更高效，而且至今也不清楚算法的创作者是谁。算法叫平方根倒数速算法。

float Q_rsqrt( float number )
{
    long i;
    float x2, y;
    const float threehalfs = 1.5F;
    x2 = number * 0.5F;
    y   = number;
    i   = * ( long * ) &y;   // evil floating point bit level hacking
    i   = 0x5f3759df - ( i >> 1 ); // what the fuck?
    y   = * ( float * ) &i;
    y   = y * ( threehalfs - ( x2 * y * y ) ); // 1st iteration
    // y   = y * ( threehalfs - ( x2 * y * y ) ); // 2nd iteration, this can be removed

    #ifndef Q3_VM
    #ifdef __linux__
         assert( !isnan(y) ); // bk010122 - FPE?
    #endif
    #endif
    return y;
}

这个求出的是平法根的倒数，“what the fack?”那里有一个非常神奇的数字，哈哈我是看不懂了，有兴趣的朋友去研究研究吧。