数学
配合oiwiki:
位运算
int __builtin_ffs(int x)：返回 x 的二进制末尾最后一个 1 的位置，位置的编号从 1 开始（最低位编号为 1 ）。当 x 为 0 时返回 0 。int __builtin_clz(unsigned int x)：返回 x 的二进制的前导 0 的个数。当 x 为 0 时，结果未定义。int __builtin_ctz(unsigned int x)：返回 x 的二进制末尾连续 0 的个数。当 x 为 0 时，结果未定义。int __builtin_clrsb(int x)：当 x 的符号位为 0 时返回 x 的二进制的前导 0 的个数减一，否则返回 x 的二进制的前导 1 的个数减一。int __builtin_popcount(unsigned int x)：返回 x 的二进制中 1 的个数。int __builtin_parity(unsigned int x)：判断 x 的二进制中 1 的个数的奇偶性。
可以在函数名末尾添加ll如__builtin_popcountll以操作long long
如果需要操作的集合非常大，可以使用 bitset
bitset
可视为多位二进制数
n位bitset执行一次位运算的时间复杂度可视为n/32
同样支持 ~ | ^ & >> << == != 等操作符,可以用 [] 查询
s.count() 返回二进制串中有多少个1；
s.set()把s所有位变为1；
s.set(k,v)把s的第k位改为v,即s[k]=v；
s.reset()把s的所有位变为0.
s.reset(k)把s的第k位改为0,即s[k]=0；
s.flip()把s所有位取反.即s=~s；
s.flip(k)把s的第k位取反，即s[k]^=1；
高精度
https://oi-wiki.org/math/bignum/
const int LEN = 1004;
void clear(int a[]) {
for (int i = 0; i < LEN; ++i) a[i] = 0;
}
void read(int a[]) {
string s;
cin>>s;
clear(a);
int len = s.size();
for (int i = 0; i < len; ++i) a[len - i - 1] = s[i] - '0';
}
void print(int a[]) {
int i;
for (i = LEN - 1; i >= 1; --i)
if (a[i] != 0) break;//忽略前导0
for (; i >= 0; --i) putchar(a[i] + '0');
}
void add(int a[], int b[], int c[]) {
clear(c);
for (int i = 0; i < LEN - 1; ++i) {
c[i] += a[i] + b[i];
if (c[i] >= 10) {
c[i + 1] += 1;
c[i] -= 10;
}
}
}
void sub(int a[], int b[], int c[]) {
clear(c);
for (int i = 0; i < LEN - 1; ++i) {
c[i] += a[i] - b[i];
if (c[i] < 0) {
c[i + 1] -= 1;
c[i] += 10;
}
}
}
void mul(int a[], int b[], int c[]) {
clear(c);
for (int i = 0; i < LEN - 1; ++i) {
for (int j = 0; j <= i; ++j) c[i] += a[j] * b[i - j];
if (c[i] >= 10) {
c[i + 1] += c[i] / 10;
c[i] %= 10;
}
}
}
inline bool greater_eq(int a[], int b[], int last_dg, int len) {
if (a[last_dg + len] != 0) return true;
for (int i = len - 1; i >= 0; --i) {
if (a[last_dg + i] > b[i]) return true;
if (a[last_dg + i] < b[i]) return false;
}
return true;
}
void div(int a[], int b[], int c[], int d[]) {//a/b=c,a%b=d
clear(c);
clear(d);
int la, lb;
for (la = LEN - 1; la > 0; --la)
if (a[la - 1] != 0) break;
for (lb = LEN - 1; lb > 0; --lb)
if (b[lb - 1] != 0) break;
if (lb == 0) {
puts("> <");
return;
}
for (int i = 0; i < la; ++i) d[i] = a[i];
for (int i = la - lb; i >= 0; --i) {
while (greater_eq(d, b, i, lb)) {
for (int j = 0; j < lb; ++j) {
d[i + j] -= b[j];
if (d[i + j] < 0) {
d[i + j + 1] -= 1;
d[i + j] += 10;
}
}
c[i] += 1;
}
}
}
数论基础
平凡约数（平凡因数）：对于整数 \(b\ne0\)，\(\pm1\)、\(\pm b\) 是 \(b\) 的平凡约数。当 \(b=\pm1\) 时，\(b\) 只有两个平凡约数
对于整数 \(b\ne 0\)，\(b\) 的其他约数称为真约数（真因数、非平凡约数、非平凡因数）
同余
- 自反性：\(a\equiv a\pmod m\)
- 对称性：若 \(a\equiv b\pmod m\)，则 \(b\equiv a\pmod m\)
- 传递性：若 \(a\equiv b\pmod m\),\(b\equiv c\pmod m\)，则 \(a\equiv c\pmod m\)
- 线性运算：若 \(a,b,c,d\in\mathbf{Z}\),\(m\in\mathbf{N}^*\),\(a\equiv b\pmod m\),\(c\equiv d\pmod m\) 则有：
- \(a\pm c\equiv b\pm d\pmod m\)
- \(a\times c\equiv b\times d\pmod m\)
- 若 \(a,b\in\mathbf{Z},k,m\in\mathbf{N}^*\),\(a\equiv b\pmod m\), 则 \(ak\equiv bk\pmod{mk}\)
- 若 \(a,b\in\mathbf{Z},d,m\in\mathbf{N}^*,d\mid a,d\mid b,d\mid m\)，则当 \(a\equiv b\pmod m\) 成立时，有\(\dfrac{a}{d}\equiv\dfrac{b}{d}\left(\bmod\;{\dfrac{m}{d}}\right)\)
- 若 \(a,b\in\mathbf{Z},d,m\in\mathbf{N}^*,d\mid m，则当 a\equiv b\pmod m 成立时，有 a\equiv b\pmod d\)
- 若 \(a,b\in\mathbf{Z},d,m\in\mathbf{N}^*\)，则当 \(a\equiv b\pmod m\) 成立时，有 \(\gcd(a,m)=\gcd(b,m)\) 若
d能整除m及a,b中的一个，则d必定能整除a,b中的另一个
素数
素数个数：\(\pi(x) \sim \dfrac{x}{\ln(x)}\)
所有大于 3 的素数都可以表示为 \(6n\pm 1\) 的形式
int 范围内的素数间隔是小于 319 ， long long 范围内的素数间隔小于 1525
哥德巴赫猜想
- 关于偶数的哥德巴赫猜想：任一大于2的偶数都可写成两个素数之和。
- 关于奇数的哥德巴赫猜想：任一大于7的奇数都可写成三个质数之和的猜想。
Miller_Rabin
概率性素性测试
\(O(k \log^3n)\) \(k\)为测试次数
质因数分解
唯一分解定理
\(a={p_1}^{\alpha_1}{p_2}^{\alpha_2}\cdots{p_s}^{\alpha_s},p_1<p_2<\cdots<p_s\)
a 的正约数个数为 \(\prod^s_{i=1}(c_i+1)\)
a 的所有正约数和为 \(\prod^s_{i=1}(\sum^{c_i}_{j=0}p^j_i)\)
n 的正因数个数上界是 \(2\sqrt n\)
但实际上这个边界很宽松， \(10^9\) 内的数，正因数最多有 1344 个；\(10^{18}\) 内的数，正因数最多有 103680 个。
\(O(\sqrt n)\)
void divide(int n){
for(int i=2;i*i<=n;i++){
if(n%i==0){
pri[++cnt]=i;
c[cnt]=0;
while(n%i==0){
n/=i;
c[cnt]++;
}
}
}
if(n>1){
pri[++cnt]=n;
c[cnt]=1;
}
}
Pollard Rho 算法
Pollard-Rho 算法是一种用于快速分解非平凡因数的算法（注意！非平凡因子不是素因子）
\(O(n^{\frac{1}{4}})\) 通过倍增可以优化求 \(gcd\) 的用时
P4718
对于每个数字检验是否是质数，是质数就输出 Prime；如果不是质数，输出它最大的质因子是哪个。
#include<bits/stdc++.h>
using namespace std;
#define ll long long
inline ll gcd(ll a,ll b){
return b?gcd(b,a%b):a;
}
inline ll qmul(ll a,ll b,ll P){//不用int128非常耗时
return (__int128)a*b%P;
// ll res=0;
// while(b){
// if(b&1) res=(res+a)%P;
// a=(a+a)%P;
// b>>=1;
// }
// return res;
}
inline ll qpow(ll a,ll b,ll P){
ll res=1;
while(b){
if(b&1) res=qmul(res,a,P);
a=qmul(a,a,P);
b>>=1;
}
return res;
}
ll RandInt(ll l , ll r) {//随机数
static mt19937 Rand(time(0));
uniform_int_distribution<ll> dis(l, r);
return dis(Rand);
}
bool Miller_Rabin(ll n) {
if(n < 3 || n % 2 == 0) return n == 2;
ll a = n - 1 , b = 0;
while(a % 2 == 0) {
a /= 2;
b++;
}
for(int i = 1 , j; i <= 10; i++) {
ll x = RandInt(2 , n - 1);
ll v = qpow(x , a , n);
if(v == 1) continue;
for(j = 0; j < b; j++) {
if(v == n - 1) break;
v = qmul(v,v,n);
}
if(j == b) return 0;
}
return 1;
}
ll Pollard_Rho(ll n) {
if(Miller_Rabin(n)||n==1) return n;//特判,n为1或极大质数是都会卡
ll s = 0 , t = 0;
ll c = RandInt(1 , n - 1);
int step = 0 , goal = 1;
ll value = 1;
auto f = [=](ll x) {
return (qmul(x,x,n)+c)%n;
};
for(goal = 1;; goal <<= 1, s = t , value = 1) {
for(step = 1; step <= goal; step++) {
t = f(t);
value = qmul(value , abs(t - s) , n);
if(step % 127 == 0) {
ll d = gcd(value , n);
if(d > 1) return d;
}
}
ll d = gcd(value , n);
if(d > 1) return d;
}
}
ll Ans;
void Fac(ll n) {//找最大质因子
if(n <= Ans || n < 2) return;
if(Miller_Rabin(n)) {
Ans = max(Ans , n);
return;
}
ll p = n;
while(p == n) p = Pollard_Rho(n);
while((n % p) == 0) n /= p;
Fac(n);
Fac(p);
}
ll N,T;
int main() {
cin >> T;
while(T--) {
Ans = 0;
cin >> N;
Fac(N);
if(Ans == N) cout << "Prime" << endl;//最大质因子为自己则为质数
else cout << Ans << endl;
}
}
质因数分解
queue<ll> aria;
void find(ll n) {
if(n == 1) return;
if(Miller_Rabin(n)) {
aria.push(n);
return;
}
ll p = n;
while(p == n) p = Pollard_Rho(n);
find(p);
find(n/p);
}
int main() {
ll n;
while(~scanf("%lld", &n)) {
find(n);
cout << aria.front();
aria.pop();
while(!aria.empty()) {
cout << "*" << aria.front();
aria.pop();
}
cout << endl;
}
return 0;
}
反素数
反素数：任何小于 n 的正数的约数个数都小于 n 的约数个数，即 n 以内因子最多且最小的数。
const int N=1e6+5,inf=0x3f3f3f3f;
int a[11]={0,2,3,5,7,11,13,17,19,23,29};//打表大法好（质因子种数不超过10）
long long n,ans,tot;//tot为求到的最大的约数个数
void f(long long x,long long now,long long shu,long long num)
{
//x为当前递归的质因子，now为当前求得的数，num为now的约数个数
if(x==11)return ;//递归边界1
long long tmp=1,i;
for(i=1;i<=shu;i++)//当前递归的质因子的个数不超过shu（想不到其他变量名惹...无奈词汇量太小）
{
tmp*=a[x];//tmp暂时存储
if(now*tmp>n)return ;//递归边界2
if(num*(i+1)==tot&&now*tmp<ans)ans=now*tmp;//如果约数个数相同，并且当前得到的数now小于之前得到的数ans，那就更新ans；
if(num*(i+1)>tot)//如果now的约数个数num大于之前求到的最大的约数个数tot，那就更新tot，并且更新ans；
{
tot=num*(i+1);
ans=now*tmp;
}
f(x+1,now*tmp,i,num*(i+1));//往下递归
}
}
模型
求具有给定除数的最小自然数。请确保答案不超过 \(10^{18}\)
unsigned long long p[16] = {
2, 3, 5, 7, 11, 13, 17, 19,
23, 29, 31, 37, 41, 43, 47, 53}; // 根据数据范围可以确定使用的素数最大为53
unsigned long long ans;
unsigned long long n;
// depth: 当前在枚举第几个素数
// temp: 当前因子数量为 num的时候的数值
// num: 当前因子数
// up：上一个素数的幂，这次应该小于等于这个幂次嘛
void dfs(unsigned long long depth, unsigned long long temp,
unsigned long long num, unsigned long long up) {
if (num > n || depth >= 16) return; // 边界条件
if (num == n && ans > temp) { // 取最小的ans
ans = temp;
return;
}
for (int i = 1; i <= up; i++) {
if (temp * p[depth] > ans)
break; // 剪枝：如果加一个这个乘数的结果比ans要大，则必不是最佳方案
dfs(depth + 1, temp = temp * p[depth], num * (i + 1),
i); // 取一个该乘数，进行对下一个乘数的搜索
}
}
int main() {
scanf("%llu", &n);
ans = ~(unsigned long long)0;
dfs(0, 1, 1, 64);
printf("%llu\n", ans);
return 0;
}
约数
\(O(\log n)\)
ll gcd(ll a,ll b){
return b?gcd(b,a%b);
}
ll lcm(ll a,ll b){
return a/gcd(a,b)*b;
}
正因数集合的求法
试除法适用于求单个正整数 n 的正因数集合。
\(O(\sqrt n)\)
void get_factor(int n, vector<int> &factor) {
for (int i = 1;i * i <= n;i++) {
if (!(n % i)) {
factor.push_back(i);
if (i != n / i) factor.push_back(n / i);
}
}
}
倍数法适用于求一个区间 \([1,n]\) 的每个数的正因数集合
此法常用于一些因子相关的求和，如\(\sum^n_{i=1}\sum_{d\mid i} d\)
\(O(n\log n)\)
const int N = 1e6 + 7;
vector<int> factor[N];
void get_factor(int n) {
for (int i = 1;i <= n;i++)
for (int j = 1;i * j <= n;j++)
factor[i * j].push_back(i);
}
拓展欧几里得
常用于求 \(ax+by=\gcd(a,b)\) 的一组可行解
或求解 \(ax+by=m\) 的解，当 \(m|\gcd(a,b)\) 时原方程有整数解为上式的解乘 \(\frac{m}{\gcd(a,b)}\)
对于模线性方程 \(ax≡b (mod n)\) 可以化简为 \(ax + ny = b\)，设 \(d = gcd(a, n)\) 当且仅当 \(b % d == 0\) 时有解，且有d个解
#include<bits/stdc++.h>
using namespace std;
#define ll long long
ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){x=1,y=0;return a;}
ll d=exgcd(b,a%b,x,y);
ll t=x;x=y,y=t-(a/b)*y;
return d;
}
ll a,b,c,x,y;
signed main(){
cin>>a>>b>>c;
ll d=exgcd(a,b,x,y);
cout<<d<<endl;
if(c%d==0){
x*=c/d,y*=c/d;//得到特解
cout<<x<<" "<<y<<endl;
ll p=b/d,q=a/d,k;
if(x<0) k=ceil((1.0-x)/p),x+=p*k,y-=q*k;//x提到最小正整数
else k=(x-1)/p,x-=p*k,y+=q*k;//x降到最小正整数
if(y>0){//有在整数解组
(y-1)/q+1;//解个数
x+(y-1)/q*p;//最大正整x
x;//最小正整x
y;//最大正整y
(y-1)%q+1;//最小正整y
}else{
x;//最小正整x
y+q*(ll)ceil((1.0-y)/q);//最小正整y
}
}
}
由特解求通解
\(x'=x_0+k*\frac{b}{\gcd(a,b)}\)
\(y'=y_0-k*\frac{a}{\gcd(a,b)}\)
数论分块
一般 \(O(\sqrt n)\)
一些小结论：
-
\(\forall a,b,c\in\mathbb{Z},\left\lfloor\frac{a}{bc}\right\rfloor=\left\lfloor\frac{\left\lfloor\frac{a}{b}\right\rfloor}{c}\right\rfloor\)
-
对于常数 n，使得式子
\(\left\lfloor\dfrac ni\right\rfloor=\left\lfloor\dfrac nj\right\rfloor\)
成立的最大的满足 \(i\leq j\leq n\) 的
j的值为 \(\left\lfloor\dfrac n{\lfloor\frac ni\rfloor}\right\rfloor\)。即值 \(\left\lfloor\dfrac ni\right\rfloor\) 所在的块的右端点为 \(\left\lfloor\dfrac n{\lfloor\frac ni\rfloor}\right\rfloor\) -
\(a\%b\) 可以表示为 \(a-b*\lfloor\frac{a}{b}\rfloor\) 此结论可用于高精度取模
数论分块的过程大概如下：考虑和式
\[\sum_{i=1}^nf(i)\left\lfloor\dfrac ni\right\rfloor \]那么由于我们可以知道
\(\left\lfloor\dfrac ni\right\rfloor\) 的值成一个块状分布（就是同样的值都聚集在连续的块中），那么就可以用数论分块加速计算，降低时间复杂度。
利用上述结论，我们先求出 \(f(i)\) 的 前缀和（记作 \(s(i)=\sum_{j=1}^i f(j)\)），然后每次以 \([l,r]=[l,\left\lfloor\dfrac n{\lfloor\frac ni\rfloor}\right\rfloor]\) 为一块，分块求出贡献累加到结果中即可。
N 维数论分块
求含有 \(\left\lfloor\dfrac {a_1}i\right\rfloor\)、\(\left\lfloor\dfrac {a_2}i\right\rfloor\cdots\left\lfloor\dfrac {a_n}i\right\rfloor\) 的和式时，数论分块右端点的表达式从一维的 \(\left\lfloor\dfrac ni\right\rfloor\) 变为 \(\min\limits_{j=1}^n\{\left\lfloor\dfrac {a_j}i\right\rfloor\}\)，即对于每一个块的右端点取最小（最接近左端点）的那个作为整体的右端点。可以借助下图理解：
一般我们用的较多的是二维形式，此时可将代码中 r = n / (n / i) 替换成 r = min(n / (n / i), m / (m / i))
大数阶乘取模
似乎不太算数论分块QAQ\((1 \le N \le 1e10)\)
分块打表,每1e7个数打一个表
const int a[100]={
682498929,491101308,76479948,723816384,67347853,27368307,625544428,199888908,888050723,927880474,
281863274,661224977,623534362,970055531,261384175,195888993,66404266,547665832,109838563,933245637,724691727,
368925948,268838846,136026497,112390913,135498044,217544623,419363534,500780548,668123525,128487469,30977140,
522049725,309058615,386027524,189239124,148528617,940567523,917084264,429277690,996164327,358655417,568392357,
780072518,462639908,275105629,909210595,99199382,703397904,733333339,97830135,608823837,256141983,141827977,
696628828,637939935,811575797,848924691,131772368,724464507,272814771,326159309,456152084,903466878,92255682,
769795511,373745190,606241871,825871994,957939114,435887178,852304035,663307737,375297772,217598709,624148346,
671734977,624500515,748510389,203191898,423951674,629786193,672850561,814362881,823845496,116667533,256473217,
627655552,245795606,586445753,172114298,193781724,778983779,83868974,315103615,965785236,492741665,377329025,
847549272,698611116
};//。。。
const int MOD=1000000007;
int main(){
cin>>n>>p;
if (p==1000000007){
if (n>=p) {
cout<<"0";
return 0;
}
if(n<10000000) now=1;
else now=a[n/10000000-1];
for(int i=n/10000000*10000000+1;i<=n;i++) now=now*i%MOD;
} else{
now=1;
if (n>=p) now=0;
else for(int i=1;i<=n;i++)
now=now*i%p;
}
cout<<now<<endl;
}
P2261[CQOI2007]余数求和
给出正整数 \(n\) 和 \(k\)，请计算
\[G(n, k) = \sum_{i = 1}^n k \bmod i \]\(ans=\sum^n_{i=1} k-i*\lfloor\frac{k}{i}\rfloor=n*k-\sum^n_{i=1}i*\lfloor\frac{k}{i}\rfloor\)
ll ans=n*k;
for(ll l=1,r;l<=n;l=r+1) {
if(k/l!=0) r=min(k/(k/l),n);
else r=n;
ans-=(k/l)*(r-l+1)*(l+r)/2;
}
Calculating
若 \(x\) 分解质因数结果为 \(x=p_1^{k_1}p_2^{k_2}\cdots p_n^{k_n}\)，令\(f(x)=(k_1+1)(k_2+1)\cdots (k_n+1)\)，求 \(\sum_{i=l}^rf(i)\) 对 \(998\,244\,353\) 取模的结果。
明显 \(f(i)\) 表示的是i的因子个数
则有：\(\sum^n_{i=1}f(i)=\sum^n_{i=1}\lfloor\frac{n}{i}\rfloor\)
略证：设 \(i\) 为因子, 有 \(n/i\) 个数含有因子 \(i\)
int block(int n){
int res=0;
for(int l=1,r;l<=n;l=r+1){
r=n/(n/l);
res=(res+(n/l)*(r-l+1)%P)%P;
}
return res;
}
void solve(){
cout<<(block(r)-block(l-1)+P)%P<<endl;
}
欧拉函数
欧拉函数（Euler's totient function），即 \(\varphi(n)\)，表示的是小于等于 \(n\) 和 \(n\) 互质的数的个数
\[\varphi(n)=\sum_{i=1}^n[(i,n)=1] \]比如说 \(\varphi(1) = 1\)
当 \(n\) 是质数的时候，显然有 \(\varphi(n) = n - 1\)
性质:
-
欧拉函数是积性函数。
即如果有 \(\gcd(a, b) = 1\)，那么 \(\varphi(a \times b) = \varphi(a) \times \varphi(b)\)
特别地，当 \(n\) 是奇数时 \(\varphi(2n) = \varphi(n)\)
-
\(n = \sum_{d \mid n}{\varphi(d)}\)
-
若 \(n = p^k\)，其中 \(p\) 是质数，那么 \(\varphi(n) = p^k - p^{k - 1}\)
-
由唯一分解定理，设 \(n = \prod_{i=1}^{s}p_i^{k_i}\)，其中 \(p_i\) 是质数，有 \(\varphi(n) = n \times \prod_{i = 1}^s{\dfrac{p_i - 1}{p_i}}\)
求单个数的欧拉函数,直接根据定义质因数分解的同时求就好了,可以用Pollard Rho优化
int euler_phi(int n) {
int ans = n;
for (int i = 2; i * i <= n; i++)
if (n % i == 0) {
ans = ans / i * (i - 1);
while (n % i == 0) n /= i;
}
if (n > 1) ans = ans / n * (n - 1);
return ans;
}
欧拉定理
\(\gcd(a, m) = 1\)，则 \(a^{\varphi(m)} \equiv 1 \pmod{m}\)
\[a^b\equiv \begin{cases} a^{b\bmod\varphi(p)}, &\gcd(a,\,p)=1\\ a^b,&\gcd(a,\,p)\ne1,\,b<\varphi(p)\\ a^{b\bmod\varphi(p)+\varphi(p)},&\gcd(a,\,p)\ne1,\,b\ge\varphi(p) \end{cases} \pmod p \]扩展欧拉定理
裴蜀定理
设 \(a\),\(b\) 是不全为零的整数，则存在整数 \(x\),\(y\), 使得 \(ax+by=\gcd(a,b)\)
推论
设自然数 a、b 和整数 n。a 与 b 互素。考察不定方程：\(ax+by=n\) 其中 x 和 y 为自然数。如果方程有解，称 n 可以被 a、b 表示。
记 \(C=ab-a-b\)。由 a 与 b 互素，C 必然为奇数。则有结论：
对任意的整数 n，n 与 C-n 中有且仅有一个可以被表示。
即：可表示的数与不可表示的数在区间 [0,C] 对称（关于 C 的一半对称）。0 可被表示，C 不可被表示；负数不可被表示，大于 C 的数可被表示。
noip2017
小凯手中有两种面值的金币，两种面值均为正整数且彼此互素。每种金币小凯都有无数个。在不找零的情况下，仅凭这两种金币，有些物品他是无法准确支付的。现在小凯想知道在无法准确支付的物品中，最贵的价值是多少金币？
ans=ab-a-b
NOIP2005 过河
在河上有一座独木桥，一只青蛙想沿着独木桥从河的一侧跳到另一侧。在桥上有一些石子，青蛙很讨厌踩在这些石子上。由于桥的长度和青蛙一次跳过的距离都是正整数，我们可以把独木桥上青蛙可能到达的点看成数轴上的一串整点：0，1，……，L（其中L是桥的长度）。坐标为0的点表示桥的起点，坐标为L的点表示桥的终点。青蛙从桥的起点开始，不停的向终点方向跳跃。一次跳跃的距离是S到T之间的任意正整数（包括S,T）。当青蛙跳到或跳过坐标为L的点时，就算青蛙已经跳出了独木桥。
题目给出独木桥的长度L（1<=L<=1e9），青蛙跳跃的距离范围S,T<=10，桥上石子的位置。你的任务是确定青蛙要想过河，最少需要踩到的石子数。
路径压缩
假设每次走p或者p+1步.我们知道 \(\gcd(p,p+1)=1\)
我们需要求得一个最小的值tt使得对于所有\(s>t\) 的 \(px+(p+1)y=spx+(p+1)y=s\)一定有非负整数解。根据NOIP2017提高组D1T1的结论，我们可以知道这个数为 \(t=p(p+1)-p-(p+1)t=p(p+1)−p−(p+1)\)。由于本题的最大步长为10，因此 \(t_{max}=9\times10-9-10=71\)
但是要注意，对于 \(s=t\) 这种特殊情况，这种方法是不成立的应为在这种情况下，每次是不能够走p+1步的，因此需要另外特殊判断。
而且有可能跳过终点,所以dp的时候要循环到L+t-1
a,b不互质时不便压缩,因为必须有 \(s|\gcd(a,b)\)
欧拉定理 & 费马小定理
费马小定理:
若 \(p\) 为素数，\(\gcd(a, p) = 1\)，则 \(a^{p - 1} \equiv 1 \pmod{p}\)
欧拉定理:
若 \(\gcd(a, m) = 1\)，则 \(a^{\varphi(m)} \equiv 1 \pmod{m}\)
扩展欧拉定理:
\[a^b\equiv \begin{cases} a^{b\bmod\varphi(p)}, &\gcd(a,\,p)=1\\ a^b,&\gcd(a,\,p)\ne1,\,b<\varphi(p)\\ a^{b\bmod\varphi(p)+\varphi(p)},&\gcd(a,\,p)\ne1,\,b\ge\varphi(p) \end{cases} \pmod p \]无论是费马小定理，还是（扩展）欧拉定理，一个很重要的应用就是降幂，从而将不可能的表达式化为可能
cf Notepad
你有一个本子，你要往上面写全部的长度为\(n\)的\(b\)进制数字，每一页可以写\(c\)个。要求所有数字必须严格不含前导\(0\)。求最后一页上有多少个数字
\(2~\leq~b~<~10^{10^6}~,~1~\leq~n~<~10^{10^6}~,~1~\leq~c~\leq~10^9\)
即求 \((a-1)a^{n-1} mod p\)
#include <bits/stdc++.h>
using namespace std;
#define ll long long
const int N=1000010;
char sa[N],sn[N];
ll p,a,n,phi,q,ans;
int len1,len2;
bool flag;
ll power(ll x,ll k){
ll ans=1;
for (;k;k>>=1,x=x*x%p)
if (k&1) ans=ans*x%p;
return ans;
}
int main(){
scanf("%s %s %lld",sa+1,sn+1,&p);
phi=q=p; len1=strlen(sa+1); len2=strlen(sn+1);
for (ll i=2;i*i<=q;i++)if (!(q%i)){
phi=phi/i*(i-1);
while (!(q%i)) q/=i;
}
if (q>1) phi=phi/q*(q-1);
for (int i=1;i<=len1;i++)
a=(a*10+sa[i]-48)%p;
for (int i=len2;i>=1;i--) //n-1
if (sn[i]==48) sn[i]='9';
else{
sn[i]--;
break;
}
for (int i=1;i<=len2;i++){
n=n*10+sn[i]-48;
if (n>=phi) flag=1;
n%=phi;
}
if (flag) n+=phi;
ans=((a-1)*power(a,n)%p+p)%p; //注意这里可能为负数，所以要加p再模p，被HACK了一次
if (ans) printf("%lld",ans);
else printf("%lld",p);
return 0;
}
逆元
如果一个线性同余方程 \(ax \equiv 1 \pmod b\)，则 \(x\) 称为 \(a \bmod b\) 的逆元，记作 \(a^{-1}\)
\(b\) 为素数时 \(x=a^{b-2}\)
a,b不互质时,有公式: \(x/d \% m = x\%(d*m)/d\)
线性求逆元
inv[0] = 0;
inv[1] = 1;
for (int i = 2; i <= n; ++i) {
inv[i] = (ll)(p - p / i) * inv[p % i] % p;
}
线性求任意\(n\)个数的逆元
首先计算 \(n\) 个数的前缀积，记为 \(s_i\)，然后使用快速幂或扩展欧几里得法计算 \(s_n\) 的逆元，记为 \(sv_n\)。
因为 \(sv_n\) 是 \(n\) 个数的积的逆元，所以当我们把它乘上 \(a_n\) 时，就会和 \(a_n\) 的逆元抵消，于是就得到了 \(a_1\) 到 \(a_{n-1}\) 的积逆元，记为 \(sv_{n-1}\)
同理我们可以依次计算出所有的 \(sv_i\)，于是 \(a_i^{-1}\) 就可以用 \(s_{i-1} \times sv_i\) 求得。
s[0] = 1;
for (int i = 1; i <= n; ++i) s[i] = s[i - 1] * a[i] % p;//
sv[n] = qpow(s[n], p - 2);
for (int i = n; i >= 1; --i) sv[i - 1] = sv[i] * a[i] % p;
for (int i = 1; i <= n; ++i) inv[i] = sv[i] * s[i - 1] % p;
求整数除法小数点后的第n位开始的3位数（0<a,b,n<1000000000）
即求 \(a*10^{n+2}\%(b*1000)/b\)
线性同余方程
\(ax\equiv b\pmod n\)
等价于: \(ax+nk=b\)
有整数解的充要条件为 \(\gcd(a,n) \mid b\)
中国剩余定理
中国剩余定理 (Chinese Remainder Theorem, CRT) 可求解如下形式的一元线性同余方程组（其中 \(n_1\), \(n_2\), \(\cdots\), \(n_k\) 两两互质）：
\[\begin{cases} x &\equiv a_1 \pmod {n_1} \\ x &\equiv a_2 \pmod {n_2} \\ &\vdots \\ x &\equiv a_k \pmod {n_k} \\ \end{cases} \]过程
-
计算所有模数的积 \(n\)
-
对于第 \(i\) 个方程：
- 计算 \(m_i=\frac{n}{n_i}\)
- 计算 \(m_i\) 在模 \(n_i\) 意义下的 逆元 \(m_i^{-1}\)
- 计算 \(c_i=m_im_i^{-1}\)（不要对 \(n_i\) 取模）
-
方程组在模 \(n\) 意义下的唯一解为：
\(x=\sum_{i=1}^k a_ic_i \pmod n\)
ll crt(int k, ll* a, ll* r) {
ll n = 1, ans = 0;
for (int i = 1; i <= k; i++) n = n * r[i];
for (int i = 1; i <= k; i++) {
ll m = n / r[i], b, y;
exgcd(m, r[i], b, y); // b * m mod r[i] = 1
//r[i]为质数可以用费马，r[i]互质用exgcd
ans = (ans + a[i] * m * b % n) % n;//可能被卡,用__int128
}
return (ans % n + n) % n;
}
Garner 算法
\(O(k^2)\)
CRT 的另一个用途是用一组比较小的质数表示一个大的整数。
例如，若 \(a\) 满足如下线性方程组，且 \(a < \prod_{i=1}^k p_i\)（其中 \(p_i\) 两两互质）：
\[\begin{cases} a &\equiv a_1 \pmod {p_1} \\ a &\equiv a_2 \pmod {p_2} \\ &\vdots \\ a &\equiv a_k \pmod {p_k} \\ \end{cases} \]我们可以用以下形式的式子（称作 \(a\) 的混合基数表示）表示 \(a\) ：
\(a = x_1 + x_2 p_1 + x_3 p_1 p_2 + \ldots + x_k p_1 \ldots p_{k-1}\)
Garner算法将用来计算系数 \(x_1, \ldots, x_k\)
令 \(r_{ij}\) 为 \(p_i\) 在模 \(p_j\) 意义下的逆：\(p_i \cdot r_{i,j} \equiv 1 \pmod{p_j}\)
for (int i = 0; i < k; ++i) {
x[i] = a[i];
for (int j = 0; j < i; ++j) {
x[i] = r[j][i] * (x[i] - x[j]);
x[i] = x[i] % p[i];
if (x[i] < 0) x[i] += p[i];
}
}
扩展：模数不互质的情况
ll excrt(int k,ll a[],ll r[]){
ll n=r[1],ans=a[1];//第一个方程的解特判
for(int i=2;i<=k;i++){
ll b=r[i],c=(a[i]-ans%b+b)%b,x,y;//ax≡c(mod b)
ll d=exgcd(n,b,x,y),bg=b/d;
if(c%d!=0) return -1; //判断是否无解
x=(lll)c/d*x%bg;
ans+=x*n;//更新前k个方程组的答案
n*=bg;//M为前k个m的lcm
ans=(ans%n+n)%n;
}
return (ans%n+n)%n;
}
P2480 古代猪文
给出 \(1 \leq G,n \leq 10^9\) 求 \(G^{\sum_{k\mid n}\binom{n}{k}} \bmod 999911659\)
由欧拉定理转化为: \(G^{\sum_{k\mid n}\binom{n}{k}\bmod 999911658} \bmod 999911659\)
明显 \(999911658\) 不是质数但可以分解为: \(2*3*4679*35617\)
构造出同余方程组:
\[\begin{cases} x &\equiv \sum_{k\mid n}\binom{n}{k} \pmod {2} \\ x &\equiv \sum_{k\mid n}\binom{n}{k} \pmod {3} \\ x &\equiv \sum_{k\mid n}\binom{n}{k} \pmod {4679} \\ x &\equiv \sum_{k\mid n}\binom{n}{k} \pmod {35617} \\ \end{cases} \]发现四个数都很小而且都是质数,可以用Lucas定理求出 \(a_i\) 后用CRT求 \(x\)
#include <bits/stdc++.h>
using namespace std;
#define Mod 999911659
#define mod 999911658
#define maxn 40005
#define ll long long
ll n,g,d[maxn],tot,p[10],cnt;
inline ll qpow(ll a,ll k,ll p){
ll res=1;
while(k)
{
if(k&1) res=(res*a)%p;
a=(a*a)%p;
k>>=1;
}
return res%p;
}
ll fac[maxn],inv[maxn];
inline void init(ll p){
fac[0]=1;
for(int i=1;i<p;i++) fac[i]=fac[i-1]*i%p;
inv[p]=0;
inv[p-1]=qpow(fac[p-1],p-2,p);
for(int i=p-2;i>=0;i--) inv[i]=inv[i+1]*(i+1)%p;
}
inline ll C(ll n,ll m,ll p){
if(m>n) return 0;
return fac[n]*inv[m]%p*inv[n-m]%p;
}
inline ll Lucas(ll n,ll m,ll p){
if(m==0) return 1;
return Lucas(n/p,m/p,p)*C(n%p,m%p,p)%p;
}
ll a[10];
inline void calc(int x){
init(p[x]);
for(int i=1;i<=tot;i++)
a[x]=(a[x]+Lucas(n,d[i],p[x]))%p[x];
}
inline ll CRT(){
ll ans=0;
for(int i=1;i<=cnt;i++){
ll M=mod/p[i],t=qpow(M,p[i]-2,p[i]);
ans=(ans+a[i]%mod*t%mod*M%mod)%mod;
}
return (ans+mod)%mod;
}
int main(){
scanf("%lld%lld",&n,&g);
if(g%Mod==0){
printf("0\n");
return 0;
}
ll t=mod;
for(int i=2;i*i<=mod;i++){
if(t%i==0){
p[++cnt]=i;
while(t%i==0) t=t/i;
}
}
if(t!=1) p[++cnt]=t;
for(int i=1;i*i<=n;i++){
if(n%i==0){
d[++tot]=i;
if(i*i!=n) d[++tot]=n/i;
}
}
for(int i=1;i<=cnt;i++) calc(i);
printf("%lld",qpow(g,CRT(),Mod));
}
威尔逊定理
Wilson 定理：对于素数 \(p\) 有 \((p-1)!\equiv -1\pmod p\)
hdu 6608
给出一个质数P，找出小于P的最大的质数N，求出N的阶乘模P。（P∈[1e10,1e14]）
一个数\(n\)若是质数，则有\((n−1)!\equiv n−1\pmod n\). 于是可以先令\(ans=p−1\), 再对\(p−1\)到\(q\)的数对\(p\)求逆元。\(p\)到\(q\)之间的距离大约300以上，Miller Robin大素数判断可以找到最近的素数。
int main(){
cin>>t;
while(t--){
cin >> Q;
ll ans;
for(ll i=Q-1;;i--){
if(Miller_Rabin(i)){
ans=i;
break;
}
}
ll sum=Q-1;
for(ll i=ans+1;i<Q;i++){
sum=mult_mod(sum,pow_mod(i,Q-2,Q),Q);
}
cout<<sum<<endl;
}
}
卢卡斯定理
\(\binom{n}{m}\bmod p = \binom{\left\lfloor n/p \right\rfloor}{\left\lfloor m/p\right\rfloor}\cdot\binom{n\bmod p}{m\bmod p}\bmod p\)
要求 p 的范围不能够太大，一般在 10^5 左右。边界条件：当 m=0 的时候，返回 1。
时间复杂度为 \(O(f(p) + g(n)\log n)\)，其中 \(f(n)\) 为预处理组合数的复杂度，\(g(n)\) 为单次求组合数的复杂度。
ll Lucas(ll n,ll m,ll p) {
if (m == 0) return 1;
return (C(n % p, m % p, p) * Lucas(n / p, m / p, p)) % p;
}
exlucas
P4720
#include<bits/stdc++.h>
using namespace std;
#define ll long long
ll ksm(ll a,ll b,ll P){
ll res=1;
while(b){
if(b&1) res=res*a%P;
a=a*a%P;
b>>=1;
}
return res;
}
ll exgcd(ll a,ll b,ll &x,ll &y){
if(!b){x=1,y=0;return a;}
ll d=exgcd(b,a%b,x,y);
ll t=x;x=y,y=t-(a/b)*y;
return d;
}
ll CRT(int n, ll* a, ll* m) {
ll M = 1, p = 0;
for (int i = 1; i <= n; i++) M = M * m[i];
for (int i = 1; i <= n; i++) {
ll w = M / m[i], x, y;
exgcd(w, m[i], x, y);
p = (p + a[i] * w * x % M) % M;
}
return (p % M + M) % M;
}
ll calc(ll n, ll x, ll P) {
if (!n) return 1;
ll s = 1;
for (ll i = 1; i <= P; i++)
if (i % x) s = s * i % P;
s = ksm(s, n / P, P);
for (ll i = n / P * P + 1; i <= n; i++)
if (i % x) s = i % P * s % P;
return s * calc(n / x, x, P) % P;
}
ll inverse(ll a,ll P){//求逆元，因为P不一定为质数所以用exgcd，用费马会wa一个点
ll x,y;
exgcd(a,P,x,y);
return x;
}
ll multilucas(ll m, ll n, ll x, ll P) {
int cnt = 0;
for (ll i = m; i; i /= x) cnt += i / x;
for (ll i = n; i; i /= x) cnt -= i / x;
for (ll i = m - n; i; i /= x) cnt -= i / x;
return ksm(x, cnt, P) % P * calc(m, x, P) % P * inverse(calc(n, x, P), P) %
P * inverse(calc(m - n, x, P), P) % P;
}
ll exlucas(ll m, ll n, ll P) {
int cnt = 0;
ll p[20], a[20];
for (ll i = 2; i * i <= P; i++) {
if (P % i == 0) {
p[++cnt] = 1;
while (P % i == 0) p[cnt] = p[cnt] * i, P /= i;
a[cnt] = multilucas(m, n, i, p[cnt]);
}
}
if (P > 1) p[++cnt] = P, a[cnt] = multilucas(m, n, P, P);
return CRT(cnt, a, p);
}
int main(){
ll a,b,c;
cin>>a>>b>>c;
cout<<exlucas(a,b,c)<<'\n';
}
二次剩余
一个数 \(a\)，如果不是 \(p\) 的倍数且模 \(p\) 同余于某个数的平方，则称 \(a\) 为模 \(p\) 的 二次剩余。而一个不是 \(p\) 的倍数的数 \(b\)，不同余于任何数的平方，则称 \(b\) 为模 \(p\) 的 二次非剩余。
对二次剩余求解，也就是对常数 \(a\) 解下面的这个方程：
\[x^2 \equiv a \pmod p \]通俗一些，可以认为是求模意义下的开平方 运算。
这里只讨论 \(\boldsymbol{p}\) 为奇素数 的求解方法
对于奇素数 \(p\) 和集合 \(\left\lbrace 1,2,\dots ,p-1\right\rbrace\)，在模 \(p\) 意义下二次剩余的数量等于二次非剩余的数量,即 \(\frac{p-1}{2}\)
欧拉准则
\(n^{\frac{p-1}{2}} \equiv 1\)与 \(n\) 是二次剩余是等价的，由于 \(n^{\frac{p-1}{2}}\) 不为 \(1\) 是只能是 \(-1\) ，那么 \(n^{\frac{p-1}{2}} \equiv -1\)与 \(n\) 是非二次剩余等价。
Cipolla
给出 \(N,p\)，求解方程:
\[x^2 \equiv N \pmod{p} \]多组数据，且保证 \(p\) 是奇素数。
输出共 \(T\) 行。
对于每一行输出，若有解，则按 \(\bmod ~p\) 后递增的顺序输出在 \(\bmod~ p\) 意义下的全部解；若两解相同，只输出其中一个；若无解，则输出Hola!
#include<bits/stdc++.h>
using namespace std;
#define ll long long
struct num {
ll x;// 实部
ll y;// 虚部(即虚数单位√w的系数)
};
ll t,w,n,p;
num mul(num a,num b,ll p) {// 复数乘法
num res;
res.x=( (a.x*b.x%p+a.y*b.y%p*w%p) %p+p)%p;// x = a.x*b.x + a.y*b.y*w
res.y=( (a.x*b.y%p+a.y*b.x%p) %p+p)%p;// y = a.x*b.y + a.y*b.x
return res;
}
ll qpow_r(ll a,ll b,ll p) {// 实数快速幂
ll res=1;
while(b) {
if(b&1) res=res*a%p;
a=a*a%p;
b>>=1;
}
return res;
}
ll qpow_i(num a,ll b,ll p) {// 复数快速幂
num res={1,0};
while(b) {
if(b&1) res=mul(res,a,p);
a=mul(a,a,p);
b>>=1;
}
return res.x%p;// 只用返回实数部分，因为虚数部分没了
}
ll cipolla(ll n,ll p) {
n%=p;
if(qpow_r(n,(p-1)/2,p)==-1+p) return -1;// 据欧拉准则判定是否有解
ll a;
while(1) {// 找出一个符合条件的a
a=rand()%p;
w=( ((a*a)%p-n) %p+p)%p;// w = a^2 - n，虚数单位的平方
if(qpow_r(w,(p-1)/2,p)==-1+p) break;
}
num x={a,1};
return qpow_i(x,(p+1)/2,p);
}
int main() {
srand(time(0));
cin>>t;
while(t--) {
cin>>n>>p;
if(!n) {
printf("0\n");
continue;
}
ll ans1=cipolla(n,p),ans2=-ans1+p;// 另一个解就是其相反数
if(ans1==-1) printf("Hola!\n");
else {
if(ans1>ans2) swap(ans1,ans2);
if(ans1==ans2) printf("%lld\n",ans1);
else printf("%lld %lld\n",ans1,ans2);
}
}
}
Dirichlet 卷积
参考链接:
https://blog.bill.moe/multiplicative-function-sieves-notes/
对于两个数论函数 f(x) 和 g(x)，则它们的狄利克雷卷积得到的结果 h(x) 定义为：
\[h(x)=\sum_{d\mid x}{f(d)g\left(\dfrac xd \right)}=\sum_{ab=x}{f(a)g(b)} \]上式可以简记为：
\[h=f*g \]满足交换律，结合律，分配律
数论函数的积性，在狄利克雷生成函数中的对应具有封闭性
等式的性质: \(f=g\) 的充要条件是 \(f*h=g*h\)，其中数论函数 \(h(x)\) 要满足 \(h(1)\ne 0\)
数论函数
定义域为正整数的函数，值域为复数的函数。
积性函数
规定 \(f(1)=1\)，当 \((a,b)=1\) 时满足 \(f(ab)=f(a)f(b)\) 的函数
特别地: 满足 \(∀a,b，f(ab)=f(a)f(b)\) ,则为完全积性函数
常见的积性函数:
\[\varepsilon(n)=[n=1](完全积性)\\ id(n)=n,id_k(n)=n^k(完全积性)\\ 1(n)=1(完全积性)\\ d(n)=\sum_{i\mid n} 1\\ \sigma(n)=\sum_{i\mid n}i\\ \sigma_k(n)=\sum_{i\mid n}i^k\\ \mu(n)= \begin{cases} 1&n=1\\ 0&n\text{ 含有平方因子}\\ (-1)^k&k\text{ 为 }n\text{ 的本质不同质因子个数}\\ \end{cases}\\ \varphi(n)=\sum_{i=1}^n[(i,n)=1]\\ \]单位元
单位函数 \(\varepsilon\) 是 Dirichlet 卷积运算中的单位元，即对于任何数论函数 \(f\)，都有 \(f*\varepsilon=f\)
\[\varepsilon(x)=[x=1] \]逆元：
对于任何一个满足 \(f(x)\ne 0\) 的数论函数，如果有另一个数论函数 \(g(x)\) 满足 \(f*g=\varepsilon\)，则称 \(g(x)\) 是 \(f(x)\) 的逆元。由等式的性质可知，逆元是唯一的
常见的卷积
注意:转化是等号两侧的转化,双向箭头两侧只是表达方式不同
\[d(n)=\sum_{d\mid n}1\quad\Leftrightarrow\quad d=1*1\\ \sigma(n)=\sum_{d\mid n}d\quad\Leftrightarrow\quad \sigma=d*1\\ \varphi(n)=\sum_{d\mid n}\mu(d)\frac nd\quad\Leftrightarrow\quad\varphi=\mu*id\\ e(n)=\sum_{d\mid n}\mu(d)\quad\Leftrightarrow\quad e=\mu*1\\ \]线性筛
要求
\(f(x)\) 为积性函数
\[n=\prod_{i=1}^kp_i^{a_i} \]线性筛思想
使用最小的 \(p_1\) 去筛掉其他的数。
将 \(n\) 分为三类考虑
- \(n\) 是质数,根据定义得到\(f(i)\)的值
- \(\frac n{p_1}\bmod p_1=0\),说明\(\frac{n}{p_1}\)与\(n\)的质因子集相同，考虑其变化。
- \(\frac n{p_1}\bmod p_1\neq0\),说明\(\frac{n}{p_1}\)与\(p_1\)互质,利用积性函数的性质得到\(f(n)=f(\frac n{p_1})f(p_1)\)
素数线性筛
int vst[maxn],Prime[maxn],cnt=0; //prime
void Prime_Sieve(int n) {
for(int i=2; i<=n; i++) {
if(!vst[i])Prime[++cnt]=i;
for(int j=1; j<=cnt&&i*Prime[j]<=n; j++) {
vst[i*Prime[j]]=1;
if(i%Prime[j]==0)break;
}
}
}
欧拉函数
设\(p_1\)是\(n\)的最小质因子，\(n'=\frac{n}{p_1}\)，在线性筛中，\(n\)通过\(n'\times p1\)被筛掉。
- 当\(n'\bmod p_1=0\),即\(a_1 \gt 1\)时,\(n'\)含有\(n\)的所有质因子，则有：\[\varphi(n)=p_1 \times \varphi(n') \]
- 当\(n'\bmod p_1\neq 0\)即\(a_1=1\)时，\(n'\)与\(p_1\)互素，而欧拉函数是积性函数，故：\[\varphi(n)=(p_1-1) \varphi(n') \]
const int maxn=1000005;
int vst[maxn],Prime[maxn],cnt=0; //prime
int Phi[maxn]; //phi
void Phi_Table(int n) {
Phi[1]=1;
for(int i=2; i<=n; i++) {
if(!vst[i]) {
Prime[++cnt]=i;
Phi[i]=i-1;
}
for(int j=1; j<=cnt&&i*Prime[j]<=n; j++) {
vst[i*Prime[j]]=1;
if(i%Prime[j]==0) {
Phi[i*Prime[j]]=Phi[i]*Prime[j];
break;
}
Phi[i*Prime[j]]=Phi[i]*Phi[Prime[j]];
}
}
}
莫比乌斯函数
https://www.cnblogs.com/icyM3tra/p/16150523.html
当\(n\)为素数时,根据定义有\(\mu(n)=-1\)
设\(p_1\)是\(n\)的最小质因子，\(n'=\frac{n}{p_1}\)，在线性筛中，\(n\)通过\(n'\times p1\)被筛掉。
- 当\(n'\bmod p_1=0\),即\(a_1 \gt 1\)时,由定义可知：\[\mu(n)=0 \]
- 当\(n'\bmod p_1\neq 0\)即\(a_1=1\)时，\(n'\)与\(p_1\)互素，而莫比乌斯函数是积性函数，故：\[\mu(n)=- \mu(n') \]
const int maxn=1000005;
int vst[maxn],Prime[maxn],cnt=0; //prime
int Mobius[maxn]; //mobius
void Mobius_Table(int n) {
Mobius[1]=1;
for(int i=2; i<=n; i++) {
if(!vst[i]) {
Prime[++cnt]=i;
Mobius[i]=-1;
}
for(int j=1; j<=cnt&&i*Prime[j]<=n; j++) {
vst[i*Prime[j]]=1;
if(i%Prime[j]==0) {
Mobius[i*Prime[j]]=0;
break;
}
Mobius[i*Prime[j]]=-Mobius[i];
}
}
}
P4450 双亲数
求\(\sum^A_{i=1}\sum^B_{j=1}[(i,j)==d],1\leq A,B \leq 10^6\)
常用套路:\([(i,j)==1]=\sum_{t|(i,j)}\mu(t)=\sum_{t|i,t|j}\mu(t)\)
\[\begin{aligned} \sum^A_{i=1}\sum^B_{j=1}[(i,j)==d]&=\sum^{A/d}_{i=1}\sum^{B/d}_{j=1}[(i,j)==1]\\ &=\sum^{A/d}_{i=1}\sum^{B/d}_{j=1}\sum_{t|(i,j)}\mu(t)\\ &=\sum^{\min(A,B)/d}_{t=1}\mu(t)\sum^{A/d}_{i=1}[t|i]\sum^{B/d}_{j=1}[t|j]\\ &=\sum^{\min(A,B)/d}_{t=1}\mu(t)(\frac{A}{dt})(\frac{B}{dt}) \end{aligned} \]P1829
求\(\sum^A_{i=1}\sum^B_{j=1}lcm(i,j),1\leq n,m\leq 10^7\)
\[\begin{aligned} \sum_{i=1}^n\sum_{j=1}^m{\rm lcm}(i,j) &=\sum_{d=1}^N\sum_{i=1}^n\sum_{j=1}^m[(i,j)=d]\frac{ij}d \\ &=\sum_{d=1}^Nd\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}[(i,j)=1]ij \\ &=\sum_{d=1}^Nd\sum_{i=1}^{n/d}\sum_{j=1}^{m/d}\sum_{t|i,t|j}\mu(t)ij\\ &=\sum_{d=1}^Nd\sum_{t=1}^{N/d}\mu(t)\sum_{i=1}^{n/d}[t|i]i\sum_{j=1}^{m/d}[t|j]j\\ &=\sum_{d=1}^Nd\sum_{t=1}^{N/d}\mu(t)t^2S(\frac{n}{dt})S(\frac{m}{dt}) \end{aligned} \]其中\(S(n)=\frac{n(n+1)}{2}\)
欧拉反演
大部分莫反题都是用\(\sum_{d|n}\mu(d)代换式子中出现的\)[n=1]$
但在某些情形下，存在另一种做法：用\(\sum_{d|n}\varphi(d)\)代换式子里的\(n\)
P1447
求\(\sum^n_{i=1}\sum^m_{j=1}2(i,j)-1,1\leq n,m \leq 10^5\)
\(\sum\limits_{i=1}^n\sum\limits_{j=1}^m\gcd(i,j) =\sum\limits_{i=1}^n\sum\limits_{j=1}^m\sum\limits_{d\mid i,d\mid j}\varphi(d) =\sum\limits_{d=1}^N\varphi(d)(n/d)(m/d)\)
约数个数函数
若\(n=\prod_{i=1}^kp_i^{a_i}\),则:
\[d(n)=\prod_{i=1}^k(a_i+1) \]当\(n\)为素数时，根据定义有\(d(n)=2\)
设\(p_1\)是\(n\)的最小质因子,\(n'=\frac{n}{p_1}\)
- 当\(n'\bmod p_1=0\),即\(a_1\gt 1\)时,设\(last\)为\(n'\)中\(p_1\)的质数,由约数个数公式可知:\[d(n)=d(n')\times\frac{last+2}{last+1} \]
- 当\(n'\bmod p_1\neq 0\),即\(a_1=1\)时,\(n'\)与\(p_1\)互质,而约数个数函数是积性函数，故：\[\begin{aligned} d(n)&=d(p_1)\times d(n') \\ &=2d(n') \end{aligned} \]
const int maxn=1000005;
int vst[maxn],Prime[maxn],cnt=0; //prime
int d[maxn],Min_Divnum[maxn]; //divisors
void Divisors_Number_Table(int n) {
for(int i=2; i<=n; i++) {
if(!vst[i]) {
Prime[++cnt]=i;
Min_Divnum[i]=1;
d[i]=2;
}
for(int j=1; j<=cnt&&i*Prime[j]<=n; j++) {
vst[i*Prime[j]]=1;
if(i%Prime[j]==0) {
Min_Divnum[i*Prime[j]]=Min_Divnum[i]+1;
d[i*Prime[j]]=d[i]/(Min_Divnum[i]+1)*(Min_Divnum[i]+2);
break;
}
Min_Divnum[i*Prime[j]]=1;
d[i*Prime[j]]=d[i]*d[Prime[j]];
}
}
}
约数和函数
若\(n=\prod_{i=1}^kp_i^{a_i}\),则:
\[\begin{aligned} \sigma(n)&=(1+p_1+p_1^2+\cdots+p_1^{a_1})\times(1+p_2+p_2^2+\cdots+p_2^{a_2})\times\cdots\times(1+p_k+p_k^2+\cdots+p_k^{a_k}) \\ &=\prod_{i=1}^k\sum_{j=0}^{a_i}p_i^j \end{aligned} \]当\(n\)为素数时，根据定义有\(\sigma(n)=n+1\)
设\(p_1\)是\(n\)的最小质因子,\(n'=\frac{n}{p_1}\)
- 当\(n'\bmod p_1=0\),即\(a_1\gt 1\)时,设\(last=p^{a_1-1}_1\)也就是\(n'\)中\(p_1\)为底的最后一个幂，设\(sum=\sum^{a_1-1}_{i=0}p^i_1\),由约数个数公式可知：\[\sigma(n)=\sigma(n')\times \frac{sum+last}{sum} \]
- 当\(n'\bmod p_1\neq 0\),即\(a_1=1\)时,\(n'\)与\(p_1\)互质,而约数个数函数是积性函数，故：\[\begin{aligned} \sigma(n)=\sigma(p_1)\times \sigma(n')\\ =(p_1+1)\sigma(n') \end{aligned} \]
typedef long long LL;
const int maxn=1000005;
int vst[maxn],Prime[maxn],cnt=0; //prime
LL f[maxn],Min_Fac_last[maxn],Min_Fac_sum[maxn]; //divisors
void Divisors_Sum_Table(int n) {
f[1]=1;
for(int i=2; i<=n; i++) {
if(!vst[i]) {
Prime[++cnt]=i;
Min_Fac_last[i]=i;
f[i]=Min_Fac_sum[i]=i+1;
}
for(int j=1; j<=cnt&&i*Prime[j]<=n; j++) {
vst[i*Prime[j]]=1;
if(i%Prime[j]==0) {
Min_Fac_last[i*Prime[j]]=Min_Fac_last[i]*Prime[j];
Min_Fac_sum[i*Prime[j]]=Min_Fac_sum[i]+Min_Fac_last[i*Prime[j]];
f[i*Prime[j]]=f[i]/Min_Fac_sum[i]*Min_Fac_sum[i*Prime[j]];
break;
}
f[i*Prime[j]]=f[i]*(Prime[j]+1);
Min_Fac_last[i*Prime[j]]=Prime[j];
Min_Fac_sum[i*Prime[j]]=Prime[j]+1;
}
}
}
杜教筛
狗都不看
\(O(n^{\frac{2}{3}})\)
P3768 简单的数学题
输入一个整数 \(n\) 和一个整数 \(p\)，你需要求出：
\[\left(\sum_{i=1}^n\sum_{j=1}^n ij \gcd(i,j)\right) \bmod p \]\(n \leq 10^{10}\)，时限4s。
\(5 \times 10^8 \leq p \leq 1.1 \times 10^9\) 且 \(p\) 为质数。
https://www.luogu.com.cn/blog/Soulist/solution-p3768
沙阁筛
\(O(\frac{n^{\frac{3}{4}}}{\ln n})\)
标签：return,数论,ll,long,int,sum,模板,equiv From： https://www.cnblogs.com/EndPB/p/17125118.html