技术宅的结界

 找回密码
 立即注册→加入我们

QQ登录

只需一步,快速开始

搜索
热搜: 下载 VB C 实现 编写
查看: 1444|回复: 3
收起左侧

【C】不使用sqrt函数计算开平方

[复制链接]

85

主题

263

帖子

3427

积分

用户组: 管理员

No. 418

UID
418
精华
13
威望
52 点
宅币
1969 个
贡献
1026 次
宅之契约
0 份
在线时间
252 小时
注册时间
2014-8-9
发表于 2015-9-26 14:31:33 | 显示全部楼层 |阅读模式

欢迎访问技术宅的结界,请注册或者登录吧。

您需要 登录 才可以下载或查看,没有帐号?立即注册→加入我们

x
不使用sqrt函数计算开平方
零、介绍
@元始天尊
和天尊兄第一次见面,他就给我出了一道难题:
不使用Sqrt函数,只能使用+-*/和流程控制,实现sqrt函数。
一、开始
我们先来实验一个数:sqrt(2)
1.JPG

然后设初始值为1,将1的代入公式:
2.JPG
得到x=1.5,将1.5再次代入公式:
3.JPG
得到x = 1.416666666,再次带入
第四次,所得到的x值已经为:x=1.414213562和同等小数精度的sqrt值相等
至此我们已经得到sqrt(2)的迭代公式:
4.JPG
由此,我们来推广sqrt(A),A属于任意实数的迭代公式:
5.JPG
我们来试一下:
假设A=5,求sqrt(A)的值:
6.JPG
反复带入:
第2次sqrt(5) = 2.333333.....
第3次sqrt(5) = 2.238095238
第4次sqrt(5) = 2.236068896
第5次sqrt(5) = 2.236067978
第6次sqrt(5) = 2.236067977
第7次sqrt(5) = 2.236067977
.........
所以,对于任意实数A,其sqrt(A)的迭代公式为:
6.JPG
二、关于迭代精度问题
假设迭代初始值都从1开始,我们分别对2、3、4、5、6、7、8、9进行计算
当迭代结果精度与sqrt函数返回值小数第九位的精度相等,表示我们已经达到目的
开始测试:
对象2
Sqrt(2) ~=1.414213562
1.Sqrt(2) ~=1.5
2. Sqrt(2) ~=1.416666667
3. Sqrt(2) ~=1.414215686
4. Sqrt(2) ~=1.414213562(OK)
对象 3
Sqrt(3) ~=1.732050808
1. Sqrt(3) ~=2
2. Sqrt(3) ~=1.75
3. Sqrt(3) ~=1.732142857
4. Sqrt(3) ~=1.73205081
5. Sqrt(3) ~=1.732050808(OK)
对象 4
Sqrt(4) =2
1. Sqrt(4) ~=2.5
2. Sqrt(4) ~=2.05
3. Sqrt(4) ~=2.000609756
4. Sqrt(4) ~=2.000000093
5. Sqrt(4) ~=2.000000000(OK)
对象 5
Sqrt(5) ~=2.236067977
1. Sqrt(5) ~=3
2. Sqrt(5) ~=2.333333333
3. Sqrt(5) ~=2.238095238
4. Sqrt(5) ~=2.236068896
5. Sqrt(5) ~=2.236067978
6. Sqrt(5) ~=2.236067977(OK)
对象 6
Sqrt(6) ~=2.449489743
1. Sqrt(6) ~=3.5
2 Sqrt(6) ~=2.607142857
3. Sqrt(6) ~=2.454235636
4. Sqrt(6) ~=2.449494372
5. Sqrt(6) ~=2.449489743(OK)
对象 7
Sqrt(7) ~=2.645751311
1. Sqrt(7) ~=4
2. Sqrt(7) ~=2.875
3. Sqrt(7) ~=2.654891304
4.Sqrt(7) ~=2.645767044
5. Sqrt(7) ~=2.645751311(OK)
对象 8
Sqrt(8) ~=2.828427125
1. Sqrt(8) ~=4.5
2. Sqrt(8) ~=3.138888889
3. Sqrt(8) ~=2.843780728
4. Sqrt(8) ~=2.828468572
5. Sqrt(8) ~=2.828427125(OK)
对象 9
Sqrt(9) =3
1.Sqrt(9) ~=5
2. Sqrt(9) ~=3.4
3. Sqrt(9) ~=3.023529412
4. Sqrt(9) ~=3.000091554
5. Sqrt(9) ~=3.000000001
6. Sqrt(9) ~=3.000000000(OK)
结论,可见,要使sqrt函数精确到小数点后第九位,迭代初始数为1,迭代次数5次以上即可。
三、算法的优化
我们说,迭代初始值从1开始还是太慢,计算量太大了。因为,如果我们计算sqrt到小数点后1千位呢?所以,这一小节我们试着探索优化迭代初始值。
那么,经过猜测迭代初始值应该越接近sqrt函数的返回值(也就是x的开平方)越好:
考虑到算法的复杂性要求,我们将x迭代初始值设为(x/2)。因为:对于任一实数r有sqrt(r)恒大于等于r/2(几何原理)。
那么根据迭代公式:
7.JPG
改变迭代初始值为:
8.JPG
我们继续试验:
这里我用5,9,3来进行测试,因为上面测试结果中5和9的最大精确迭代次数都超过了平均值,而数字3则是从先前测试结果大于等于平均值的对象中随机选取的。
对象 5
迭代初始值=5/2
Sqrt(5) ~=2.236067977
1. Sqrt(5) ~=2.25
2. Sqrt(5) ~=2.236111111
3. Sqrt(5) ~=2.236067978
4. Sqrt(5) ~=2.236067977(OK)
对象 9
迭代初始值=9/2
Sqrt(9) =3
1.Sqrt(9) ~=4.5
2. Sqrt(9) ~=3.25
3. Sqrt(9) ~=3.009615385
4. Sqrt(9) ~=3.00001536
5. Sqrt(9) ~=3(OK)
对象 3
迭代初始值=1.5
Sqrt(3) ~=1.732050808
1. Sqrt(3) ~=1.75
2. Sqrt(3) ~=1.732142857
3. Sqrt(3) ~=1.732050801
4. Sqrt(3) ~=1.732050808(OK)
那么根据以上结果,三个测试数最小精确迭代次数均<=先前平均值。
从而证明:迭代初始值设为x/2能够起到一部分减小运算量的效果。
因为测试数据量不够,所以此结果仅仅作为猜想。
四、算法的C语言表示:
[C] 纯文本查看 复制代码
#include <stdio.h>
#include <math.h>

