## Do you have any efficient algorithm?

Arrangements, combinations and permutations, probability, ...
square1001
Posts: 27
Joined: Tue Mar 15, 2016 2:58 am
Location: Tokyo, Japan
Contact:

### Do you have any efficient algorithm?

I have a question:
You are given an integer n.
Find $n! \mod 1,000,000,007$.
I can calculate for $O(n)$ for a naive calculation.

But if $n \le 1,000,000,000$, I can't calculate for a few seconds.

Do you have an "efficient" algorithm?
If you have, please post.

sjhillier
Posts: 525
Joined: Sun Aug 17, 2014 4:59 pm
Location: Birmingham, UK
Contact:

### Re: Do you have any efficient algorithm?

It may be giving away too much to answer that directly, but how about trying to find $1,000,000,006!$ mod $1,000,000,007$, then think about the answer, and try to think how that might help you.

square1001
Posts: 27
Joined: Tue Mar 15, 2016 2:58 am
Location: Tokyo, Japan
Contact:

### Re: Do you have any efficient algorithm?

We can prove that $1,000,000,006! \equiv 1,000,000,006 \text{ } (\text{mod } 1,000,000,007)$.

So, if you want to calculate for big n, you can use mod inverse like this code:

Code: Select all

#include <iostream>
using namespace std;
int modpow(int a, int b, int m) {
int ret = 1;
for(int i = 30; i >= 0; i--) {
ret = 1LL * ret * ret % m;
if(b & (1 << i)) ret = 1LL * ret * a % m;
}
return ret;
}
int n; const int mod = 1000000007;
int main() {
cin >> n;
int ret = 1;
if(n < 970000000) {
for(int i = 1; i <= n; i++) ret = 1LL * ret * i % mod;
}
else {
ret = mod - 1;
for(int i = mod - 1; i > n; i--) ret = 1LL * ret * modpow(i, mod - 2, mod) % mod;
}
cout << ret << endl;
return 0;
}

The time complexity is about $O\left(\frac{n\log n}{\log n + 1}\right)$. It improved a little!

But the time is about 30-40 second in the worst case.

Do you have a more efficient algorithm?

sjhillier
Posts: 525
Joined: Sun Aug 17, 2014 4:59 pm
Location: Birmingham, UK
Contact:

### Re: Do you have any efficient algorithm?

Yes, that's the method I was thinking of, but as you say it only helps in some cases. Anything beyond that is also beyond my knowledge right now!

square1001
Posts: 27
Joined: Tue Mar 15, 2016 2:58 am
Location: Tokyo, Japan
Contact:

### Re: Do you have any efficient algorithm?

I improved this solution and it is the fast $O(n)$ algorithm.

This is my C++ code:

Code: Select all

#include <iostream>
using namespace std;
int modpow(int a, int b, int m) {
int ret = 1;
for(int i = 30; i >= 0; i--) {
ret = 1LL * ret * ret % m;
if(b & (1 << i)) ret = 1LL * ret * a % m;
}
return ret;
}
int modinv(int x, int m) {
return modpow(x, m - 2, m);
}
int n; const int mod = 1000000007;
int main() {
cin >> n;
int ret = 1;
if(n < mod / 2) {
for(int i = 1; i <= n; i++) ret = 1LL * ret * i % mod;
}
else {
ret = mod - 1;
int mul = 1;
for(int i = mod - 1; i > n; i--) mul = 1LL * mul * i % mod;
ret = 1LL * ret * modinv(mul, mod) % mod;
}
cout << ret << endl;
return 0;
}

But it takes about 20 seconds in the worst case. Can you improve more?

v6ph1
Posts: 120
Joined: Mon Aug 25, 2014 7:14 pm

### Re: Do you have any efficient algorithm?

With a list of primes, it might be improved a little bit more:
e.g. 20! = 2^(10+5+2+1) * 3^(6+2) * 5^4 * 7^2 * 11 * 13 * 17 *19

