Modular operation
In programming competitions, there are quite a few problems with too large data results, which often need to be modeled. For example, (a*b)% P, if the result of a*b can't be stored, then remove the module, the result is obviously wrong. In order to prevent overflow, we can take the module of a, B and then quadrature respectively.
Modular operation formula:
Add: (a +% B)% P = (a% P + B% P)% p
Subtraction: (a - b)% P = ((a% p - B% P) + P)% p
Multiplication: (a * b)% P = (a% P) * (B% P)% p
Power operation: A ^ B% P = ((a% P) ^ b)% p
Division: (A / b)% P = (a% p * inv% P)% p
Note: (A / b)% P cannot be directly / b. you need to find a number inv so that inv * B% P = 1 We call inv the inverse element. Now we will introduce the method of finding the inverse element.
In the competition, we usually give the result% 1e9+7, so method 1 is recommended. Of course, if a happens to be a multiple of b, it can be calculated directly.
/*==================================================*\
| Application scope of method 1 / method 2: mod is prime number
| According to Fermat's small theorem, (a/b)%mod = a*pow_mod(b,mod-2)%mod
| pow_mod For fast power functions
\*==================================================*/
ll getInv(ll b) {
return pow_mod(b, mod-2);
}
ll getInv(ll b) {
if(b==1) return 1;
return (mod-mod/b) * getInv(mod%b)%mod;
}
/*==================================================*\
| Method 3: the inverse element of each number can be obtained by solving the inverse element linearly
\*==================================================*/
const int N = 100005;
long long inv[N];
void inv_init() {
inv[0] = inv[1] = 1;
for (int i = 2; i < N; i++) {
inv[i] = ((mod - mod/ i) * inv[mod % i]) % mod;
}
}
/*==================================================*\
| Find the inverse element of b under mod, no inverse element returns - 1
| exgcd To extend Euclid algorithm
| Scope of application: as long as there is inverse element, it can be solved. It is applicable to the case of small number but large mod
\*==================================================*/
ll getInv(int b) {
ll x, y;
ll d = exgcd(b,mod,x,y);
return d==1 ? (x+mod)%mod : -1;
}Modular multiplication of large numbers
/*==================================================*\
| Find a * B% Mod
| Algorithm understanding: mul_mod (100,10) = > 10 binary 1010
| 100 * 10 == 100*0 + 200*1 + 400*0 + 800*1
\*==================================================*/
ll mul_mod(int a,int b) {
ll res = 0;
while( b ) {
if(b & 1) res = (res + a) % mod;
a = 2 * a % mod;
b >>= 1;
}
return res;
}Fast power modulus
/*==================================================*\
| Algorithm immediate: 11 binary 1011 = = > 11 = 2 ^ 0 + 2 ^ 1 + 2 ^ 3
| 2^11 = 2^(2^0) * 2^(2^1) * 2^(2^3)
| If the multiplication overflows, use the mul mod function
\*==================================================*/
ll pow_mod(ll x,ll n) {
ll res=1;
while(n > 0) {
if(n & 1) res = res*x % mod;
x = x*x % mod;
n >>= 1;
}
return res;
}Euler power down formula
/*==================================================*\
| Find A^n mod m
| A^n ≡ A^(n%ϕ(m)+ϕ(m)) % m , n>ϕ(m)
| Pow ﹣ mod is a fast power function, and euler is a euler function, which will be discussed later
\*==================================================*/
ll euler_pow(ll a, char* n) {
ll phi, _n = 0;
phi = euler(mod);
int len = strlen(n);
for(int i=0; i < len; i++)
_n = (_n*10+(n[i]-'0')) % phi; // Get n% ϕ (m)
_n += phi; // Get n% ϕ (m) + ϕ (m)
return pow_mod(a, _n);
}Square root of large number
/*==================================================*\
| Square root of large number (represented by string array)
\*==================================================*/
void Sqrt(char *str) {
double i, r, n;
int j, l, size, num, x[1000];
size = strlen(str);
if( size == 1 && str[0] == '0' ) {
printf("0\n"); return;
}
if( size%2 == 1 ) {
n = str[0]-48; l = -1;
} else {
n = (str[0]-48)*10+str[1]-48;
l = 0;
}
r = 0; num = 0;
while(true) {
i = 0;
while( i*(i+20*r) <= n ) ++i; --i;
n -= i*(i+20*r);
r = r*10+i;
x[num] = (int)i;
++num;
l += 2;
if( l >= size ) break;
n = n*100+(double)(str[l]-48)*10+(double)(str[l+1]-48);
}
for(j = 0;j < num; j++)
printf("%d", x[j]);
printf("\n");
} To find the number of factorials
#include <bits/stdc++.h>
using namespace std;
const double PI = acos(-1.0);
const double e = exp(1.0);
int main() {
int n ;
cin >> n;
//Stirling's theorem for the number of N
double res = log10(sqrt(2*PI*n)) + log10(n/e)*n + 1;
cout << res << endl;
// //log for N
double ans = 1;
for(int i=1;i <= n; i++)
ans += log10((double)i);
cout << ans << endl;
return 0;
}The last nonzero position of factorial
/*==================================================*\
| Factorial last nonzero, complexity O(nlogn)
\*==================================================*/
//Return this bit, n passed in as a string
#include <string.h>
#define MAXN 10000
const int mod[20]={1,1,2,6,4,2,2,4,2,8,4, 4,8,4,6,8,8,6,8,2};
int lastdigit(char* buf) {
int len=strlen(buf), a[MAXN], i, c, ret=1;
if (len==1) return mod[buf[0]-'0'];
for (i=0;i<len;i++) a[i]=buf[len-1-i]-'0';
for (; len; len-=!a[len-1]){
ret = ret * mod[a[1] % 2*10+a[0]] % 5;
for (c=0,i=len-1; i>=0; i--)
c = c*10+a[i], a[i]=c/5, c%=5;
}
return ret+ret%2*5;
} Euler function
/*---------------------Euler's theorem------------------------
Euler function: number of elements in 1~n and n-coprime φ (n)
Euler Theorem: if gcd(a, n)=1, then a ^ φ (n)=1 (mod n)
When B is large, a^b = a^(b mod φ (n))(mod n)
When n is an odd number, φ (2n) = φ (n)
For prime p, φ (p) = p-1, φ (p^n) = p^n-p^(n-1)
When (m,n)=1, φ (mn) = φ (m) × φ (n).
-----------------------------------------------------*/
// Calculate separately the Euler function phi (n), 12 = 2 * 2 * 3
unsigned euler(unsigned n) {
unsigned res = n;
for (int i = 2; i*i <= n; i++) {
if(n % i == 0) {
res -= res / i; // res = res / i * (i - 1);
while (n % i == 0) n /= i; // Remove the same factors
}
}
if (n > 1) // If it is not 1, it means that the current n is a quality factor of res
res -= res / n;
return res;
}
// Recurrence of Euler function phi(i)
void initEuler() {
for (i = 1; i <= maxn; i++) phi[i] = i; // Initialize to itself
for (i = 2; i <= maxn; i += 2) phi[i] /= 2;
for (i = 3; i <= maxn; i += 2) if(phi[i] == i) {
for (j = i; j <= maxn; j += i)
phi[j] = phi[j] - phi[j]/i; // phi[j] = phi[j] / i * (i-1)
}
}Finding prime number by screening method
const int N= 1e6; // Prime table range
bool isp[N+5];
int len = 0, p[N+5];
void prime() { //Recommend this, faster
isp[0] = isp[1] = true;
for (int i = 2; i < maxn; i++) {
if(!isp[i]) p[++len] = i;
for (int j = 1; j <= len && p[j]*i < maxn; j++) {
isp[i*p[j]] = true;
if (i%p[j] == 0) break;
}
}
}
void prime() {
for(int i=2; i <= N; i++)
isp[i] = true;
int n = sqrt(N);
for(int i=2; i <= n; i++) if(isp[i]) {
for(int j=2; i*j <= N; j++)
isp[i*j] = false;
}
}primality testing
/*==================================================*\
| Use when n is a large number
| pow_mod Modular multiplication of large numbers
\*==================================================*/
const int S=20; // The greater the number of decision, the smaller the probability of error
bool Miller_Rabin(ll n) {
if (n == 2)
return true;
if (n < 2 || !(n & 1))
return false;
ll m = n - 1, k = 0;
while (!(m & 1)) {
k++;
m >>= 1;
}
for (int i = 1; i <= S; i++) {
ll a = rand() % (n - 1) + 1;
ll x = pow_mod(a, m, n);
ll y;
for (int j = 1; j <= k; j++) {
y = mul_mod(x, x, n);
if (y == 1 && x != 1 && x != n - 1)
return false;
x = y;
}
if (y != 1)
return false;
}
return true;
}Maximum common divisor (GCD)
/*==================================================*\
| GCD greatest common divisor
\*==================================================*/
int gcd(int x, int y){
if (!x || !y) return x > y ? x : y;
for (int t; t = x % y; x = y, y = t);
return y;
}
/*==================================================*\
| Fast GCD
\*==================================================*/
int kgcd(int a, int b){
if (a == 0) return b;
if (b == 0) return a;
if (!(a & 1) && !(b & 1)) return kgcd(a>>1, b>>1) << 1;
else if (!(b & 1)) return kgcd(a, b>>1);
else if (!(a & 1)) return kgcd(a>>1, b);
else return kgcd(abs(a - b), min(a, b));
}
/*==================================================*\
| Extended GCD algorithm
| Find x, y to make gcd(a, b) = a * x + b * y = d;
| Application: solution of indefinite equation, solution of linear congruence equation, inverse element of seeking touch
\*==================================================*/
int exgcd(int a, int b, int & x, int & y){
if (b == 0) { x=1; y=0; return a; }
int d = exgcd(b, a % b, x, y);
int t = x; x = y; y = t - a / b * y;
return d;
} Expanding China's surplus theorem
/*==================================================*\
| x ≡ a1 (mod m1),x ≡ a2 (mod m2),x ≡ a3 (mod m3)Seeking x.
| exgcd Extended Euclid, Mul Mod
\*==================================================*/
ll excrt() {
ll x, y;
ll M = mn[1], ans = an[1];
for(int i=2; i <= N; i++) {
ll a=M, b=mn[i], c=(an[i]-ans%b+b)%b;
ll d = exgcd(a, b, x, y);
if(c%d != 0) return -1;
x = mul_mod(x, c/d, b/d);
ans += x*M;
M *= b/d;
ans = (ans%M+M) % M;
}
return (ans%M+M) % M;
}
combination
/*==================================================*\
| Combination number C(n, r)
\*==================================================*/
int com(int n, int r){// return C(n, r)
if( n-r > r ) r = n-r; // C(n, r) = C(n, n-r)
int s = 1;
for(int i=0, j=1; i < r; ++i ) {
s *= (n-i);
for( ; j <= r && s%j == 0; ++j ) s /= j;
}
return s;
}
/*==================================================*\
| Combination number C(n, r)
\*==================================================*/
void init() {
ll i;
fac[0]=1;
for (ll i = 1; i < N; i++)
fac[i] = fac[i - 1] * i % MOD;
}
ll inv(ll a, ll n) {
ll x, y;
exgcd(a, n, x, y);
return (x + n) % n;
}
ll C(ll n, ll r) {
return fac[n] * inv(fac[r] * fac[n - r] % mod, mod) % mod;
}
/*==================================================*\
| O(n) Find the combination number and recurs from left to right
\*==================================================*/
C[0] = 1;
for(int i = 1; i <= n; i++)
C[i] = C[i - 1] * (n - i + 1) / i;
polya count
/*==================================================*\
| c Beads of different colors make up a necklace of s. the necklace has no direction and starting position;
\*==================================================*/
int gcd (int a, int b) { return b ? gcd(b,a%b) : a; }
int main (void){
int c, s;
while (scanf("%d%d", &c, &s)==2) {
int k;
long long p[64]; p[0] = 1; // power of c
for (k=0 ; k<s ; k++) p[k+1] = p[k] * c;
// reflection part
long long count = s&1 ? s * p[s/2 + 1] : (s/2) * (p[s/2] + p[s/2 + 1]);
// rotation part
for (k=1 ; k<=s ; k++) count += p[gcd(k, s)];
count /= 2 * s;
cout << count << '\n';
}
return 0;
}Staggered formula
/*------------------------------------------------------------
Question: ten different books are on the shelf.
Now it is rearranged so that each book is not in the original place. How many pendulums are there?
To generalize this problem, it is the problem of permutation, which is one of the problems in combinatorics.
Consider an arrangement with n elements. If all elements in an arrangement are not in their original positions,
Then such an arrangement is called a wrong arrangement of the original arrangement The number of misplaced elements is recorded as D(n).
The study of the number of permutations is called the permutation problem or the replacement problem.
------------------------------------------------------------*/
//dp[i] = (i - 1)*(dp[i - 1] + dp[i - 2]); i > 2 dp[1]=0 dp[2]=1
ll D(int n) {
ll a = 0, b = 1, c;
if(n < 3) return n-1;
for (int i = 3; i <= n; i++) {
c = ((i - 1) * 1ll * (a + b)) % mod;
a = b;
b = c;
}
return c;
}Common generating function
/*--------------------------------------------------
Definition: G(x) = a0 · x^0 + a1 · x^1 +...+ an · x^n
Then G(x) is the generating function of a0, a1..., an
Classic example: there are 1, 3 and 2 weights with a mass of 1, 2 and 4
Q: how many different qualities can be weighed out? How many ways to call quality 3?
Solution: constructor G(x)=(1+x)*(1+x^2+x^4+x^6)*(1+x^4+x^8)
In the second bracket, X^6 means: use the second item with 3 pieces of quality, and so on
The code simulates the calculation process of the constructor, such as:
(1+x) * (1+x^2+x^4+x^6) * (1+x^4+x^8)
= (1+x^2+x^4+x^6 + x+x^3+x^5+x^7) * (1+x^4+x^8)
= ...
--------------------------------------------------*/
// v[i] the ith unknown value, s/e: the maximum / minimum number of items
int v[MAX], s[MAX], e[MAX];
// a is the final result and b is the intermediate result.
int a[MAX], b[MAX];
// P is the maximum possible index
int P;
void mu(int n) { // n is the number of species
memset(a,0,sizeof(a));
a[0]=1;
for(int i=1; i <= n; i++) { //Cycle each factor
memset(b, 0, sizeof(b));
// J is the possible quantity of the ith item, j*v[i] is the coefficient of the j item in the I bracket
for (int j=n1[i]; j<=n2[i] && j*v[i]<=P; j++) {
// Calculate the k-th coefficient of the result of multiplying the former i-1
for (int k=0; k+j*v[i] <= P; k++)
b[k+j*v[i]] += a[k];
}
memcpy(a, b, sizeof(b)); // Save b array back to a array
}
}Exponential generating function
/*--------------------------------------------------
Definition: G (x) = a0 + (a1/1!) · x^1 + (a2/2!) · x^2 +
The exponential generating function corresponding to the series a0, a1, a2,..., an
Used to find the number of multiple permutations, as follows:
Classic example: suppose there are multiple elements, in which a1 repeats n1 times, a2 repeats n2 times,
a3 Repeat n3 times. Take r permutations from them and find their permutations.
Constructor G(x) = (1+x/1!+x^2/2!+...+x^n1/n1!)*
(1+x/1!+...+x^n2/n2!) *(1+x/1!+...+x^n3/n13)
It's not hard to understand the meaning of each of them according to the common generating function
--------------------------------------------------*/
ll fac[N]; // fac[i] represents the ith factorial value
int a[N]; //1~n number of each type
double c1[N],c2[N]; // double
void cal(int n,int m) { //There are n kinds, m
memset(c1, 0, sizeof(c1));
memset(c2, 0, sizeof(c2));
c1[0] = 1.0/fac[0];
for(int i=1; i<=n; i++) {
for(int j=0; j<=m; j++) {
for(int k=0; k+j<=m && k<=a[i]; k++)
c2[k+j] += c1[j]/fac[k];
}
for(int j=0; j<=m; j++) {
c1[j] = c2[j];
c2[j]=0;
}
}
}
ans=c1[m]*fac[m]; //The number of multiple permutations when m is takenMobius function
// Recursive screen method n*logn
void getMu(){
for(int i=1; i<=MAXN; i++) {
int target = i==1?1:0;
int delta = target - mu[i];
mu[i]=delta;
for(int j=i*2; j<=MAXN; j+=i)
mu[j]+=delta;
}
}
//Seeking Mobius function by linear sieve method
bool check[MAXN+10];
int prime[MAXN+10];
int mu[MAXN+10];
void Moblus() {
memset(check,false,sizeof(check));
mu[1] = 1;
int tot = 0;
for(int i = 2; i <= MAXN; i++) {
if( !check[i] ) {
prime[tot++] = i;
mu[i] = -1;
}
for(int j = 0; j < tot; j++) {
if(i * prime[j] > MAXN) break;
check[i * prime[j]] = true;
if( i % prime[j] == 0){
mu[i * prime[j]] = 0;
break;
}else{
mu[i * prime[j]] = -mu[i];
}
}
}
}
// Find sum(d) = ∑ p|d μ (dp)
const int maxn = 10000010;
int prime[maxn],mu[maxn],sum[maxn];
bool check[maxn];
void Mobius(){
memset(check,false,sizeof(check));
mu[1] = 1;
prime[0] = 0;
for(int i=2;i<maxn;i++){
if(!check[i]){
mu[i] = -1;
sum[i] = 1;
prime[++prime[0]] = i;
}
for(int j=1;j<=prime[0];j++){
if(i*prime[j] >= maxn) break;
check[i*prime[j]] = true;
if(i % prime[j]){
mu[i*prime[j]] = -mu[i];
sum[i*prime[j]] = mu[i] - sum[i];
}
else{
mu[i*prime[j]] = 0;
sum[i*prime[j]] = mu[i];
break;
}
}
}
}
Inverse number
#include <bits/stdc++.h>
using namespace std;
const int maxn=500005;
int n;
int aa[maxn]; //Array after discretization
int c[maxn]; //Tree array
struct Node {
int v;
int order;
}in[maxn];
int lowbit(int x) {
return x&(-x);
}
void update(int t,int value) {
int i;
for(i=t;i<=n;i+=lowbit(i)) {
c[i]+=value;
}
}
int getsum(int x) {
int i;
int temp=0;
for(i=x;i>=1;i-=lowbit(i)) {
temp+=c[i];
}
return temp;
}
bool cmp(Node a ,Node b) {
return a.v<b.v;
}
int main() {
int i,j;
while(scanf("%d",&n)==1 && n) {
//Discretization
for(i=1;i<=n;i++) {
scanf("%d",&in[i].v);
in[i].order=i;
}
sort(in+1,in+n+1,cmp);
for(i=1;i<=n;i++) aa[in[i].order]=i;
//Reverse order of tree array
memset(c,0,sizeof(c));
long long ans=0;
for(i=1;i<=n;i++) {
update(aa[i],1);
ans+=i-getsum(aa[i]);
}
cout<<ans<<endl;
}
return 0;
}
discrete logarithm
ll exBSGS(ll a, ll b, ll p) {
a %= p; b %= p;
ll ret = 1;
for(ll i = 0; i <= 50; ++i) {
if(ret == b) return i;
ret = (ret*a) % p;
} //Enumerate smaller i
ll x,y,d, v = 1, cnt = 0;
while((d = gcd(a, p)) != 1) {
if(b % d) return -1;
b /= d, p /= d;
v = (v * (a/d)) % p;
++cnt;
}//Approximately until (a, p) == 1
map<ll, ll> h;
ll m = ceil(sqrt(p)), t = 1;
for(ll i = 0; i < m; ++i) {
if(h.count(t)) h[t] = min(h[t], i);
else h[t] = i;
t = (t*a) % p;
}
for(ll i = 0; i < m; ++i) {
d = extgcd(v, p, x, y);
x = (x* (b/d) % p + p) % p;
if(h.count(x)) return i*m + h[x] + cnt;
v = (v*t) % p;
}
return -1;
}
Primitive root
#include <bits/stdc++.h>
using namespace std;
const int N = 1000010;
bitset<N> prime;
int p[N],pri[N];
int k,cnt;
void isprime() {
prime.set();
for(int i=2; i<N; i++) {
if(prime[i]) {
p[k++] = i;
for(int j=i+i; j<N; j+=i)
prime[j] = false;
}
}
}
void Divide(int n) {
cnt = 0;
int t = (int)sqrt(1.0*n);
for(int i=0; p[i]<=t; i++) {
if(n%p[i]==0) {
pri[cnt++] = p[i];
while(n%p[i]==0) n /= p[i];
}
}
if(n > 1)
pri[cnt++] = n;
}
int main() {
int P;
isprime();
while(cin>>P) {
Divide(P-1);
for(int g=2; g<P; g++) {
bool flag = true;
for(int i=0; i<cnt; i++) {
int t = (P - 1) / pri[i];
if(pow_mod(g,t,P) == 1) {
flag = false;
break;
}
}
if(flag) {
int root = g;
cout<<root<<endl;
}
}
}
return 0;
}
Note: most codes come from the network