Project Euler+ #97 modular exponentation not working

Question

I am trying to solve Project Euler+ #97 from Hackerrank. The problem asks to calculate the last 12 digits of A x B ** C + D. My attempt was to use the modular exponentiation mod 10 ** 12 from Wikipedia in order to efficiently calculate the last 12 digits and avoid overflow. However, for all cases aside from the sample 2 x 3 ** 4 + 5 I am getting wrong. According to the constraints there should be no overflow for unsigned long long.

---

The problem:

Now we want to learn how to calculate some last digits of such big numbers. Let's assume we have a lot of numbers A x B ** C + D and we want to know last 12 digits of these numbers.

Constraints:

• 1 ≤ T ≤ 500000
• 1 ≤ A, B, C, D ≤ 10 ** 9

Input: First line contains one integer T - the number of tests. T lines follow containing 4 integers (A, B, C and D) each.

Output: Output exactly one line containing exactly 12 digits - the last 12 digits of the sum of all results. If the sum is less than 10 ** 12 print corresponding number of leading zeroes then.

---

My attempt in C

#include <stdio.h>

int main() {
    const unsigned long long limit = 1000000000000;
    int cases;
    for (scanf("%d", &cases); cases; --cases) {
        // mult = A, base = B, exp = C, add = D
        unsigned long long mult, base, exp, add;
        scanf("%llu %llu %llu %llu", &mult, &base, &exp, &add);
        base = base % limit;
        while (exp) {
            if (exp & 1) {
               mult = (mult * base) % limit;
            }
            exp >>= 1;
            base = (base * base) % limit;
        }
        printf("%012llu\n", (mult + add) % limit);
    }
    return 0;
}

Answer 1

I think you can overflow unsigned long long math (e.g. - modulo 2^64) because your computation of base in your inner loop can get as high as (10^12 - 1)^2 ~= 10^24 ~= 2^79.726, which is much more than 2^64. For example, think about B = 10^6 - 1 and C = 4.

On my MacBook Pro running a 64b version of Mac OS X with clang 8.1.0:

#include <stdio.h>

int main()
{
  fprintf(stdout, "sizeof(unsigned long long) = %u\n", (unsigned) sizeof(unsigned long long));
  fprintf(stdout, "sizeof(__uint128_t) = %u\n", (unsigned) sizeof(__uint128_t));
  fprintf(stdout, "sizeof(long double) = %u\n", (unsigned) sizeof(long double));

  return 0;
}

Says:

sizeof(unsigned long long) = 8
sizeof(__uint128_t) = 16
sizeof(long double) = 16

If your platform says 16 or 10 instead for long long, then I think you are in the clear. If it says 8 like mine does, then you need to rework your answer to either force 128b (or 80b) integer math natively or mimic it some other way.

You can try __uint128_t, which is supported by gcc and clang. Otherwise, you'd need to resort to something like long double and fmodl(), which might have enough mantissa bits but might not give exact answers like you want.

Also, you don't accumulate multiple results like the task says. Here's my shot at it, based on your program, but using __uint128_t.

#include <stdio.h>
#include <stdlib.h>

#define BILLION     1000000000
#define TRILLION 1000000000000

int main() 
{
  const __uint128_t limit = TRILLION;
  unsigned long     cases = 0;
  __uint128_t       acc   = 0;

  if (scanf("%lu", &cases) != 1 || cases == 0 || cases > 500000)
    abort();

  while (cases-- > 0)
  {            
    unsigned long a, b, c, d;
    __uint128_t b2c = 1, bbase;

    if (scanf("%lu %lu %lu %lu", &a, &b, &c, &d) != 4 ||
        a == 0 || a > BILLION || b == 0 || b > BILLION || 
        c == 0 || c > BILLION || d == 0 || d > BILLION)
      abort();

    for (bbase = b; c > 0; c >>= 1) 
    {
      if ((c & 0x1) != 0)
        b2c = (b2c * bbase) % limit;    // 64b overflow: ~10^12 * ~10^12 ~= 10^24 > 2^64

      bbase = (bbase * bbase) % limit;  // same overflow issue as above
    }

    // can do modulus on acc only once at end of program instead because
    // 5 * 10^5 * (10^9 * (10^12 - 1) + 10^9) = 5 * 10^26 < 2^128

    acc += a * b2c + d;  
  }

  acc %= limit;

  printf("%012llu\n", (unsigned long long) acc);

  return 0;
}