5

I need to calculate the total number of divisors of a number N (without caring about what the values of the divisors are) and do so within 40-80 operations for all such numbers N. How can I do it? This is not a homework question. I tried out Pollard's Rho algorithm but even that turned out too slow for my purposes. Here is my code in python. How can I improve its performance, if possible?

def is_prime(n):    
    if n < 2:
        return False
    ps = [2,3,5,7,11,13,17,19,23,29,31,37,41,
         43,47,53,59,61,67,71,73,79,83,89,97]
    def is_spsp(n, a):
        d, s = n-1, 0
        while d%2 == 0:
            d /= 2; s += 1
        t = pow(int(a),int(d),int(n))
        if t == 1:
            return True
        while s > 0:
            if t == n-1:
                return True
            t = (t*t) % n
            s -= 1
        return False
    if n in ps: return True
    for p in ps:
        if not is_spsp(n,p):
            return False
    return True

def gcd(a,b):
        while b: a, b = b, a%b
        return abs(a)

def rho_factors(n, limit=100):
    def gcd(a,b):
        while b: a, b = b, a%b
        return abs(a)
    def rho_factor(n, c, limit):
        f = lambda x:    (x*x+c) % n
        t, h, d = 2, 2, 1
        while d == 1:
            if limit == 0:
                raise OverflowError('limit exceeded')
            t = f(t); h = f(f(h)); d = gcd(t-h, n)
        if d == n:
            return rho_factor(n, c+1, limit)
        if is_prime(d):
            return d
        return rho_factor(d, c+1, limit)
    if -1 <= n <= 1: return [n]
    if n < -1: return [-1] + rho_factors(-n, limit)
    fs = []
    while n % 2 == 0:
        n = n // 2; fs = fs + [2]
    if n == 1: return fs
    while not is_prime(n):
        f = rho_factor(n, 1, limit)
        n = int(n / f)
        fs = fs + [f]
    return sorted(fs + [n])

def divs(n):
    if(n==1):
        return 1
    ndiv=1
    f=rho_factors(n)
    l=len(f)
    #print(f)
    c=1
    for x in range(1,l):
        #print(f[x])
        if(f[x]==f[x-1]):
            c=c+1
        else:
            ndiv=ndiv*(c+1)
            c=1
       # print ("C",c,"ndiv",ndiv)
    ndiv=ndiv*(c+1)
    return ndiv
Kevin
  • 74,910
  • 12
  • 133
  • 166
