代码之家  ›  专栏  ›  技术社区  ›  marc lincoln

计算谐波级数的python程序

  •  5
  • marc lincoln  · 技术社区  · 16 年前

    有谁知道如何用python编写一个程序来计算调和级数的加法。即1+1/2+1/3+1/4…

    10 回复  |  直到 8 年前
        1
  •  21
  •   Community CDub    8 年前

    @Kiv's answer 是正确的,但如果不需要无限精度,则对于大n来说速度很慢。最好使用 asymptotic formula 在这种情况下:

    asymptotic expansion for harmonic number

    #!/usr/bin/env python
    from math import log
    
    def H(n):
        """Returns an approximate value of n-th harmonic number.
    
           http://en.wikipedia.org/wiki/Harmonic_number
        """
        # Euler-Mascheroni constant
        gamma = 0.57721566490153286060651209008240243104215933593992
        return gamma + log(n) + 0.5/n - 1./(12*n**2) + 1./(120*n**4)
    

    @基夫的回答 对于Python2.6:

    from fractions import Fraction
    
    harmonic_number = lambda n: sum(Fraction(1, d) for d in xrange(1, n+1))
    

    例子:

    >>> N = 100
    >>> h_exact = harmonic_number(N)
    >>> h = H(N)
    >>> rel_err = (abs(h - h_exact) / h_exact)
    >>> print n, "%r" % h, "%.2g" % rel_err
    100 5.1873775176396242 6.8e-16
    

    AT N = 100 相对误差小于 1e-15 .

        2
  •  12
  •   Community CDub    8 年前

    @recursive's solution 对于浮点近似值是正确的。如果愿意,可以使用分数模块在Python3.0中获得确切的答案:

    >>> from fractions import Fraction
    >>> def calc_harmonic(n):
    ...   return sum(Fraction(1, d) for d in range(1, n + 1))
    ...
    >>> calc_harmonic(20) # sum of the first 20 terms
    Fraction(55835135, 15519504)
    

    请注意,数字的数量增长很快,因此这将需要大量的内存来存储大n。如果您真的想得到更多信息,还可以使用生成器来查看部分和的序列。

        3
  •  5
  •   joel.neely    16 年前

    只是在其他使用浮点的答案上的一个脚注;从最大除数开始并迭代 向下地 (朝向最大值的倒数)将尽可能推迟累积的舍入误差。

        4
  •  4
  •   dancavallaro    16 年前

    谐波级数发散,即其和是无穷大。

    编辑:除非你想要部分和,但你不太清楚。

        5
  •  4
  •   FutureNerd    10 年前

    一个快速、准确、平滑、复值的h函数可以用digamma函数来计算,如前所述 here . euler-mascheroni(gamma)常数和digamma函数分别在numpy和scipy库中可用。

    from numpy import euler_gamma
    from scipy.special import digamma
    
    def digamma_H(s):
        """ If s is complex the result becomes complex. """
        return digamma(s + 1) + euler_gamma
    
    from fractions import Fraction
    
    def Kiv_H(n):
        return sum(Fraction(1, d) for d in xrange(1, n + 1))
    
    def J_F_Sebastian_H(n):
        return euler_gamma + log(n) + 0.5/n - 1./(12*n**2) + 1./(120*n**4)
    
    
    

    下面是速度和精度的三种方法的比较(以kiv_h为参考):

    Kiv_H(x) J_F_Sebastian_H(x) digamma_H(x) x seconds bits seconds bits seconds bits 1 5.06e-05 exact 2.47e-06 8.8 1.16e-05 exact 10 4.45e-04 exact 3.25e-06 29.5 1.17e-05 52.6 100 7.64e-03 exact 3.65e-06 50.4 1.17e-05 exact 1000 7.62e-01 exact 5.92e-06 52.9 1.19e-05 exact

        6
  •  2
  •   recursive    16 年前

    这应该能起到作用。

    def calc_harmonic(n):
        return sum(1.0/d for d in range(2,n+1))
    
        7
  •  0
  •   zenazn    16 年前

    这个怎么样:

    partialsum = 0
    for i in xrange(1,1000000):
        partialsum += 1.0 / i
    print partialsum
    

    其中1000000是上限。

        8
  •  0
  •   duffymo    16 年前

    作业?

    这是一个发散级数,所以不可能对所有项求和。

    我不知道Python,但我知道如何用Java编写它。

    public class Harmonic
    {
        private static final int DEFAULT_NUM_TERMS = 10;
    
        public static void main(String[] args)
        {
            int numTerms = ((args.length > 0) ? Integer.parseInt(args[0]) : DEFAULT_NUM_TERMS);
    
            System.out.println("sum of " + numTerms + " terms=" + sum(numTerms));
         }
    
         public static double sum(int numTerms)
         {
             double sum = 0.0;
    
             if (numTerms > 0)
             {
                 for (int k = 1; k <= numTerms; ++k)
                 {
                     sum += 1.0/k;
                 }
             }
    
             return sum;
         }
     }
    
        9
  •  0
  •   georgegkas    8 年前

    我添加了另一个解决方案,这次使用递归,以找到第n个谐波数。

    一般实施细节

    功能原型: harmonic_recursive(n)

    功能参数: n -第n次谐波数

    基本情况: 如果 n 等于 1 返回1。

    重复步骤: 如果不是基本情况,请致电 harmonic_recursive 对于 n-1 并将结果加上 1/n . 这样我们每次都把调和级数的第i项加上之前所有项的和,直到那一点。

    伪码

    (这个解决方案也可以用其他语言轻松实现。)

    harmonic_recursive(n):
        if n == 1:
            return 1
        else:
            return 1/n + harmonic_recursive(n-1)
    

    python代码

    def harmonic_recursive(n):
        if n == 1:
            return 1
        else:
            return 1.0/n + harmonic_recursive(n-1)
    
        10
  •  -1
  •   srikar saggurthi    8 年前

    使用simple for循环

    def harmonicNumber(n):
    x=0
    for i in range (0,n):
        x=x+ 1/(i+1)
    return x