Sum of series: 1^1 + 2^2 + 3^3 + ... + n^n (mod m)

Question

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);

}

Aviator: Efficient algorithms are usually independent of the language. Shouldn't really matter whether this is Java or C (except maybe a linear factor in runtime). — Joey, Oct 01 '09 at 10:08
@Johannes: I understand. I thought of suggesting BigInteger. Thats why asked — vpram86, Oct 01 '09 at 10:09
You say you want something fast for big n (10^10), but you don't say whether m is similarly big, or if it stays around 200k. It might matter, because if m is small then you can try pre-calculating/caching some terms. If you already know a^m and a^a for all a less than m, then when you come to calculate (m+2)^(m+2) then it's just 2^(m+2) = 2^m*2^2. Then (m+3)^(m+3) = 3^m*3^3 and so on. You can probably arrange things so that you always access your stored values sequentially, not sure. — Steve Jessop, Oct 01 '09 at 11:59
Thinking about it, you might also want to cache 1^2m ... (m-1)^2m as well, while you're calculating the 2m+1 ... 3m-1 terms. Then use these values to calculate 1^3m ... (m-1)^3m, and replace the stored value with the new value for use in calculating 1^4m ... (m-1)^4m. Without writing the code I've no idea whether this will actually be faster than Mehrdad's solutino, but unless I've missed something fatal, it's O(n) instead of O(n log n). Obviously requires O(m) memory though. — Steve Jessop, Oct 01 '09 at 12:04
Oh, and O(m log m) time for calculating the cached values using Mehrdad's code. So if m grows with n, it's probably no help at all, but if m is bounded it's probably an improvement. — Steve Jessop, Oct 01 '09 at 12:09

Mehrdad Afshari · Accepted Answer · 2009-10-01T10:21:50.350

24

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).

edited Oct 01 '09 at 10:21

answered Oct 01 '09 at 09:59

Mehrdad Afshari

414,610
91
852
789

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 – Oct 01 '09 at 10:01
3

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'. – vgru Oct 01 '09 at 10:06
Groo: Yep, did that in code, missed in equations. Thanks. Fixed. – Mehrdad Afshari Oct 01 '09 at 10:08
4

@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. – dave4420 Oct 01 '09 at 10:10
@rahul. In case of n=1000 and m=2000, result 1000^1000%2000 to (1000^2%2000)^500%2000. – user172818 Oct 01 '09 at 10:11
It is getting killed with expmod function also for n= 1000000000 – Oct 01 '09 at 10:44
rahul: No exceptions for me. It just takes a long time. – Mehrdad Afshari Oct 01 '09 at 10:58

score 5 · Answer 2 · 2009-10-01T10:47:19.973

5

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.

edited Oct 01 '09 at 10:47

answered Oct 01 '09 at 10:21

This does ring some bells, but I'm afraid you'll have to explain. Also, iirc, computing phi isn't trivial. – Kobi Oct 01 '09 at 10:26
2

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' – Oct 01 '09 at 10:35
Jaska: It's irrelevant here. What if `(a,b) != 1`? – Mehrdad Afshari Oct 01 '09 at 10:43
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. – Oct 01 '09 at 10:49
2

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 Oct 01 '09 at 11:04
2

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 Oct 01 '09 at 22:52

score 3 · Answer 3 · edited May 23 '17 at 12:24

3

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.

edited May 23 '17 at 12:24

Community

1
1

answered Oct 01 '09 at 10:09

user172818

4,518
1
18
20

score 1 · Answer 4 · answered Oct 11 '12 at 06:26

I encountered similar question recently: my 'n' is 1435, 'm' is 10^10. Here is my solution (C#):

ulong n = 1435, s = 0, mod = 0;
mod = ulong.Parse(Math.Pow(10, 10).ToString());
for (ulong i = 1; i <= n; 
{
     ulong summand = i;
     for (ulong j = 2; j <= i; j++)
     {
         summand *= i;
         summand = summand % mod;
     }
     s += summand;
     s = s % mod;
}

At the end 's' is equal to required number.

score 0 · Answer 5 · answered Oct 01 '09 at 12:57

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;
}

score 0 · Answer 6 · answered Oct 01 '09 at 18:29

0

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

answered Oct 01 '09 at 18:29

Jaska

775
1
8
13

Sum of series: 1^1 + 2^2 + 3^3 + ... + n^n (mod m)

6 Answers6

Linked

Related