But this does not reduce the complexitivity:
for each prime (n/log(n) times):
calc the power (log(n))
powermod (log(log(n))

But you could reduce the number divisions.

Btw. is there a formula to calc ((n-1)/2)! mod n?
((n-1)/2)! = (n-1)!/((n-1)/2)! * (-1)^((n-1)/2) (trivial to proof)

(((n-1)/2)!)^2 = (-1)^((n+1)/2)

But the squareroot function has 2 solutions

Min_25
Posts: 2
Joined: Sun Sep 28, 2014 8:59 am
Location: Japan

### Re: Do you have any efficient algorithm?

($n$! modulo $p$) can be computed in $O(p^{1/2} \log^{2+}{p})$ time, but it requires a polynomial library and it should be well-optimized to beat a naive algorithm.

Here is an algorithm: http://fredrikj.net/blog/2012/03/factor ... s-theorem/.

My program can compute 500,000,003 mod 1,000,000,007 in 0.3 seconds (0.87 seconds on ideone), but its source code has hundreds of lines .

As mentioned by v6ph1, we can also use the equality
$$n! = \prod _{p_i} p^{v(n, p_i)} ,$$
where $v(n, p) = \sum_{e=1}^{\infty}\left\lfloor\frac{n}{p^e}\right\rfloor$.
This solution runs in $O(n \log\log{n})$ time, but a cache-friendly segmented sieve of Eratosthenes can sieve primes less than $5 \times 10^8$ within 0.5 seconds (you can challenge the task in http://www.spoj.com/problems/PRIMES2/ ).

In addition, we need not to compute the power of primes $p > \sqrt{N}$ (for example, $11^5 \cdot 13^4 \cdot 17^3 \cdot 19^2 \cdot 23$ can be computed as $11 \cdot (11 \cdot 13) \cdot (11 \cdot 13 \cdot 17) \cdot (11 \cdot 13 \cdot 17 \cdot 19) \cdot (11 \cdot 13 \cdot 17 \cdot 19 \cdot 23$) if we keep the partial product). Probably this algorithm is easy to code and can compute $n! \pmod {10^9+7}$ within 1 seconds (3 seconds on ideone).

square1001
Posts: 27
Joined: Tue Mar 15, 2016 2:58 am
Location: Tokyo, Japan
Contact:

### Re: Do you have any efficient algorithm?

Thank you for reply, Min_25!

It's able to calculate factorials and binomial coefficients in hundreds of milliseconds with this algorithm.

Min_25
Posts: 2
Joined: Sun Sep 28, 2014 8:59 am
Location: Japan

### Re: Do you have any efficient algorithm?

I have added a simpler source code for calculating ($n!$ mod $p$).

The following code can compute $500,000,003 \pmod{10^9 + 7}$ in 0.91 sec on my (64-bit) computer or in 5 second on ideone, which uses an old 32-bit computer.

(The template parameter is used for the optimization on 64-bit computer, which converts the 32-bit constant remainder operation to two 64-bit x 64-bit multiplications.)

Code: Select all

#include <cstdio>
#include <vector>
#include <algorithm>
#include <cassert>

using namespace std;

using ll = long long;

int pow_mod(int b, int e, int mod) {
int ret = 1;
for (; e; e >>= 1, b = ll(b) * b % mod) {
if (e & 1) ret = ll(ret) * b % mod;
}
return ret;
}

int fact_valuation(int n, int p) {
int ret = n /= p;
while (n >= p) {
n /= p;
ret += n;
}
return ret;
}

template <int p>
int fact_mod(int n) {
if (n >= p) return 0;
if (n * 2 > p) {
int ret = fact_mod<p>(p - n - 1);
ret = (p - n - 1) & 1 ? ret : p - ret;
return pow_mod(ret, p - 2, p);
}
vector<int> ns;
for (int i = n; i > 1; i /= 5)
for (int j = i; j > 1; j /= 3)
for (int k = j; k > 1; k >>= 1)
ns.push_back(k);

std::sort(ns.begin(), ns.end());
const int d[8] = {6, 4, 2, 4, 2, 4, 6, 2};
int ret = 1, prod = 1, i = 1, r = 0;
for (auto end : ns) {
for (; i <= end; i += d[r & 7], ++r) {
prod = ll(prod) * i % p;
}
ret = ll(ret) * prod % p;
}
ret = ll(ret) * pow_mod(2, fact_valuation(n, 2), p) % p;
ret = ll(ret) * pow_mod(3, fact_valuation(n, 3), p) % p;
ret = ll(ret) * pow_mod(5, fact_valuation(n, 5), p) % p;
return ret;
}

int main() {
const int mod = 1e9 + 7;
printf("%d\n", fact_mod<mod>((mod - 1) / 2));
}