Modulo C (n, m) % (10 ^ 10)

Source: Internet
Author: User

Test Data

1515151 1213
.... 0836060000
151144 2002
.... 3558733440
10000000000 11411
.... 0000000000
115123131 1210
.... 2126393000
54546161515130 121231321
.... 6496000000


Because I am too watery in number theory, I spent more than a day on this topic, and finally found it a floating point Precision problem ......

Evaluate C (n, m) % (10 ^ 10), n <= 10 ^ 18 0 <= m <= n

There must be some number theory knowledge in such a large scope ~

Note that 10 ^ 10 can be split into 2 ^ 10*5 ^ 10

Can we modulo and merge the two items separately !!!!

So the current task is to find

C (n, m) % 2 ^ 10;

And C (n, m) % 5 ^ 10;

We can see that C (n, m) % (p ^ a) p is a prime number.

After reading the ultimate version of the AC blog, you will definitely do this.

Example C (n, m) % 1024

That is

 


Then convert it to n! % 1024

The denominator of the numerator is obtained separately, and then the inverse is obtained.

1*2*3 * .. * 1024

1025 *... * 2048

...

It is obvious that the result is the same for a period of 1024 (cyclical !!!), There are still a few remaining parts, and the brute force is enough.

Note: if we want to obtain the inverse element of the denominator, the denominator must have mutual quality with 2. Therefore, we first remove the multiples of 2. Obviously, the periodicity is still satisfied.

For example, 11 is required now! % (2 ^ 10)

Assume that calc (n, mod) can calculate n! Mod remainder of all odd-number products in

1*2*3*4*5*6*7*8*9*10*11

1*3*5*7*9*11*2*4*6*8*10 = 1*3*5*7*9*11*2 ^ 5*(1*2*3*4*5 ), the problem can be solved recursively !!!!

Therefore, you can write the formula calc (11/2) * 2 ^ (11/4) * calc (5) * (2 ^ 11/8) * calc (2) * 2 ^) * calc (1)

Complexity is log-level

This question is special. There are only two qualitative factors. When P is a general number, it is better to separate it and combine it with the Chinese Remainder Theorem.

However, if p ^ a is large after prime factor decomposition, the above method cannot be used because the array cannot be opened and the cycle cannot be used ....

So when n m k is a 10 ^ 9 value, it seems that it cannot be done ????

I think K is 10 ^ 9, and m must be 10 ^ 6.

If K is 10 ^ 5, n m can be used at will.

Determine whether C (n, m) is greater than 10 ^ 10. If double is used in the middle, WA 9 is used at the beginning. If it is not checked out, use eps to determine whether it will pass.

Another way is to refer to the code area as follows:

Paste code


