Friday, 12 September 2008

1+2+...+n的约数(5) - 再回首:一切竟是那么简单

嗯,优化的算法找好了,让我们替换原来的求质函数吧...
等等,我们需要计算n(n+1)/2个质数吗?
汗,一滴滴流了下来。
不需要!答案是勿庸置疑的!

写计算n个质数的算法...
等等,我们需要计算n个质数吗?
又是不需要!
求解sqrt(n)就够了...

天啊,我竟然用了这么多无用的计算...

把新的求质函数放在一旁,改进现有算法:

代码:
class Divisor
  # prime: [2, 3, 5, 7, ...]
  #attr_accessor :primes, :offset
  
  def initialize
    @primes = Array.new(1, 2)
    @offset = 3
  end
  
  def n(num)
    return 1 if (num == 1)
    
    n = 2
    n += 1 while (num_of_divs(n) < num)
    
    return n
  end
  
  # add n if it is a prime
  def add_prime(n)
    sqrt = Math.sqrt(n).to_i    
    @primes.each do |prime|
      break if (prime > sqrt)
      return if (n % prime == 0)
    end
    
    @primes << n
  end
  
  def num_of_divs(n)
    n1 = n
    n2 = n +1
    sqrt = Math.sqrt(n2).to_i
    
    # if not added (into primes' list), try to add it,
    if (sqrt >= @offset)
      add_prime(sqrt) 
      @offset += 1  #increase prime offset
    end
    
    # f(n) = n(n+1)/2 could be as multiplication of two numbers, either (n/2) * (n+1) or n * [(n+1)/2]
    if (n1 % 2 == 0)
      n1 /= 2
    else
      n2 /= 2
    end
    count = 1
    @primes.each do |prime|
      # exit the loop when both n1 & n2 are 1
      break if (n1 == 1 && n2 == 1)
      
      # try to use current prime divide n1 & n2
      tmp = 0
      while (n1 % prime == 0)
        tmp += 1
        n1 /= prime
      end
      while (n2 % prime == 0)
        tmp += 1
        n2 /= prime
      end
      
      count *= (tmp + 1) if (tmp > 0)
    end
    
    #if n1 or n2 is not 1, multiple the count by 2 (as current n1 or n2 is a prime)
    count *= 2 if (not n1 == 1)
    count *= 2 if (not n2 == 1)
    
    return count
  end
end
start = Time.now
div = Divisor.new
print '->' + div.n(500).to_s + '<-'
final = Time.now
cost =final - start
print 'Cost: ' + cost.to_s
100个约数需要0.016s,500个约数需要0.828s...
原来...一切都是这么简单...

No comments:

Post a Comment