高精度计算的若干问题

高精度计算的若干问题•目录：•一、高精度计算的若干问题•1.高精度乘法与分治算法•2.牛顿迭代法•3.高精度除法•4.高精度开方•5.π值计算1.高精度乘法与分治算法高精度计算的核心运算是高精度乘法。加减法的计算复杂度是O(n)阶的，不会有阶上的改进，除法、开方等运算可转换为乘法。两个高精度数的乘法一般是O(n2)阶的，采用快速Fourier变换(FFT)或快速数论变换(FNT)可降至O(n*logn)阶,其算法较复杂。这里介绍一种分治算法，可将计算复杂度降至O(n1.58)阶。•设A,B是两个N(N=2n)位数，且•A=a110N/2+a2,B=b110N/2+b2•则：AB=(a110N/2+a2)(b110N/2+b2)•=a1b110N+(a1b2+a2b1)10N/2+a2b2•=p10N+(r–p–q)10N/2+q•其中，p=a1b1,q=a2b2,r=(a1+a2)(b1+b2)•这样，计算两个N位数的乘法AB，可转换为计算３个N/2位数的乘法(p,q,r)。•设Tn为N=2n时所做的１位乘法的次数，则Tn=3Tn-1=32Tn-2=…=3nT0=3n,•由于N=2n，n=log2N,•Tn=3logN=Nlog3≈N1.58•当N很大时，其效果还是很明显的。•具体处理见附录５。2.牛顿迭代法•给定方程:f(x)=0,已知x0是方程的一个近似根，f’(x)是f(x)的导函数,•如果f(x)是多项式：f(x)=anxn+an-1xn-1+…+a1x+a0,则•f’(x)=nanxn-1+(n-1)an-1xn-2+…+a2x+a0•求该方程的更精确的根的牛顿迭代公式是:•牛顿迭代公式在方程求根和高精度计算中具有重要的应用,当f’(xn)不为零时，对于适当的初值x0，牛顿迭代公式收敛，且具有二阶的收敛速度。即，如果设准确根为α,则误差|xn+1-α|≈C|xn-α|2,或者说，如果xn有k位有效数字，则xn+1有2k位有效数字。•．设x=1/b,b是高精度数,取•f(x)=bx2-x=0,f’(x)=2bx-1≈1,利用牛顿迭代公式，可得计算x的迭代公式为:•xn+1=xn-(bxn-xn)=xn(2-bxn)•说明：•1.若a,b都是高精度数,则a/b=a*(1/b),先按上述公式求出1/b,再计算a*(1/b)。•如果b不是高精度数,可归纳出有关算法，直接相除，不宜先计算1/b。3.高精度除法•2.牛顿法对收敛条件要求较高，建议事先乘以适当的因子,使b满足：0.1=b1。收敛结束的条件为|xn+1-xn|10-N，其中N为要计算的位数。•假定b,xn等均用长度为N+1的一维数组存放，b[0]存放b的整数部分，b[N]存放最后一位小数，xn+1，xn类似，为容纳舍入误差的积累，将数组体积延长到N+k,（例如，k=2~10)。为检验|xn+1-xn|10-N，只需对数组xn+1，xn的第N-r至第N位进行检验。(r：20-100)。•3.利用牛顿法二阶收敛的特性，可大致确定迭代次数，设初值x0有r位有效数字，迭代k次后，xk应有r*2k位有效数字。对于N位数，k可取为log2(N/r),迭代k次后,再检验条件|xn+1-xn|10-N。•在最初的几次迭代中，乘法不一定算满k位。为此可编制一个计算A(M)*B(N)=C(L)的子程序，要求M=N=L即可。•参考程序4.高精度开方•设x=1/，b是高精度数,取•f(x)=bx3-x=0,f’(x)=3bx2-1≈2,利用牛顿迭代公式，可得计算x的迭代公式为:•xn+1=xn-(bxn2-xn)/2=xn(3-bxn2)*0.5•于是，√b=b*(1/√b)=b*xb•附注：以上两个公式摘自D.H.Bailey1988年的一篇论文(Math.Comp.,Vol.50),该论文介绍了作者将π计算到29,360,000位小数的相关问题.5.π值计算•π这个数渗透了整个数学！•-------------陈省身•π值计算的世界纪录：•日本东京大学教授金田康正(Kanada)•2002年11月16日•计算出圆周率小数点后12411亿位数•(上次的纪录是他在1999年9月创造的，为2061亿位数）。π值计算公式•1.Machin公式•参看附录1239145116arctgarctg12)1...(75312753nxxxxxarctgxnnπ值计算公式（2）•下面两个著名公式不能用于π的高精度计算：)1665,(7566534431222)1671,(71513114WallisGregory•谢谢大家！