1

Forgive me if I am being a bit silly, but I have only very recently started programming, and am maybe a little out of my depth doing Problem 160 on Project Euler. I have made some attempts at solving it but it seems that going through 1tn numbers will take too long on any personal computer, so I guess I should be looking into the mathematics to find some short-cuts.

Project Euler Problem 160:

For any N, let f(N) be the last five digits before the trailing zeroes in N!. For example,

9! = 362880 so f(9)=36288 10! = 3628800 so f(10)=36288 20! = 2432902008176640000 so f(20)=17664

Find f(1,000,000,000,000)

New attempt:

#include <stdio.h>


main()
{
    //I have used long long ints everywhere to avoid possible multiplication errors
    long long f; //f is f(1,000,000,000,000)
    f = 1;
    for (long long i = 1; i <= 1000000000000; i = ++i){
        long long p;
        for (p = i; (p % 10) == 0; p = p / 10) //p is i without proceeding zeros
            ;
        p = (p % 1000000); //p is last six nontrivial digits of i
        for (f = f * p; (f % 10) == 0; f = f / 10)
            ;
        f = (f % 1000000); 
    }
    f = (f % 100000);
    printf("f(1,000,000,000,000) = %d\n", f);
}

Old attempt:

#include <stdio.h>

main()
{
    //This part of the programme removes the zeros in factorials by dividing by 10 for each factor of 5, and finds f(1,000,000,000,000) inductively
    long long int f, m; //f is f(n), m is 10^k for each multiple of 5
    short k; //Stores multiplicity of 5 for each multiple of 5
    f = 1;
    for (long long i = 1; i <= 100000000000; ++i){
        if ((i % 5) == 0){
            k = 1;
            for ((m = i / 5); (m % 5) == 0; m = m / 5) //Computes multiplicity of 5 in factorisation of i
                ++k;
            m = 1;
            for (short j = 1; j <= k; ++j) //Computes 10^k
                m = 10 * m;
            f = (((f * i) / m) % 100000);
        }
        else f = ((f * i) % 100000);
    }
    printf("f(1,000,000,000,000) = %d\n", f);
}
Martijn Pieters
  • 1,048,767
  • 296
  • 4,058
  • 3,343
