求平方根问题 (C++ 实现)

下面是用二分法和牛顿迭代法求一个正数的平方根。

二分法

这里的题目稍微宽了一点点，包含了整数和小数的情况，这里二分法就不用多说了，如果中间值的平方与目标值在误差范围内，则返回，否则根据大小情况改变左/右区间的端点。

#include<iostream>
using namespace std;
int main()
{
	double eps = 1e-6;
	double k = 0.8;  # 输入正数
	double l = 0.0,r,mid;
	if (k>=1) r = k;  # 若输入正数大于1，则右端点设为 k
    if (k<1)  r = 1;  # 若输入整数小于1，则右端点设为 1

	while (fabs(l-k/l)>eps)
	{
		mid = l + (r - l) /2 ;   
		if (mid<k/mid)
		{
			l = mid;
		}
		else {
			r = mid;
		}
	}
	cout << l << endl;
	return 0;
}

上面主要是注意两点小细节，一个是 mid 的计算方式：

mid = (l+r)/2;
mid = l + (k-l)/2;  # 推荐

还有判断条件：

if(mid*mid<k){}
if(mid < k/mid){}  # 推荐

这里推荐的原因是推荐使用的这两种方式都可以防止计算的数值溢出的情况。

牛顿迭代法

这里主要是使用了方程的泰勒展开式，将求平方根问题看作求 $f(x)=x^2-k$ 与 $x$ 轴的交点。.

如上图所示，就是求出曲线上某一点（初始时为二分右端点）切线与x轴的交点，然后将此交点代入原方程中，继续求切线，求切线与x轴交点，直到满足条件，即使平方根的值。 由泰勒方程得切线 L 为：

.

$L:y=f(x_0)+f^{\prime }(x_0)(x-x_0)$

所以切线交点 x ₁ 为：
.

$x_1=x_0-\frac{f(x_0)}{f^{\prime }(x_0)}$

以此类推：
.

$x_{n+1}=x_n-\frac{f(x_n)}{f^{\prime }(x_n)}$

以上类推来看：
.

$\begin{array}{l}f\left(x\right)=x^m-k\\ f^{\prime }\left(x\right)=m{x}^{m-1}\\ x_{n+1}=x_n-\frac{x_n^m-k}{m{x}^{m-1}}=\left(1-\frac{1}{m}\right)x_n+\frac{k}{m{x}_n^{m-1}}\end{array}$

若是求平方根可以将 m=2 直接带入得到：
.

$x_{n+1}=\frac{1}{2}\left(x_n+\frac{k}{x_n}\right)$

由上可以得到：

#include<iostream>
using namespace std;

double mysqrt(double x)
{
	if (x == 1 || x == 0)
		return x;
	double temp = x / 2;
	while (1)
	{
		double a = temp;
		temp = (temp + x / temp) / 2;
		if (fabs(a - temp) < 0.000001)
		{
			return temp;
		}
	}
}

int main()
{
	double k = 0.5;
	double result= 0.0;
	
	lresult= mysqrt(k);

	cout << result << endl;
	return 0;
}

进一步，我们将内容拓展到 正数的 3次方根，4次方根…，根据上面的迭代式可以得到：

#include<iostream>
using namespace std;

double func(double x, double target,int n)  # F(x)
{
	return pow(x, n)-target;
}

double funct(double x, int n)  # F‘(x)
{
	return n*pow(x, n - 1);
}

double myRootEq(double target,int n)
{
	double x1, x0;
	x0 = target;
	if (x0 == 1 || x0 == 0)
		return x0;
	while (1)
	{
		x1 = x0 - func(x0, target,n) / funct(x0, n);  # 递推式
		if (fabs(x1 - x0)<0.000001)
		{
			return x1;
		}
		else
		{
			x0 = x1;
		}
	}
}

int main()
{
	double k = 0.5;
	double result= 0.0;
	
	result= myRootEq(k,3);

	cout << result << endl;
	return 0;
}

参考博客：https://blog.csdn.net/ycf74514/article/details/48996383