#include using namespace std; using i64 = int64_t; using u64 = uint64_t; using base = complex; void fft(vector &a, bool inv){ int n = a.size(), j = 0; for(int i=1; i> 1); while(j >= bit){ j -= bit; bit >>= 1; } j += bit; if(i < j) swap(a[i], a[j]); } double ang = 2 * acos(-1) / n * (inv ? -1 : 1); vector roots(n/2); for(int i=0; i multiply(vector &v, vector &w, i64 mod){ int n = 2; while(n < v.size() + w.size()) n <<= 1; vector v1(n), v2(n), r1(n), r2(n); for(int i=0; i> 15, v[i] & 32767); for(int i=0; i> 15, w[i] & 32767); fft(v1, 0); fft(v2, 0); for(int i=0; i ret(n); for(int i=0; i lagrange(vector h, i64 P) { int d = (int)h.size()-1; assert(d >= 0); auto mul = [P](i64 a, i64 b){ return a*b%P; }; auto ipow = [&mul](i64 a, i64 b) { i64 res = 1; while(b) { if(b&1) res = mul(res, a); a = mul(a, a); b >>= 1; } return res; }; auto modInv = [&ipow, P](i64 a) { return ipow(a, P-2); }; vector fact(4*d+2), invfact(4*d+2); fact[0] = 1; for(int i=1; i<=4*d+1; ++i) fact[i] = mul(i,fact[i-1]); invfact[4*d+1] = modInv(fact[4*d+1]); for(int i=4*d; i>=0; --i) invfact[i] = mul(invfact[i+1], i+1); vector f(d+1); for(int i=0; i<=d; ++i) { f[i] = h[i]; f[i] = mul(invfact[i], f[i]); f[i] = mul(invfact[d-i], f[i]); if((d-i)%2==1) f[i] = P-f[i]; if(f[i] == P) f[i] = 0; } vector inv(4*d+2); for(int i=1; i<4*d+2; ++i) inv[i] = mul(fact[i-1], invfact[i]); vector g(4*d+2); g[d+1] = 1; for(int j=0; j<=d; ++j) g[d+1] = mul(g[d+1], d+1-j); for(int i=d+2; i<4*d+2; ++i) g[i] = mul(g[i-1], mul(i, inv[i-d-1])); vector conv = multiply(f, inv, P); vector ret(4*d+2); for(int i=0; i<=d; ++i) ret[i] = h[i]; for(int i=d+1; i<4*d+2; ++i) ret[i] = mul(g[i], conv[i]); return ret; } vector squarepoly(vector poly, i64 P) { vector ss = lagrange(poly, P); vector ret(ss.size()/2); auto mul = [P](i64 a, i64 b){ return a*b%P; }; for(int i=0; i<(int)ss.size()/2; ++i) ret[i] = mul(ss[2*i], ss[2*i+1]); return ret; } i64 factorial(i64 N,i64 P) { i64 d = 1; vector fact_part = {1%P, 2%P}; while(N > d*(d+1)) { fact_part = squarepoly(fact_part, P); d *= 2; } auto mul = [P](i64 a, i64 b) { return a*b%P; }; i64 ans = 1, bucket = N/d; for(int i=0; i> N; P = 1e9 + 7; if(N>=P) { cout << 0 << '\n'; return 0; } cout << factorial(N,P) << '\n'; return 0; }