最新要闻

广告

手机

iphone11大小尺寸是多少？苹果iPhone11和iPhone13的区别是什么？

iphone11大小尺寸是多少？苹果iPhone11和iPhone13的区别是什么？

 警方通报辅警执法直播中被撞飞：犯罪嫌疑人已投案

警方通报辅警执法直播中被撞飞：犯罪嫌疑人已投案

家电

  1. 首页
  2. 移动通信
  3. 内容

数论模板

来源：博客园

数学

配合oiwiki:


(资料图)

https://oi-wiki.org/math/

位运算

  1. int __builtin_ffs(int x)：返回 x 的二进制末尾最后一个 1 的位置，位置的编号从 1 开始（最低位编号为 1 ）。当 x 为 0 时返回 0 。
  2. int __builtin_clz(unsigned int x)：返回 x 的二进制的前导 0 的个数。当 x 为 0 时，结果未定义。
  3. int __builtin_ctz(unsigned int x)：返回 x 的二进制末尾连续 0 的个数。当 x 为 0 时，结果未定义。
  4. int __builtin_clrsb(int x)：当 x 的符号位为 0 时返回 x 的二进制的前导 0 的个数减一，否则返回 x 的二进制的前导 1 的个数减一。
  5. int __builtin_popcount(unsigned int x)：返回 x 的二进制中 1 的个数。
  6. 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

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；如果不是质数，输出它最大的质因子是哪个。

#includeusing namespace std;#define ll long longinline 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 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 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*tmptot)//如果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};  // 根据数据范围可以确定使用的素数最大为53unsigned 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 &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 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个解

#includeusing namespace std;#define ll long longll 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<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<

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<

欧拉函数

欧拉函数（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 using namespace std;#define ll long longconst 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*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}\]

过程

  1. 计算所有模数的积 \(n\)

  2. 对于第 \(i\) 个方程：

    1. 计算 \(m_i=\frac{n}{n_i}\)
    2. 计算 \(m_i\) 在模 \(n_i\) 意义下的 逆元 \(m_i^{-1}\)
    3. 计算 \(c_i=m_im_i^{-1}\)（不要对 \(n_i\) 取模）

  3. 方程组在模 \(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 using namespace std;#define Mod 999911659#define mod 999911658#define maxn 40005#define ll long longll 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=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

卢卡斯定理

\(\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;}

exlucasP4720

#includeusing namespace std;#define ll long longll 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<

二次剩余

一个数 \(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!

#includeusing namespace std;#define ll long longstruct 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\) 分为三类考虑

  1. \(n\) 是质数,根据定义得到\(f(i)\)的值
  2. \(\frac n{p_1}\bmod p_1=0\),说明\(\frac{n}{p_1}\)与\(n\)的质因子集相同，考虑其变化。
  3. \(\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; //primevoid 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; //primeint Phi[maxn]; //phivoid 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; //primeint Mobius[maxn]; //mobiusvoid 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; //primeint d[maxn],Min_Divnum[maxn]; //divisorsvoid 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; //primeLL f[maxn],Min_Fac_last[maxn],Min_Fac_sum[maxn]; //divisorsvoid 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})\)

关键词： 欧拉定理 数论函数 二次剩余