thedarkknight
  • 147
  • 2
  • 8
  • 4
    Please clean up your code formatting. – Steven Rumbalski Sep 06 '12 at 17:43
  • Possibly related: [How can I write a fast function to calculate total divisors of a number?](http://stackoverflow.com/questions/12288671/how-can-i-write-a-fast-function-to-calculate-total-divisors-of-a-number) – Kevin Sep 06 '12 at 17:44
  • How big do you expect N to be? – andredor Sep 06 '12 at 17:55
  • 4
    A good start would be to code it in C :P – Haile Sep 06 '12 at 17:56
  • 3
    If this is a problem with some real-world code, can you explain why you need to do this? If this is an online programming challenge, can you provide a link to the problem description? If this is an interview question, does your interviewer know you're looking for help online? :-) – Kevin Sep 06 '12 at 17:56
  • @Kevin: I have looked through it seems the answers focused on finding the prime factors of the number ,I only want the total number of divisors of a number I dont want values of the divisors – thedarkknight Sep 06 '12 at 17:57
  • 1
    Also related: [What is the best way to get all the divisors of a number?](http://stackoverflow.com/q/171765) – Martijn Pieters Sep 06 '12 at 17:58
  • And: [Algorithm to calculate the number of divisors of a given number](http://stackoverflow.com/q/110344) – Martijn Pieters Sep 06 '12 at 17:59
  • @Kevin:It is a part of SPOJ problem I am solving ,it is no contest problem .It is not for any employment . – thedarkknight Sep 06 '12 at 17:59
  • @MartijnPieters:Please I dont need the divisors...I only need the total number of divisors – thedarkknight Sep 06 '12 at 18:00
  • Can you provide a link to the problem description? – Kevin Sep 06 '12 at 18:01
  • 1
    @user1652233: I know, that's why I am linking it as a related post, not marking this one as a duplicate. :-) – Martijn Pieters Sep 06 '12 at 18:01
  • @MartijnPieters:I have already tried out the answers written there ..none of the answers proved fast enough..they had too many operations for my solution to be fast enough ,I need 80 operations per test case on average for any N ... – thedarkknight Sep 06 '12 at 18:02
  • 3
    I'm pretty sure finding the number of divisors of any number N is very very hard to do efficiently. It's strictly harder than testing if any number N is prime or not, and that's the kind of problem that institutes offer million dollar prizes for solving. – Kevin Sep 06 '12 at 18:04
  • @Kevin:I hope you understand that I am a new programmer ,I am not looking for any prizes...Also isn't there any trick to reduce the number of operations for my code? – thedarkknight Sep 06 '12 at 18:09
  • 1
    There has to be a better way than brute force, right? (right?) I couldn't find anything, but it would somehow violate my faith in the order of the universe if you needed to factorize a number just to count the divisors. – harold Sep 06 '12 at 18:11
  • @Kevin:also I am storing the prime divisors ,then finding the total number of divisors using the divs(n) functions.I have a feeling it is not required as I basically just need the total number of divisors..these 2 above operations can somehow be merged.Just that i am not able to figure out how – thedarkknight Sep 06 '12 at 18:11
  • @user1652233 keep in mind that some SPOJ problems are VERY HARD to solve in slow languages like Python. You will not only need a good algorithm: for some problem you'll need a VERY DEEP knowledge on how to write fast and efficent python code. – Haile Sep 06 '12 at 18:17
  • @Haile:If you can suggest a great algorithm which reduces the number of operations,I will implement in C/C++/JAVA/Python no problem with that. – thedarkknight Sep 06 '12 at 18:21
  • Can you provide a link to the problem description? This isn't a trick question, and I'm not asking so I can say, "ah ha, it's a programming contest, I won't help you now". Reading the actual original problem would give us a deeper understanding of the problem and might explain why your current solution isn't working. – Kevin Sep 06 '12 at 18:21
  • @Haile:Do you think the solution could be faster in c++ ,if yes how much faster ..can you put an estimate to it,depending on that I will make the code in c/c++ – thedarkknight Sep 06 '12 at 18:22
  • @user1652233 It can be 10x up to 100x faster than python. (the same algorithm, if the python code is not well tuned). – Haile Sep 06 '12 at 18:23
  • Did you try to run this with PyPy? usually it's something like 5 times faster than CPython for highly numerical problems. Especially if those integers can be represented with plain ints, and not longs. – Bakuriu Sep 06 '12 at 18:31
  • What can be the maximum number of prime factors of a number upto 10^12 – thedarkknight Sep 06 '12 at 18:31
  • I think it's 12 (2*3*5*...*41 should be bigger than 10^12). By the way, why do you redefine gcd inside the rho_factors function? If you want to have the minimal speed-up of local variables you can do that without redefining it. Also because having closures it means that everytime you call rho_factors those functions are re-created[even though I think they are compiled only once]. – Bakuriu Sep 06 '12 at 18:35
  • "I need 80 operations per test case on average for any N" What counts as an operation? If any elementary arithmetic operation counts, you won't manage to do it in so few operations. – Daniel Fischer Sep 06 '12 at 18:39
  • @user1652233 The maximum number of prime divisors a number `<= 10^12` can have is 11, but that's pretty irrelevant, the numbers that take long are those with few (large) prime factors. – Daniel Fischer Sep 06 '12 at 18:42
  • @DanielFischer: I am basically trying to say that a complexity like o(log N) is required with a small constant factor. – thedarkknight Sep 06 '12 at 18:43
  • 1
    Perhaps there's something about the problem that you overlooked, that would reduce the problem space and therefore allow for a more efficient solution. _do you have a link to the original problem?_ – Kevin Sep 06 '12 at 18:48
  • take the prime list creation out of the is_prime function, this makes a new list every time this function runs! – unddoch Sep 06 '12 at 18:50
  • 2
    Well, bad news. Nobody has so far found a way to factorise numbers that is even polynomial in the number of digits (that'd be `O((log N)^k)` for some `k`). – Daniel Fischer Sep 06 '12 at 18:50
  • @DanielFischer:what I was trying to say that number of iterations for a particular number should not be upto square root n ....and a lot faster than that. – thedarkknight Sep 06 '12 at 18:58
  • @Kevin..I am not posting the link as it will divert the attention from this question itself ,I am quite sure that divisors need to be found and in the worst case the number whose divisors are to be found can be as large as 10^12 – thedarkknight Sep 06 '12 at 19:01
  • 4
    The usual method is to first remove small prime factors by trial division, and then switch to a method like Pollard's rho or elliptic curve factorisation, and hope for the best. For numbers in the given range, the rho algorithm finds the larger prime factors reasonably quickly, so if that's too slow after it's been optimised, you must be missing a trick. Could you please link to the problem statement? – Daniel Fischer Sep 06 '12 at 19:05
  • 1
    Do the divisors have to be unique? If not, 2^39 would have the most with 39 divisors. – Mark Ransom Sep 06 '12 at 19:06
  • 1
    @user1652233 go to the problem page, click on "best solutions", click on "python 2.7". How many solutions are there? – Haile Sep 06 '12 at 19:13
  • @DanielFischer:Is there any fast implementation of the Pollard rho algo for c++ which I can refer to. – thedarkknight Sep 06 '12 at 19:19
  • 1
    Undoubtedly, there will be several fast implementations, in C++ and some other languages. But I don't know what libraries would provide one. Google might turn some up. – Daniel Fischer Sep 06 '12 at 19:30
  • Wonder anyone here has experienced life as a low IQ person like me......I try hard but limited success – thedarkknight Sep 06 '12 at 19:32
  • 1
    @user1652233: if you're experimenting with Pollard's rho in your spare time, I wouldn't worry about low IQ. I might worry about low self-esteem, though! For my part, even though I have lots of letters after my name, I do very poorly on these kinds of problems. – DSM Sep 06 '12 at 19:43
  • may be http://stackoverflow.com/questions/12251962/prime-factorization-of-large-numbers helps – Wayne Rooney Sep 06 '12 at 21:27

3 Answers3

2

I remember solving this problem before on SPOJ, but I don't remember the exact method I used (it'd be great if you provide the problem ID). Did you try the naive method here? It has a complexity of O(sqrt n), which is about O(10 ^ 6) in the worst case. The modulo operator might be a bit slow, but it's worth giving it a try. Here's how it should look when done in C++:

int cntdiv = 0;
for(int i = 2; i * i <= x; i ++) {
    if(x % i == 0) {
        cntdiv += 1 + (i * i != x);
    }
}
//cntdiv is now your count
Chris
  • 26,544
  • 5
  • 58
  • 71
2

First, do you mean find the total number of divisors, the number of primes in its factorization, or the number of distinct prime divisors? For instance, 12 = 2 * 2 * 3 has 6 divisors (1,2,3,4,6,12), 3 primes in its factorization (2,2,3), and 2 distinct prime divisors (2,3). Do you want 6, 3, or 2 as your result? I'm going to assume you want the second of these for the rest of this, but I don't think anything materially changes if you were interested in one of the others...

Second, you're going to have to fully factor your number. There is no known shortcut that can find the number of prime factors without finding the factors themselves. (With the notable exception that you can test whether the number of factors is ==1 or >=2 quickly.)

10^12 is not that big. You only need to test divisors up to the square root of the number, which is at most 10^6. Say a divide takes 20 cycles on a modern CPU at 2GHz, that's only 10 milliseconds to test a million divisors.

#include <stdio.h>
#include <stdlib.h>
int main(int argc, char *argv[]) {
  long long n = atoll(argv[1]);
  for (int i = 2; i < 1000000; i++) {
    while (n % i == 0) { printf("%d\n", i); n /= i; }
  }
  if (n > 1) printf("%lld\n", n);
}

Takes 23 milliseconds on my machine. Wonder where that other 13 milliseconds went?

Python is about 10x slower, as this code still takes only 0.23 seconds on my machine:

import sys
n = int(sys.argv[1])
for i in xrange(2, 1000000):
  while n%i==0: print i; n/=i
if n>1: print n

How fast do you want it?

Keith Randall
  • 22,985
  • 2
  • 35
  • 54
  • :I only want total number of divisors of a number(both prime and composite) – thedarkknight Sep 07 '12 at 19:37
  • 1
    Then once you have the prime factorization, use http://stackoverflow.com/questions/110344/algorithm-to-calculate-the-number-of-divisors-of-a-given-number to calculate the total number of divisors. – Keith Randall Sep 07 '12 at 20:29
  • @keith Randall.....ok for example take n to be 99, then (99)^(1/2) gives 9 or 10(if u consider smallest integer function then calculate prime no.s upto 10 then what abt prime no. 11 which is also the factor of 99??? – Sumit Kumar Saha Sep 17 '13 at 17:31
  • @SumitKumarSaha: If a number is not prime, it will have *a* factor less than or equal to the square root of the number. It is not guaranteed that *all* factors are <= the square root, but you don't need that. 99=11*3*3, so you'll find both factors of 3 by checking divisibility by 2..9. Compute 99/3/3=11 as the final prime factor (if it is not 1). – Keith Randall Sep 18 '13 at 15:27
-1

I remember there is a solution based on sum of digits in a number and other features
As an example, 3411 is divisible by 9, because 3+4+1+1 = 9, sum of digits is divisible by 9, than a number is also divisible by 9. with other numbers rules are similar.

oleg.foreigner
  • 993
  • 2
  • 13
  • 28
  • @RobertHarvey actually there are variants in different bases (such as 16 which is more convenient) and for different divisors. They don't work nearly as well as just taking the mod with that number and testing whether it's zero, though. – harold Sep 07 '12 at 09:19
  • 1
    OK, come up with 78498 rules for me... (78498 = number of primes < 1000000) – Keith Randall Sep 07 '12 at 20:34