und
  • 21
  • 3
  • Hmmm ... maybe try outputting o 0?? "f(1,000,000,000,000) = 09376" – pmg May 15 '14 at 18:06
  • 1
    @pmg: No, 09376 is not the correct answer :-) – Martin R May 15 '14 at 18:07
  • 3
    *"Since only the last five digits make any difference, f(1,000,000,000,000) is the same as f(100,000)^(10,000,000) mod 100,000"* ... That is probably not correct. – Martin R May 15 '14 at 18:10
  • how does that for loop raise f to the power 10m? – Red Alert May 15 '14 at 18:10
  • This does not relate to the power part (that bit doesn't as I'd hoped although I have not reasoned why), but the first part of my code (for finding f(n)) also seems to fail, particularly for large integers, like say 1,000,000,000. Putting in a printf in the for loop, I find that a lot of zeros are introduced to the ends of the f's, which is definitely not expected as there should be no fives in there at all. – und May 15 '14 at 19:22
  • Solve a simpler problem. Forget the last five nonzero digits; can you find the *last nonzero digit*? Now you can try your algorithm on very small cases that you can verify by hand. Try your algorithm *by hand* on these small cases; if your algorithm is flawed then you will find out that the problem is the algorithm and not some detail of its implementation. – Eric Lippert May 15 '14 at 21:19
  • Your theory is that `f(x^y)` is the same as `f(x)^y mod 100000`. Is that true for `f(2^3)` ? `f(2^3)` is 4032. `f(2)^3 mod 100000` is 8. So I think we have disconfirmed your mathematics here. – Eric Lippert May 15 '14 at 22:34
  • Indeed I realised after checking small cases that the second part of the programme is wrong. The reason is that I confused n! with n! minus the proceeding 0s. However, I am confused that the first part of the programme seems to work for all small integers I have tried, but not when I go to, say, 1bn. It computes but putting a printf in the loop shows that at some point it starts introducing zeros. Do you see any issues in the code that would cause this behaviour? – und May 16 '14 at 00:00
  • Suggested new approach: rather `%10==0`, count `%2==0` occurrences and each time `%5==0` happens, decrement the `%2==0` count. e.g. `for (term = 1; term <= 1000000000000ULL; term++) { uint64_t t = term; while (m2 < 64 && t%2==0) { t /= 2; m2++; } while (t%5==0) { t /= 5; m2--; } product *= t; product %= Ten5; } while (m2 > 0) { m2--; product = (product*2)%Ten5; }` – chux - Reinstate Monica May 16 '14 at 19:17
  • Here is a thing to ponder. If you loop 10^12 times on a machine which run at a few GHz, you expect to be going for well over 1000 seconds. I suspect that you are terminating too soon. Print `i` after your loop. – Floris May 17 '14 at 08:13
  • @Floris Sorry, I have changed the OP a couple of times so it does not make sense now. The algorithms above all take several hours to compute by me estimate, which is why I am still looking for a quicker way to do the problem. – und May 17 '14 at 13:50
  • Just to clarify: you expect the solution to the problem with N=1,000,000,000,000 (12 zeros) to take just a few seconds - so you can't loop over every number there is simply not enough time. That puts this in the math domain... – Floris May 17 '14 at 17:01
  • @Floris I see, thanks for letting me know. I realise that it's probably quite a basic fact that looping through 1 tn numbers will take too long, but I did not know this. In that case I will be looking for some mathematical short-cuts to reduce the computational element of the solution. – und May 17 '14 at 17:35
  • @und - just create a loop that does nothing but count to 1 trillion, and see how long it takes. I don't know if problem 160 has a time limit, but if it does I'm pretty sure a brute force solution will exceed it. That also suggests you might want to post the problem on math.stackexchange.com or other sites. This is not, at the heart of it, a programming problem. At least that's my hunch. – Floris May 17 '14 at 18:09
  • @und The approach I mentioned above took 8100 seconds with the result `16576`. Your results may vary. – chux - Reinstate Monica May 19 '14 at 15:54
  • If you want code to be 100000 times faster: This problem is really to count the number of times the numbers (1 to 1,000,000,000,000) have a factor of 2 or 5 in it. The number of factors in xxx00002 has the same number of 2s,5s in it regardless of xxx. The number of factors in xxx00006 has the same number of 2s,5s in it regardless of xxx. This and other shortcuts can be applied to most of the numbers so code need not iterate for every 5 digits ending. – chux - Reinstate Monica May 19 '14 at 16:11

5 Answers5

1

The problem is:

For any N, let f(N) be the last five digits before the trailing zeroes in N!. Find f(1,000,000,000,000)

Let's rephrase the question:

For any N, let g(N) be the last five digits before the trailing zeroes in N. For any N, let f(N) be g(N!). Find f(1,000,000,000,000).

Now, before you write the code, prove this assertion mathematically:

  • For any N > 1, f(N) is equal to g(f(N-1) * g(N))

