探索c++的函数pow()的实现方法·数学与程序设计的结合

 

计算t的m次幂的方法：(探索c++的函数pow()的实现方法，数学与程序设计的结合)

欢迎看看另一篇文章《实现log()和exp()函数的方法,并以此计算pow() 》增加了对Taylor展开后的进一步处理。
源程序下载：http://pjy.studio.googlepages.com/powP.cpp
或到我的收藏下载。   
    由于c++刚学完函数一章，而练习需要用pow()这个函数，于是就特发奇想，想自己能否写一个能实现pow()功能的函数，经过一段努力，算有了一些结果。
    众所周知，pow(double t,double m)是c++提供计算x的y次幂的函数，虽然系统提供了这个pow(),但我还是想自己写一个自己的pow()。不过要写出这个pow可能不太容易，因为 指数m要求是double的，即可以是小数，那就不是简单的循环可以做出的。
    刚开始想这个问题时，第一个想法是把y的整数部分和小数部分分离，例如若m＝12.34，则把x^12.34转化为(x^12)*(x^34)/ (x^100),这样就把指数为double的情况转化为三个指数为int的情况。而取其整数，最简单的方法可能是int(y)了，而去小数则用m- int(m)即可，然后把小数按例子那样分离。不过这里可能存在两个问题：
    1.int(y)可能会造成意想不到的结果，例如数据掉失。
    2.把一个浮点数分成三部份的工作比较麻烦。
    以上两点，再加上我估计c++不大可能会用这种方法，而应该从数学方法去入手，所以我就另觅途径了。说来也是比较巧合，我的数学分析课刚好学完了幂级数这一章，我就有了用把t^m展开成幂级数这一方向入手的想法，而刚好书上有这个例子： 
由幂级数的知识推出当-1<x<1，n->无穷(n为整数)时恒有：(1+x)^m=1+mx+....[m*(m-1)*..*(m-n+1)*(x^n)]/(n!)+...
（具体证明不列出了）
并把通项记为tmpm，即tmpm(n)＝[m*(m-1)*..*(m-n+1)*(x^n)]/(n!)
    由于-1<x<1，那么用这个方法只能算出0～2的m次幂，不过如果1+x>=2,那么可以把它转为求1/(1+x)的-m次幂即可。 综上述，也就可以求出正数t的m次幂了（把t^m看为(1+x)^m）。估计c++的pow可能也是用类似的方法实现的。

    由公式，首先我们需要求x^n,即一个数的整数次幂，这个也不难：
double pow_i(double num,int n)//计算num的n次幂，其中n为整数 
{
   double powint=1;
   int i;
   for(i=1;i<=n;i++) powint*=num;
   return powint;
}
    而后，我们还应该需要一个求n!的函数。
最后就可以写出计算展开式的函数了：
double pow_f(double num,double m)//计算num的m次幂，num和m可为双精度，num大于零 
{
    int i,j;
    double powf=0,x,tmpm=1;
    x=num-1;
    for(i=1;tmpm>1e-12 || tmpm<-1e-12;i++)//当tmpm不在次范围时，停止循环,范围可改 
           {
           for(j=1,tmpm=1;j<=i;j++) 
                tmpm*=(m-j+1)*x/j;
                powf+=tmpm;
           }
    return powf+1;
}
注意到，我在实现求每一个n对应的tmpm时，并没有用到pow_i()，也没有用到一个函数来求n!。而是直接写成tmpm*=(m-j+1)*x/j，并通过循环此式来求n对应的tmpm。这是有原因的：
    1.这样写更加简洁。
    2.如果用函数返回n!时，会出现循环只会进行171次的后果，因为171!已超过double的范围。这样有可能造成因为循环不足而导致结果精确度差，而且有可能运行时出错。
    
    到这里，利用pow_f已能求出底数不大于2，指数为实数的幂了。容易想到，对于底数t>=2的情况，可以把它转为求1/(1+x)的-m次幂即 可。理论上是如此，但实际上机却显示一个问题，就是当t比较小时（<1），而m比较大时，误差十分大，这可能是因为循环不足而造成的。也就是说，这 个pow_f不能求出t<1，m>1时比较精确的值，那样的话对与t>2的情况也就求不出了，因为作1/t的转化后，又变成求底数小于 1的幂问题了。
    怎样解决呢，我们说，是因为t较小而m较大时，才出问题，这样的话，不妨把m分开为整数和小数（分别记为mi，mf），各自对t求幂，再把两个结果相乘不 就行了吗，这样对t^mi,可以用pow_i()，对t^mf，可以用pow_f,此时mf<1,用pow_f求t^mf，结果很准确。这个问题解 决了，那么对于t>=2的问题也就解决了。
    最后，还有一个问题，就是t<0怎么办，在这里，我的处理是，当t<0是，若指数m为整数，则用pow_i求，而m不是整数，则返回0（c++的pow()好似是用错误处理的，不过我才是个初学者，没学到，所以就返回个0吧）。
    
    综上述，我们还需要有一个函数来作为入口，在调用pow_i和pow_f前处理这些情况：
double pow_ff(double num,double m)//调用pow_f()和pow_i(),计算num的m次幂,是计算幂的入口
{
    if(num==0 && m!=0) return 0;//若num为0，则返回0
    else if(num==0 && m==0) return 1;// 若num和m都为0，则返回1
    else if(num<0 && m-int(m)!=0) return 0;//若num为负，且m不为整数数，则出错，返回0 
    if(num>2)//把底数大于2的情况转为(1/num)^-m计算
    {
        num=1/num;
        m=-m;
        }
    if(m<0) return 1/pow_ff(num,-m);//把指数小于0的情况转为1/num^-m计算
    if(m-int(m)==0) return pow_i(num,m);/*当指数为浮点数是，分成整数和小数分别求
                                        幂，这是因为但底数较小式，用pow_f直接求幂
                                        误差大，所以分为指数的整数部分用pow_i,小
                                        数部分用pow_f求.*/ 
    else return pow_f(num,m-int(m))*pow_i(num,int(m));
    return pow_f(num,m);
}

这样，这个pow_ff基本上就和c++的pow差不多了，至于系统的pow究竟是怎么实现的，是否也是这样做的，就不太清楚了。而这此的练习，使我感到数学与编程结合的美妙，^_^。