#define ITERATIONTIMES (128) //最大精确迭代次数

double CustomSquare(double number)
{
	int i;
	if (number == 0.0)
		return 0.0;
	double result = number / 2.0; //迭代初始数
	for (i = 1; i <= ITERATIONTIMES; i++)
	{
		result = 0.5 * (result + number / result);
	}
	return result;
}

int main(void)
{
	double di;
	printf("Input number:");
	scanf("%lf", &di);
	printf("Result by sqrt        function:%1.9f\n", sqrtf(di));
	printf("Result by custom sqrt function:%1.9f\n", CustomSquare(di));
}

不使用sqrt函数计算开平方.docx (20.05 KB, 下载次数: 10)

本帖被以下淘专辑推荐:

In the beginning I was not the best.
And the world was also not the best.
But I still know that I am who I am.
Because I think that it is good.
I have been working hard.
I have been keeping growth with the world.
And it was so.

2

主题

30

帖子

135

积分

用户组: 小·技术宅

UID
2285
精华
0
威望
0 点
宅币
105 个
贡献
0 次
宅之契约
0 份
在线时间
18 小时
注册时间
2017-2-25
发表于 2017-9-28 15:38:45 | 显示全部楼层
本帖最后由 乘简 于 2017-9-28 15:40 编辑

网上抄了个代码,不用循环,也可以计算sqrt,精度与库函数一样,号称比库函数更快的算法。。。
[C] 纯文本查看 复制代码
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  
  
    return y;  
}

19

主题

51

帖子

1160

积分

用户组: 上·技术宅

UID
1043
精华
5
威望
21 点
宅币
1027 个
贡献
15 次
宅之契约
0 份
在线时间
202 小时
注册时间
2015-8-15
发表于 2017-10-16 10:17:56 | 显示全部楼层
其实可以参考下牛顿迭代法
Only via chaos and excellence can beauty and success be respectively created and achieved.

0

主题

1

帖子

15

积分

用户组: 初·技术宅

UID
3846
精华
0
威望
2 点
宅币
10 个
贡献
0 次
宅之契约
0 份
在线时间
0 小时
注册时间
2018-5-15
发表于 2018-5-15 00:44:42 | 显示全部楼层
乘简 发表于 2017-9-28 15:38
网上抄了个代码,不用循环,也可以计算sqrt,精度与库函数一样,号称比库函数更快的算法。。。
[mw_shl_cod ...


这个代码中的 0x5f3759df 大有来头,有兴趣可以搜索一下,建议加个关键词 约翰卡马克

本版积分规则

QQ|申请友链|Archiver|手机版|小黑屋|技术宅的结界 ( 滇ICP备16008837号|网站地图

GMT+8, 2018-5-24 04:20 , Processed in 0.103603 second(s), 18 queries , Gzip On, Memcache On.

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

快速回复 返回顶部 返回列表