找回密码
 欢迎注册
查看: 12221|回复: 3

[原创] 高精度乘方

[复制链接]
发表于 2016-6-17 12:48:14 | 显示全部楼层 |阅读模式

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

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

×
本帖最后由 落叶 于 2016-6-17 12:56 编辑

先说说整指数乘方,百度上搜到一个程序,效率很高,代码简洁:
int pow(int x,int n)
{
int temp(x),res(1);
while(n)
{
if(n&1)
{
res *= temp;
}
temp *= temp;
n>>=1;
}
return res;
}
把2^61代入运算,内存变化图是这样的:61的二进制形式是111101
res    =       2      2        32      8192        536870912       2.305843009213693952E18                           
temp=        4    16       256   65536     4294967296        1.8446744073709551616E19
2.305843009213693952E18 是最终结果 除红色部分外有9次是关键乘法,前两个红色 2  2近似不占运算时间,后面红色1.8446744073709551616E19是冗余运算。
我的源程序和上述算法效率差不多,因构思不一样,也贴出来
Dim a(50) As Integer
i = 0
        Do
            yCzs3 = yCzs.Zx(1) Mod 2                                         ‘yCzs.Zx(1)中存放的是指数
            If yCzs3 = 0 Then
                yCzs.Zx(1) = yCzs.Zx(1) / 2
                a(i) = 2
            Else
                yCzs.Zx(1) = yCzs.Zx(1) - 1
                a(i) = 1
            End If
            yCzs3 = yCzs.Zx(1)
            i = i + 1
        Loop While yCzs.Zx(1) <> 1                '等于一退出                 ’循环完成后a()的内存情况为:1 2 2 1 2 1 2 1 2

        For j = i - 1 To 0 Step -1                                                     ‘一共需要9次乘法                  
            If a(j) = 2 Then
                zCzs1 = myMULTIPLY(zCzs1, zCzs1)                          ’ zCzs1的初始值是2,myMULTIPLY()函数是大数乘
           Else
                zCzs1 = myMULTIPLY(zCzs1, zCzs2)                            ’ zCzs2的初始值是2
            End If
        Next
我的思路是这样的:((((((((2^2)*2)^2)*2)^2)*2)^2)^2)*2


毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-6-17 13:57:24 | 显示全部楼层
本帖最后由 落叶 于 2016-6-17 19:15 编辑

上面只是说到整指数,如果指数是小数,情况就复杂多了,
有一种方法是先把小数化成分数
如2的1.2次方就等于2的5分之6次方,等价于2的6次方,再开5次方,yroot(5,2^6)  =2.2973967099940700135972538935559
这个方法缺点很大,如2^0.1111111111111111等于2的10000000000000000次方,然后再开1111111111111111次方
这么大的乘方,开方,效率低下,而且难以实现,没有这么大的高精开N次方程序,
不过有了这个公式 :a^x= 1+ax+a*(a-1)*x^2/2!+a*(a-1)*(a-2)*x^3/3!+...+a*(a-1)*(a-2)*...(a-n+1)x^n/n!
这个公式有两个难点,a的大小要接近1,x的大小要接近0
我的方法如下:例:123.45^67.891
123.45^67.891可分为(1.2345^67.891)*(10^2)^67.891=(1.2345^67.891)*(10)^(2*67.891)
对1.2345,和10分别开2^80次方,此时底数a 的形式变为1.0000.......xxx,当然结果返回后针对前面的开方要做后期处理。
对指数处理:67.891=67+0.891,指数67可以用前面的整指数乘方程序算,
然后把分解后的a 和x(这里是0.891)代入泰勒公式计算,
计算完成后要针对前面的开方作后期处理
开方后期处理是:

              string1 = myCf(string1, 2^80),    myCf是整数乘方函数

因程序太复杂,中间相互调用多,源程序也难以理解,网友可以了解了解大致思路,有好的意见请指点。
上述思路已通过程序实现,证明是正确的。              
上面的泰勒公式化简为:a^x= 1+ ax*(2*3*4*5+ (a-1)*x*(3*4*5 + (a-2)*x*(4*5+(a-3 )*x*(5+(a-4)*x....))))/n!   这种方式的化简对公式的运算提速是明显的。                 
               
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2016-6-18 09:56:43 | 显示全部楼层
乘方容易,若一个数>=2^b 而<2^(b+1),则之多需要2b次乘法。
对于开方,我觉得只需实现开平方和开立方。 其他次的开方和非整数次幂的乘法,用指数和对数来实现即可。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2016-6-19 20:25:10 | 显示全部楼层
本帖最后由 落叶 于 2016-6-19 22:35 编辑

若一个数>=2^b 而<2^(b+1),则之多需要2b次乘法,这个说法很精确,超过了这个数,代表指数分解不彻底.
对于如何用开平方,和开立方,分解或化简开N次方,我曾试过多次,一直没成功.
我查了百度,大都在线计算器开N次方算法是把,例如:n√x换成x^(1/n),可以理解为是通过小数指数的乘方算法(这个算法应该不包括调用开方程序,否则形成了相互调用)实现开N次方.
我的程序实现了两种开N次方程序,一种是用迭代法实现,优点是速度较快,但N的取值不能太大,我的程序精度一万时,只能支持最大开25500次方,
第二种方法:就是通过转换成上面所述的小数指数乘方实现,速度比第一种方法慢,但理论上可以开任意次方.
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

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

GMT+8, 2024-12-4 01:15 , Processed in 0.023305 second(s), 16 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

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