数学研发论坛

 找回密码
 欢迎注册
查看: 3019|回复: 23

[擂台] 高精度除法之估商法

[复制链接]
发表于 2016-3-21 14:44:21 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有帐号?欢迎注册

x
我用VB完整的实现了估商法,支持正负,小数,科学计数法,之所以说完整实现,是因为这个程序是我正在写的(高精度表达式计算机)的子程序,能够和其它计算程序配合,实现高精度加减乘除,阶乘,乘方,开N次方,函数等混合运算.
例:能够运行万位精度,开10000次方根正确运算,开N次方根是用的迭代法,为什么举这个例子,因为参入运算的子程序如果有问题,是计算不出正确结果的,现在这个计算机还没完善:
先说优点,国内目前找到的表达式计算机(不支持高精度),而且不完善,几乎无纠错功能,最关键是不支持函数,函数是一切功能运算的基础,这些问题在我的程序里已经建立了解决框架.可以迅速加入纠错功能,各种函数.
再说缺点,开N次方没找到好的算法,我的算法是通过利用手算公式取初值,再迭代,     三角函数是用的泰勒公式,速度也是超慢,阶乘针对硬乘优化,10000!3秒左右,准备调用别人的大数库,
乘法是硬乘和FFT乘法分工,具体如何分工没有好的方法,目前采用的是1100位*1100位以内的数用硬乘,和任一单数400位之内用硬乘,
乘方采用网上找到的快速算法,不知道名字,就是A的8次方转换成A的平方的平方,然后再平方,速度非常快.很满意.
说了这么多主要是希望坛主和其他大大们能给出好的算法和建议,
另外我再谈谈我对除法的一些看法,不能利用计算硬件的算法不是好算法,所以设计算法时要充分考虑硬件的效能,简单说一下我知道的算法,第一,用减法实现,从高位向低位减,我最早用的是这个方法并实现,    第二个方法是用的估商(不是试商,试商非常慢),估出一个不大于正确商的数,做乘法,再按减法算法继续,一直觉得这个算法不错,在位数少时,速度还满意,
第三个方法是一个论文上说的方法,用的是十进制,,提前算出除数和1-9相乘的得数,其它步骤和方法二差不多,低位时速度比第二个方法稍慢,高位时甩开了我,另外这个论文中说到了迭代法,说很慢,
第四个方法是迭代法,我实现后就没管它了,和上面说的一致,超慢!
第五个方法就是我现在用的方法,为什么我一直没实现第三种方法,就是因为它采用的是十进制,不能充分利用CPU的乘法功能,需要说明的一点是现在乘法指令和加法一样快,你不相信?,自已写个代码去试试,所以不要再想一些用加减代乘法,或除法的算法,整数除法指令也只比乘法慢一点,所以为了提高性能,32位机上要用万进制乘,用万进制的难点就是商更加难估了,解决要点就是通过找一个参数使除数接近10000,当然通用除法的调试是相当麻烦的,所以上网上你只能看到人讲原理,却看不到实现,其实我们论坛看到一个人用C写估商法,因为他的程序运行困难,我没有调试,不过他在采用万进制的情况下,第一位除数采用五位数估商,应该有问题,我曾经有过类似的经历,不过没有成功,也可能他是对的,必竟我没调试过他的程序.
希望坛友们给予指导,建议!


补充内容 (2016-3-24 14:36):
下面有我回别人的贴子,在里面有对估商法较详细的看法,有不同意见,看法发出来一起计论!

评分

参与人数 1金币 +20 收起 理由
gxqcn + 20 首贴奖励,欢迎常来。

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2016-3-21 15:57:29 | 显示全部楼层
楼主,加油
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
回复

使用道具 举报

发表于 2016-3-21 18:17:31 | 显示全部楼层
嗯。我的大数库实现除法,主要就用浮点数除法(牛顿迭代法)和估商两种算法。在位数比较小的时候(几千位之内),估商比较快,再往上就是牛顿迭代比较快了,你觉得慢肯定是乘法没实现好(实现好后是O(n ln n )的),或者控制精度没做到位。
说估商。我现在用的是2^64进制,计算机硬件支持128bit整数除以64bit整数,只要商不大于64位即可。于是算P/Q的时候,我把Q的最高64位(记作 H )拿出来,加上1,除2^127。得到 R=floor(2^127/(H+1)),满足 2^63<=R<2^64。用这个R去乘P的高位,就能得到商的高位。具体的位对齐需要仔细算。实现后我做了实验,这样估商,一次能估出64个二进制位(19.3个十进制位),估的值总是小于实际值,且差不超过5.

乘方那个算法名称是“快速幂”。

