views:

820

answers:

6

Can someone give me an idea of an efficient algorithm for large n (say 10^10) to find the sum of above series?

Mycode is getting klilled for n= 100000 and m=200000

#include<stdio.h>

int main() {
    int n,m,i,j,sum,t;
    scanf("%d%d",&n,&m);
    sum=0;
    for(i=1;i<=n;i++) {
     t=1;
     for(j=1;j<=i;j++)
      t=((long long)t*i)%m;
     sum=(sum+t)%m;
    }
    printf("%d\n",sum);

}
+19  A: 

Two notes:

(a + b + c) % m

is equivalent to

(a % m + b % m + c % m) % m

and

(a * b * c) % m

is equivalent to

((a % m) * (b % m) * (c % m)) % m

As a result, you can calculate each term using a recursive function in O(log p):

int expmod(int n, int p, int m) {
   if (p == 0) return 1;
   int nm = n % m;
   long long r = expmod(nm, p / 2, m);
   r = (r * r) % m;
   if (p % 2 == 0) return r;
   return (r * nm) % m;
}

And sum elements using a for loop:

long long r = 0;
for (int i = 1; i <= n; ++i)
    r = (r + expmod(i, i, m)) % m;

This algorithm is O(n log n).

Mehrdad Afshari
For small numbers I am doing the same thing. But it is getting killed for n = 1000000 and m =200000 . I have included the code
I believe you need one more 'mod' in those equations: '(a % m + b % m + c % m) % m', and '(a % m) * (b % m) * (c % m) % m'.
Groo
Groo: Yep, did that in code, missed in equations. Thanks. Fixed.
Mehrdad Afshari
@rahul: Your inner j-loop runs in O(i). Mehrdad's function runs in O(log i). Replace your inner loop by a call to Mehrdad's function and you will get a big speed-up.
Dave Hinton
@rahul. In case of n=1000 and m=2000, result 1000^1000%2000 to (1000^2%2000)^500%2000.
It is getting killed with expmod function also for n= 1000000000
rahul: No exceptions for me. It just takes a long time.
Mehrdad Afshari
+2  A: 

You may have a look at my answer to this post. The implementation there is slightly buggy, but the idea is there. The key strategy is to find x such that n^(x-1)<m and n^x>m and repeatedly reduce n^n%m to (n^x%m)^(n/x)*n^(n%x)%m. I am sure this strategy works.

+1  A: 

I think you can use Euler's theorem to avoid some exponentation, as phi(200000)=80000. Chinese remainder theorem might also help as it reduces the modulo.

Jaska
This does ring some bells, but I'm afraid you'll have to explain. Also, iirc, computing phi isn't trivial.
Kobi
You need to compute phi only once. Euler's theorem says that a^phi(b)=1 mod b if (a,b)=1. Then you can simplify a^c mod b to the form a^c' mod b where c'<phi(b).
Jaska
Jaska: It's irrelevant here. What if `(a,b) != 1`?
Mehrdad Afshari
It's not irrelevant as we can reduce some of the terms. Is (a,b)>1 then you should do the whole exponentiation. It is easy to check is a number is divisible by 2 or 5.
Jaska
Tip - try to edit your answer. Elaborate. Describe and explain your suggested algorithm. Try to post some code. Link to Wikipedia. Also, isn't the Chinese remainder theorem used for a set of equations?
Kobi
Euler's theorem and Chinese reminder theorem are easy to look up, and they are both (in conjunction) perfectly relevant here — use Euler's theorem to compute the sum mod each prime power in m, and use CRT to put them together.
ShreevatsaR
A: 

Are you getting killed here:

for(j=1;j<=i;j++)
    t=((long long)t*i)%m;

Exponentials mod m could be implemented using the sum of squares method.

n = 10000;
m = 20000;
sqr = n;
bit = n;
sum = 0;

while(bit > 0)
{
    if(bit % 2 == 1)
    {
        sum += sqr;
    }
    sqr = (sqr * sqr) % m;
    bit >>= 2;
}
Calyth
A: 

How about

n*(n+1)*(2*n+1)/6 (mod m)

This is the formula for the sum of squares.

tom10
@tom10, this isn't a sum of squares, it's a sum of n^ns.
DivineWolfwood
Oops... thanks. I was wondering why no one just wrote out the formula.
tom10
A: 

I can't add comment, but for the Chinese remainder theorem, see http://mathworld.wolfram.com/ChineseRemainderTheorem.html formulas (4)-(6).

Jaska
I cannot see how this helps
fortran