Note that I have not proved this myself; I might be making a mistake here. (UPDATE: It appears to be wrong! We'll have to give this more thought.) Prove it to your satisfaction. You might want to start by proving some intermediate results, like:

  • g(x * y) = g(g(x) * g(y))

And so on.

Once you have obtained a proof of this result, now you have a recurrence relation that you can use to find any f(N), and the numbers you have to deal with don't ever get much larger than N.

Eric Lippert
  • 647,829
  • 179
  • 1,238
  • 2,067
  • Thank you for answering. I have not convinced myself that the relation holds. For it to hold it would have to rely on some property of factorials which is not obvious to me. For suppose we had 499!=100002 (obviously false, but I am just convincing myself that the proof cannot be general). Then under the assumption we have f(500)=00001, but 500*100002=50001000 so f(500) should be 50001. However, I did implement the algorithm in C and it seems to check out for small values, albeit quite slowly when I try 1bn (slower than my one above, but that one fails for reasons I have not understood). – und May 16 '14 at 01:00
  • @und: Good point; my recurrence relationship seems flawed. Is there a way to modify it to make it correct? – Eric Lippert May 16 '14 at 13:52
  • Indeed, I believe it should work using the last 6 digits in the worst case. However, my implementation of this algorithm unfortunately takes too long. I will edit the code into the OP so that someone can tell me if the code is non-optimal, or whether I should be looking for a different method of solving the problem. – und May 16 '14 at 14:48
1
Prod(n->k)(k*a+c) mod a <=> c^k mod a

For example

prod[ 3, 1000003, 2000003,... , 999999000003 ] mod 1000000

equals

3^(1,000,000,000,000/1,000,000) mod 1000000

And number of trailing 0 in N! equals to number of 5 in factorisation of N!

Luka Rahne
  • 10,336
  • 3
  • 34
  • 56
  • Nice approach but need to add that you need to get rid of the trailing zeros before mod !!! which is a little bit tricky but doable. – Spektre May 22 '14 at 06:57
  • there is problem with this ... try prod[10,1000010,2000010,...,999999000010] so you need to handle every 10th item differently – Spektre Jun 03 '14 at 15:57
  • only problem with this is that is each number divisible by 5. – Luka Rahne Jun 03 '14 at 16:26
0

I would compute the whole thing and then separate first nonzero digits from LSB ... but for you I think is better this:

1.use bigger base

  • any number can be rewrite as sum of multiplies of powers of the same number (base)

  • like 1234560004587786542 can be rewrite to base b=1000 000 000 like this:

    1*b^2 + 234560004*b^1 + 587786542*b^0

2.when you multiply then lower digit is dependent only on lowest digits of multiplied numbers

A*B = (a0*b^0+a1*b^1+...)*(b0*b^0+b1*b^1+...)
     = (a0*b0*b^0)+ (...*b^1) + (...*b^2)+ ...

3.put it together

for (f=1,i=1;i<=N;i++)
 {
 j=i%base;
 // here remove ending zeroes from j
 f*=j;
 // here remove ending zeroes from f
 f%=base;
 }
  • do not forget that variable f has to be big enough for base^2
  • and base has to be at least 2 digits bigger then 100000 to cover 5 digits and overflows to zero
  • base must be power of 10 to preserve decimal digits

[edit1] implementation

uint<2> f,i,j,n,base;   // mine 64bit unsigned ints (i use 32bit compiler/app)
base="10000000000";     // base >= 100000^2 ... must be as string to avoid 32bit trunc
n="20";                 // f(n) ... must be as string to avoid 32bit trunc
for (f=1,i=1;i<=n;i++)
    {
    j=i%base;
    for (;(j)&&((j%10).iszero());j/=10);
    f*=j;
    for (;(f)&&((f%10).iszero());f/=10);
    f%=base;
    }
f%=100000;
int s=f.a[1];           // export low 32bit part of 64bit uint (s is the result)

It is too slow :(

f(1000000)=12544 [17769.414 ms]
f(     20)=17664 [     0.122 ms]
f(     10)=36288 [     0.045 ms]

for more speed or use any fast factorial implementation

[edit2] just few more 32bit n! factorials for testing

this statement is not valid :(

//You could attempt to exploit that 
//f(n) = ( f(n%base) * (f(base)^floor(n/base)) )%base
//do not forget that this is true only if base fulfill the conditions above

luckily this one seems to be true :) but only if (a is much much bigger then b and a%base=0)

g((a+b)!)=g(g(a!)*g(b!))
// g mod base without last zeroes...
// this can speed up things a lot

f(                1)=00001
f(               10)=36288
f(              100)=16864
f(            1,000)=53472
f(           10,000)=79008
f(          100,000)=56096
f(        1,000,000)=12544
f(       10,000,000)=28125

f(        1,000,100)=42016
f(        1,000,100)=g(??????12544*??????16864)=g(??????42016)->42016
  • the more is a closer to b the less valid digits there are!!!
  • that is why f(1001000) will not work ...
Spektre
  • 49,595
  • 11
  • 110
  • 380
  • 1
    Forgive me, but what is the difference between using a bigger base to reduce computations, or using modular arithmetic to do what seems like the same thing? (since in this case only say the last six non-trivial digits matter at any point, there does not seem to be any need for more generality) – und May 16 '14 at 23:18
  • @und yes it is the same (+1) but you have base 10 times lower then it should be to avoid overflow errors. btw I had an typo in the loop f*=i instead of f*=j. – Spektre May 17 '14 at 06:33
  • @und btw you need 64bit uint variables but without optimization is this also too slow. try to implement fast factorial because I just tried this and I was expecting more speed ... f(1000000)=12544 [17769.414 ms], f(20)=17664 [0.122 ms], f(10)= 36288 [0.045 ms]. also make sure that your initialization of variables is OK and not truncated to 32 bit before load to 64 bit !!! look here: http://stackoverflow.com/a/18333853/2521214 – Spektre May 17 '14 at 07:00
  • @und had just tried mine fast factorial implementation with modulo 10,base and uint and it works. but still relatively slow (becuase uint<2> is template and can be much more optimized...) also you need all primes<=n for that and mine are just 32bit so i am limited to n<2^32 (I am too lazy to recode it to uint). runtimes: f(1000000) took arount 1.8sec and f(4000000) took around 7.8sec. (the primes are the slowest part of this but can be precomputed, ... factorial is ~O(n.log(n))) – Spektre May 19 '14 at 09:48
  • Your factorials are correct only up to 10000!. I know because I got a "Right answer" the brutal way (6 hours), but sadly failed to find a sweet solution. Even the brutal way was tricky, until I started discounting `5` (which included `10`) in the multiplier. My mistake was to tinker with the product, except for (and because of) the modulus. – Weather Vane Mar 27 '16 at 18:37
0

I'm not an expert project Euler solver, but some general advice for all Euler problems.

1 - Start by solving the problem in the most obvious way first. This may lead to insights for later attempts

2 - Work the problem for a smaller range. Euler usually give an answer for the smaller range that you can use to check your algorithm

3 - Scale up the problem and work out how the problem will scale, time-wise, as the problem gets bigger

4 - If the solution is going to take longer than a few minutes, it's time to check the algorithm and come up with a better way

5 - Remember that Euler problems always have an answer and rely on a combination of clever programming and clever mathematics

6 - A problem that has been solved by many people cannot be wrong, it's you that's wrong!

I recently solved the phidigital number problem (Euler's site is down, can't look up the number, it's quite recent at time of posting) using exactly these steps. My initial brute-force algorithm was going to take 60 hours, I took a look at the patterns solving to 1,000,000 showed and got the insight to find a solution that took 1.25s.

Tony van der Peet
  • 824
  • 6
  • 7
0

It might be an idea to deal with numbers ending 2,4,5,6,8,0 separately. Numbers ending 1,3,7,9 can not contribute to a trailing zeros. Let 

A(n) = 1 * 3 * 7 * 9 * 11 * 13 * 17 * 19 * ... * (n-1).
B(n) = 2 * 4 * 5 * 6 * 8 * 10 * 12 * 14 * 15 * 16 * 18 * 20 * ... * n.

The factorial of n is A(n)*B(n). We can find the last five digits of A(n) quite easily. First find A(100,000) MOD 100,000 we can make this easier by just doing multiplications mod 100,000. Note that A(200,000) MOD 100,000 is just A(100,000)*A(100,000) MOD 100,000 as 100,001 = 1 MOD 100,000 etc. So A(1,000,000,000,000) is just A(100,000)^10,000,000 MOD 100,000.

More care is needed with 2,4,5,6,8,0 you'll need to track when these add a trailing zero. Obviously whenever we multiply by numbers ending 2 or 5 we will end up with a zero. However there are cases when you can get two zeros 25*4 = 100.

Salix alba
  • 7,536
  • 2
  • 32
  • 38