三角函数/幂函数/对数函数,用泰勒公式是没问题的,有加速方法。binary spliting是一种通用方法,还有各种公式也能用来加速,比如sin(a+b)=sinacosb+sinbcosa。这方面中文的资料比较少。

给个建议,能用C/C++就尽量用……如果考虑到效率的话。而且指针是非常强大的工具。换语言很耗时间,所以只是建议。

评分

参与人数 1威望 +2 金币 +2 贡献 +2 鲜花 +2 收起 理由
落叶 + 2 + 2 + 2 + 2 很给力!

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2016-3-22 10:51:47 | 显示全部楼层
我正在设计一个全新的高精度对数算法,期望其性能达到一个新的高度,甚至超过目前这个星球上最快的软件或者库。请静候佳音,到时我会公布我的结果的。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-3-23 13:39:54 | 显示全部楼层
本帖最后由 落叶 于 2016-3-23 21:34 编辑

前天用手机回贴,乱码
今天才回,很感谢三楼给予的详细解答,我的FFT乘法也是前几天才完成,一直想自已写FFT转换,一直没搞懂细节,没办法,正好网上有C源码,修改了一下做成DLL,前几天才搞好,
因为原先的测试偏见,没想起从新测试一下迭代除法。
仔细看了一下你的估商实现,你把除数直接提出来的方法,相当于左移操作,算完结果后你应该还位(如改小数点位置,改指数,右移等等),我还是用万进制举个例子比较直观:
例:123456换成万进制后的形式是0012  3456,乘100后变成1234 5600,而你估商时是直接用1234,相当于除数乘以100,做完结果后得数要乘100,当然你也可以用其它方法,
对于你用除数加一估商,我在方法二时,也是这样用的,不过我用的是十进制,被除法的前两位除以除数的前一位加一,但我当时的想法是加在除数的最高位,而你是加在最低位,
虽然十进制中,我所取数的最高位就是最低位。
Knuth(高德纳)对估商做出了一个数学推导,公式一大堆,一句话,就是左移除数,乘2乘2...,例:0011 0101 左移成1101 0100,当时硬件差,他以左移替代乘法,估商最大误差2,
平均误差1,这种算法最佳实现是一次单精对多精乘法,一次多精减法,一次多精加法修正(我是直接加余数),书上说的竟然是用积减除数修正,哎,不说了,代沟啊!
  其实不管是你的方法,还是老K的,还有我的估商,总的思想就是整理除数,使除数尽量接近你所采用的进制的最大值,2^32进制时,就是65535,万进制就是10000,是不是不好理解,
举例:万进制时:下面的例子中因被除数第一位小于除数第一位,合并被除数第一位和第二位,整除除数第一位估商,如下:
0999 5678 1000/1000 9999 ,估商9995,正确商为9985,把被除数,除数乘9,8996 1102 9000/9008 9991,估商9986,另外不要问我被除数变了有没有影响,几乎没影响,经测试除数整理
好后,估商误差0.3左右,不超过0.5,我也是不明白,老K的主思路和我一样,为什么误差要大些,最后发现原因就在于左移,只有乘2,乘4,乘8,除数没有调整到最佳状态,为什么这么说呢?
如果一个除数需要调整,乘以3最合适,按老K的算法你就只能乘2了,或乘4了。我觉得老K的左移操作现在对高精度数没多大用处,高精度数本身是有数据头的,向这种左移右移的操作,
在数据头里把指数改一下就可以了,如我的万进制:1234*0100,只需要把123的指数位加2就行了
再说我的算法:
0999 5678 1000/0010 0099,除数乘100(不是真乘,只是调整数的位长度,和小数点位置),jw=2(这个用于后面调整小数位),变成0999 5678 1000/1000 9900,用(1000 9900)*1-9
反复试,找出一个最大的数,使它俩的积不超过100000000,这里是9,把被除数,除数乘9,好吧,进入主程序。
我的所有代码除了总控程序,其它都准备公开源码,没有很高的效率,就是比较完善,调试通过,支持正负,小数,科学计数。
但是现在没注释好,最关键是不好运行,因为没总控程序,只有等把总控做成DLL再说吧,当然如果你们不在意,现在就可贴出来。

关于你的程序我有一点小的看法:你用的进制是满进制,也就是2^64*2^64后,寄存器已经满了,后面的进位你不能一起处理,单独处理很耗时,你可能觉得满进制乘2,4,8,16不耗时,其实高精度数本身
是有数据头的,向这种左移右移的操作,在数据头里把指数改一下就可以了,如我的万进制:1234*0100,只需要把1234的指数位加2就行了。另外你的除数加1 ,你试过没加和加过后的区别吗?