[Cpp]
# Include <cstdio>
# Include <cstring>
# Include <cmath>
Typedef long lld;
Const lld mod5 = 9765625;
Lld n, m;
Lld fac [10000000];
Lld two [64];
Lld five [64];
Lld Pow (lld a, lld B, lld mod)
{
Lld ans = 1;
While (B ){
If (B & 1) ans = ans * a % mod;
A = a * a % mod; B >>= 1;
}
Return ans;
}
Lld exgcd (lld a, lld B, lld & x, lld & y)
{
If (B = 0) {x = 1; y = 0; return ;}
Else {
Lld d = exgcd (B, a % B, x, y );
Lld t = x;
X = y;
Y = t-a/B * y;
Return d;
}
}
Lld inverse (lld num, lld mod)
{
Lld x, y;
Exgcd (num, mod, x, y );
While (x <0) x + = mod, y-= num;
Return x;
}
Void init1 ()
{
Two [0] = 1;
For (lld I = 1; I <64; I ++) two [I] = two [I-1] * 2;
Fac [1] = 1; fac [0] = 1;
For (lld I = 2; I <= 1024; I ++ ){
Lld num = I;
Fac [I] = fac [I-1];
If (I & 1) fac [I] = fac [I] * num % 1024;
}
}
Lld calc (lld n, lld mod) // n! Remove the multiples of 2 or 5 and then remainder the mod.
{
Lld t = n/mod;
Lld ans = 1;
If (t> = 1) ans = Pow (fac [mod], t, mod );
Int lim = n % mod;
Return ans * fac [lim] % mod;
}
Lld gao (lld B [], lld mod, lld num) // calculate c (n, m) % mod
{
Long fenzi _ = 0, fenmu _ = 0;
Lld tn = n, tm = m, = n-m;
For (lld I = 1; B [I] <= tn; fenzi _ + = (tn/B [I]), I ++ );
For (lld I = 1; B [I] <= tm; fenmu _ + = (tm/B [I]), I ++ );
For (lld I = 1; B [I] <=; fenmu _ + = (/B [I]), I ++ );
Lld cnt = fenzi _-fenmu _;
Lld fenzi = calc (tn, mod );
For (lld I = 1; B [I] <= tn; I ++) fenzi = fenzi * calc (tn/B [I], mod) % mod;
Lld fenmu = calc (OLAP, mod );
For (lld I = 1; B [I] <=; I ++) fenmu = fenmu * calc (/B [I], mod) % mod;
Fenmu = fenmu * calc (tm, mod) % mod;
For (lld I = 1; B [I] <= tm; I ++) fenmu = fenmu * calc (tm/B [I], mod) % mod;
Return fenzi * inverse (fenmu, mod) % mod * Pow (num, cnt, mod) % mod;
}
Void init2 ()
{
Five [0] = 1;
For (lld I = 1; I <= 50; I ++) five [I] = five [I-1] * 5;
Fac [1] = 1; fac [0] = 1;
For (lld I = 2; I <= mod5; I ++ ){
Lld num = I;
Fac [I] = fac [I-1];
If (I % 5! = 0) fac [I] = fac [I] * num % mod5;
}
}
// X % 1024 =;
// X % mod5 = B;
// Combine the two formulas with the Chinese Remainder Theorem
Lld CRT (lld a, lld B)
{
Lld mod = 0000000000ll;
Lld x, y;
Exgcd (mod5, 1024, x, y );
X = x * a % mod;
Lld p = x, q;
Exgcd (1024, mod5, x, y );
X = x * B % mod;
Q = x;
Lld ans = (mod5 * p % mod + 1024LL * q % mod) % mod;
Return (ans + mod) % mod;
}
Void print (lld ans)
{
If (n-m <m) m = n-m;
Double t = 1;
Lld a = 1;
Bool flag = false;
For (lld I = 1; I <= m; I ++)
{
T = t * (n-I + 1)/I;
A = a * (n-I + 1)/I;
If (a> = 0000000000ll | t> = 1e10)
{
Flag = true;
Break;
}
}
If (flag) printf ("... % 010I64d \ n", ans );
Else printf ("% I64d \ n", ans );
}
Int main ()
{
Freopen ("combi. in", "r", stdin );
Freopen ("combi. out", "w", stdout );
While (scanf ("% I64d % I64d", & n, & m )! = EOF ){
Init1 ();
Lld a = gao (two, 2 );
Init2 ();
Lld B = gao (five, mod5, 5 );
Lld ans = CRT (a, B );
Print (ans );
}
Return 0;
}
/*
1515151 1213
.... 0836060000
151144 2002
.... 3558733440
10000000000 11411
.... 0000000000
115123131 1210
.... 2126393000
54546161515130 121231321
.... 6496000000
*/