关于三楼坛友对我的建议,我很感谢,其时我最早看过c,c++,mfc,汇编,还有很多很杂的计算机方面的书,就是没怎么看VB,当然因为自已的学习方法错误,一直是动眼不动手,所以都只是混了个脸熟,
知道个大概意思,已经N多年没有看计算机方面的书了,这次是下决心通过写个程序整合以往的知识,并增补新知识,因为英语0,又因VB的中文参考较好,心里其实已是准备两个语言都看看,我只是觉得语言
就是个工具,大同小异,大部分都相通。MFC虽不同,真心觉得不好用,但是你要想看类的应用,简直太棒了,请看它吧?说个笑话:你是喜欢开5档的车,还是喜欢开一百档的车,编程不是编程语言,只是
它们联系很紧而已。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-3-23 13:55:54 | 显示全部楼层
再来一贴,我刚才忘了问,为什么国内几乎看不到像样的表达式计算机,高精度更是没有,难道那些工程师,科学家,计算数据是算一项记一项,或M存等,真心觉得麻烦,
有熟知这方面的坛友回一下。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2016-3-23 16:04:32 | 显示全部楼层
落叶 发表于 2016-3-23 13:55
再来一贴,我刚才忘了问,为什么国内几乎看不到像样的表达式计算机,高精度更是没有,难道那些工程师,科学 ...


我觉得到处都是车轮子,您怎会找不到软件呢, 能否具体说一下你的需求?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-3-23 16:21:48 | 显示全部楼层
本帖最后由 落叶 于 2016-3-23 17:10 编辑

我说的是这种类型的表达式计算:
-(-1)(-2)-(-3)(-4)(-5) (-6)(-7)
(-8)(-9)cuberoot(27)*3.39
+(27/3.39^2-3.39)/3  =1.8452424531466833737959119743128E6
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-3-23 16:31:11 | 显示全部楼层
本帖最后由 落叶 于 2016-3-23 21:38 编辑

不好意思,没说清:这个是一百精度的,目前我设计的是一万精度,各种函数可以方便加上:下面是我的计算器运行效果,
(-(-1)(-2)-(-3) (-4)(-5)(-6)(-7)(-8)(-9)cuberoot(27)*3.39+(27/ 3.39^2-3.39)/3+sin(5))^999   =6.1321082202579006232825733921632263575771187476670971901178450674213340042091847280844090084002659067E6259
如果你需要处理比较复杂的算式,你怎么办,一个一个算吗?如果是高精度数呢?
比较先进的计算器是这样的吗?,x=****,y=***,然后x+Y=***.还有其它类似编程的功能,虽然工程师,科学家都很聪明,但是他们都会使用里面类似于程序员编程的功能吗?
这个程序可以加我们日常省略的运算符,去空格,去回车,你可能说这不难呀,上述做的意思,就是我已经建立了一个处理语法错误的程序框架,实现函数
的意义是:大多的运算可以通过函数表达,如果为了直观可以做一个格式转换函数,把一些公式等转换成内在的函数,必竟计算机并不认识公式,但可以
方便的识别函数。

具体你问我的需求,我没需求,只是我开始学写程序时(去年11月份左右,之前只看过书,几乎没动过手),准备写一个计算机程序练练手时,当时突发奇想做一个别人没做过的,
因为我用网上的计算器感到很不习惯,觉得要是能一起算多好啊,分开计算太麻烦,计算出一个得数还要想办法存,当时查到一个贴子,日本的一个手持计算器不能识别括号,
我想也许是手持计算器不好编程吧,然后上网找电脑版计算器,不过也是很失望,最后在excel的功能里面找到了表达式混合运算,并支持函数,但不支持高精度,
然后我就开始构思这个计算器,并加上了多精度功能。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2016-3-23 20:18:59 | 显示全部楼层
"(-(-1)(-2)-(-3) (-4)(-5)(-6)(-7)(-8)(-9)cuberoot(27)*3.39+(27/  3.39^2-3.39)/3+sin(5))^999   =6.1321082202579006232825733921632263575771187476670971901178450674213340042091847280844090084002659067E6259"

你上面的这个表达式语法很奇怪,不太符合一般的语法规则。
一般说来,计算机中的表达式和通常我们数学中的表达式还是有些区别的。在计算机表达式中,乘号是不能省略的。将上面的表达式补上乘号,并将cuberoot(27),替换为sqrtn(27,3)。在PARI/GP可以很容易算出来。
其结果为-3.288935247083505843230773751 E6560
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2019-12-10 03:28 , Processed in 0.073908 second(s), 21 queries .

Powered by Discuz! X3.4

© 2001-2017 Comsenz Inc.

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