# Include <cstdio>
# Include <cstring>
# Include <cmath>
Typedef long lld;
Const lld mod5 = 9765625;
Lld n, m;
Lld fac [10000000];
Lld two [64];
Lld five [64];
Lld Pow (lld a, lld B, lld mod)
{
Lld ans = 1;
While (B ){
If (B & 1) ans = ans * a % mod;
A = a * a % mod; B >>= 1;
}
Return ans;
}
Lld exgcd (lld a, lld B, lld & x, lld & y)
{
If (B = 0) {x = 1; y = 0; return ;}
Else {
Lld d = exgcd (B, a % B, x, y );
Lld t = x;
X = y;
Y = t-a/B * y;
Return d;
}
}
Lld inverse (lld num, lld mod)
{
Lld x, y;
Exgcd (num, mod, x, y );
While (x <0) x + = mod, y-= num;
Return x;
}
Void init1 ()
{
Two [0] = 1;
For (lld I = 1; I <64; I ++) two [I] = two [I-1] * 2;
Fac [1] = 1; fac [0] = 1;
For (lld I = 2; I <= 1024; I ++ ){
Lld num = I;
Fac [I] = fac [I-1];
If (I & 1) fac [I] = fac [I] * num % 1024;
}
}
Lld calc (lld n, lld mod) // n! Remove the multiples of 2 or 5 and then remainder the mod.
{
Lld t = n/mod;
Lld ans = 1;
If (t> = 1) ans = Pow (fac [mod], t, mod );
Int lim = n % mod;
Return ans * fac [lim] % mod;
}
Lld gao (lld B [], lld mod, lld num) // calculate c (n, m) % mod
{
Long fenzi _ = 0, fenmu _ = 0;
Lld tn = n, tm = m, = n-m;
For (lld I = 1; B [I] <= tn; fenzi _ + = (tn/B [I]), I ++ );
For (lld I = 1; B [I] <= tm; fenmu _ + = (tm/B [I]), I ++ );
For (lld I = 1; B [I] <=; fenmu _ + = (/B [I]), I ++ );
Lld cnt = fenzi _-fenmu _;
Lld fenzi = calc (tn, mod );
For (lld I = 1; B [I] <= tn; I ++) fenzi = fenzi * calc (tn/B [I], mod) % mod;
Lld fenmu = calc (OLAP, mod );
For (lld I = 1; B [I] <=; I ++) fenmu = fenmu * calc (/B [I], mod) % mod;
Fenmu = fenmu * calc (tm, mod) % mod;
For (lld I = 1; B [I] <= tm; I ++) fenmu = fenmu * calc (tm/B [I], mod) % mod;
Return fenzi * inverse (fenmu, mod) % mod * Pow (num, cnt, mod) % mod;
}
Void init2 ()
{
Five [0] = 1;
For (lld I = 1; I <= 50; I ++) five [I] = five [I-1] * 5;
Fac [1] = 1; fac [0] = 1;
For (lld I = 2; I <= mod5; I ++ ){
Lld num = I;
Fac [I] = fac [I-1];
If (I % 5! = 0) fac [I] = fac [I] * num % mod5;
}
}
// X % 1024 =;
// X % mod5 = B;
// Combine the two formulas with the Chinese Remainder Theorem
Lld CRT (lld a, lld B)
{
Lld mod = 0000000000ll;
Lld x, y;
Exgcd (mod5, 1024, x, y );
X = x * a % mod;
Lld p = x, q;
Exgcd (1024, mod5, x, y );
X = x * B % mod;
Q = x;
Lld ans = (mod5 * p % mod + 1024LL * q % mod) % mod;
Return (ans + mod) % mod;
}
Void print (lld ans)
{
If (n-m <m) m = n-m;
Double t = 1;
Lld a = 1;
Bool flag = false;
For (lld I = 1; I <= m; I ++)
{
T = t * (n-I + 1)/I;
A = a * (n-I + 1)/I;
If (a> = 0000000000ll | t> = 1e10)
{
Flag = true;
Break;
}
}
If (flag) printf ("... % 010I64d \ n", ans );
Else printf ("% I64d \ n", ans );
}
Int main ()
{
Freopen ("combi. in", "r", stdin );
Freopen ("combi. out", "w", stdout );
While (scanf ("% I64d % I64d", & n, & m )! = EOF ){
Init1 ();
Lld a = gao (two, 2 );
Init2 ();
Lld B = gao (five, mod5, 5 );
Lld ans = CRT (a, B );
Print (ans );
}
Return 0;
}
/*
1515151 1213
.... 0836060000
151144 2002
.... 3558733440
10000000000 11411
.... 0000000000
115123131 1210
.... 2126393000
54546161515130 121231321
.... 6496000000
*/

 

Related Article

Contact Us

The content source of this page is from Internet, which doesn't represent Alibaba Cloud's opinion; products and services mentioned on that page don't have any relationship with Alibaba Cloud. If the content of the page makes you feel confusing, please write us an email, we will handle the problem within 5 days after receiving your email.

If you find any instances of plagiarism from the community, please send an email to: info-contact@alibabacloud.com and provide relevant evidence. A staff member will contact you within 5 working days.

A Free Trial That Lets You Build Big!

Start building with 50+ products and up to 12 months usage for Elastic Compute Service

  • Sales Support

    1 on 1 presale consultation

  • After-Sales Support

    24/7 Technical Support 6 Free Tickets per Quarter Faster Response

  • Alibaba Cloud offers highly flexible support services tailored to meet your exact needs.