前置知识

初等数论，简单数学知识。

群论

很多数论知识都需要群论证明，我们先简单了解一下群论。

设 $G$ 为一个非空集合，定义运算 $\cdot$，$(G,\cdot)$ 为一个群当且仅当：

1. 封闭性：$\forall a,b \in G,a \cdot b \in G$。
2. 结合律：$\forall a,b,c \in G,(a \cdot b) \cdot c=a \cdot (b \cdot c)$。
3. 单位元：$\exist e \in G,\forall a \in G,e \cdot a=a \cdot e=a$。
4. 逆元：$\forall a \in G,\exist b \in G,a \cdot b=b \cdot a=e$，记作 $b=a^{-1}$。

如果一个群 $G$ 的元素数量是有限的，那么它就是有限群，它的元素数量叫做它的阶，记作 $\lvert G\rvert$。

对于有限群 $G$，如果 $H \subseteq G$，并且 $G,H$ 的运算相同，则称 $H$ 为 $G$ 的子群。

对于 $a \in G$，使 $a^k=e$ 成立的最小正整数 $k$ 称作 $a$ 的阶，记作 $\operatorname{ord}(a)$。

拉格朗日定理

对于有限群 $G$，若 $H$ 为 $G$ 的子群，则 $\lvert H\rvert \mid \lvert G\rvert$，即 $H$ 的阶整除 $G$ 的阶。

以下证明中的 $\cdot$ 运算符号省略。

$\forall a \in G$，定义左陪集：

$$aH=\left\{ah \mid h\in H\right\}$$

$aH$ 的元素与 $H$ 一一对应，因此 $\lvert aH\rvert=\lvert H\rvert$。

取两个不同的陪集 $aH,bH$，假设它们有公共元素，即 $\exist h_1,h_2 \in H,ah_1=bh_2$，则 $a=b h_2h_1^{-1}$。

又因为 $h_2h_1^{-1}\in H$，所以 $a\in bH$。

$\forall ah \in aH,ah=(bh_2h_1^{-1})h=b(h_2h_1^{-1}h)\in bH$，所以 $aH \subseteq bH$。

同理可证：$bH\subseteq aH$，因此 $aH = bH$。与假设“不同的陪集”矛盾，因此 $aH,bH$ 不存在公共元素。

设不同的陪集有 $k$ 个，可得：$\lvert G\rvert=k\lvert H\rvert$。即 $\lvert H\rvert \mid \lvert G\rvert$。

推论：设 $\lvert G \rvert=m$，$\forall a \in G,a^m=e$。

取 $a$ 生成的循环子群 $\langle a\rangle=\left\{a^x \mid x\in \mathbb{N}\right\}$，容易知道其阶等于 $a$ 的阶，且为 $G$ 的子群。由拉格朗日定理得：$\operatorname{ord}(a) \mid m$。因此：

$$a^m=(a^{\operatorname{ord}(a)})^{\frac{m}{\operatorname{ord}(a)}}=e^{\frac{m}{\operatorname{ord}(a)}}=e$$

费马小定理

若 $p$ 为素数，且 $p \nmid a$，则：

$$a^{p-1}\equiv1\pmod p$$

令 $G=(\mathbb{Z}/p\mathbb{Z})^\times$（即所有与 $p$ 互质的剩余类），易知 $\lvert G\rvert=p-1$。

因为 $p \nmid a$，所以 $a \in G$，根据拉格朗日定理的推论，$a^{\lvert G\rvert}=e$，即 $a^{p-1}\equiv 1\pmod p$。

欧拉定理

若 $\gcd(a,n)=1$，则：

$$a^{\varphi(n)}\equiv 1\pmod n$$

令 $G=(\mathbb{Z}/n\mathbb{Z})^\times$，易知 $\lvert G\rvert=\varphi(n)$。

因为 $\gcd(a,n)=1$，所以 $a \in G$，$a^{\lvert G\rvert}=e$，即 $a^{\varphi(n)}\equiv 1\pmod n$。

扩展欧拉定理

$$a^b \equiv \begin{cases} a^{b\ \bmod\ \varphi(p)} & \gcd(a,p)=1 \\ a^b & \gcd(a,p) \neq 1,b<\varphi(p) \\ a^{(b\ \bmod\ \varphi(p))+\varphi(p)} & \gcd(a,p) \neq 1,b \geq \varphi(p) \end{cases} \pmod p$$

$p$ 不一定是质数。

第一项就是欧拉定理，当 $\gcd(a,p) \neq 1$ 的时候，把 $p$ 质因数分解，对于每一个 $p^k$ 证明原式成立，然后用 CRT 合并（CRT 见后文）即可。懒得证了。

裴蜀定理

若 $a,b,c$ 均为整数，$ax+by=c$ 存在整数解 $(x,y)$ 当且仅当 $\gcd(a,b) \mid c$。

设 $S=\left\{ax+by \mid x,y \in \mathbb{Z},ax+by>0\right\}$，易知 $S$ 为 $\mathbb{N^+}$ 的非空子集。设 $S$ 中元素的最小值为 $d$。容易知道必然存在整数 $x,y$ 使得 $d=ax+by$。

可以知道存在整数 $q,r$ 满足 $a=qd+r$，其中 $0 \leq r<d$。

把 $d$ 带入：$r=a-qd=a-q(ax+by)=a(1-qx)+b(-qy)$。

如果 $r>0$，那么 $r \in S$ 并且 $r<d$，矛盾，因此 $r=0$，即 $d \mid a$。

同理，可证：$d \mid b$。因此 $d$ 为 $a,b$ 的公因数。

设 $k$ 是 $a,b$ 的一个公因数，可得 $k \mid d$，即 $k \leq d$，因此 $d$ 为 $a,b$ 的最大公因数，即 $d=\gcd(a,b)$。

设 $m=\frac{c}{\gcd(a,b)}$，显然把 $x,y$ 均乘上 $m$ 就可以得到 $ax+by=c$ 的解，所以 $\gcd(a,b) \mid c$ 的充分性得证；

因为 $d=\gcd(a,b)$，所以 $d \mid (ax+by)$，即 $d \mid c$，$\gcd(a,b) \mid c$ 的必要性得证。

所以，$\gcd(a,b)\mid c$ 是存在整数解的充要条件。

扩展欧几里得算法（exGCD）

求 $ax+by=\gcd(a,b)$ 的一组特解。

由裴蜀定理容易知道一定有解。

如果 $b=0$，$\gcd(a,b)=a$，此时 $x=1,y=0$。

考虑 $a$ 除以 $b$，设 $a=qb+r$，其中 $0 \leq r<b$。可以推出 $\gcd(a,b)=\gcd(b,a\bmod b)=\gcd(b,r)$。

设 $\gcd(b,r)$ 这一步返回了 $x_0,y_0$，得到：

$$bx_0+ry_0=\gcd(a,b) \\ bx_0+(a-qb)y_0=\gcd(a,b) \\ ay_0+b(x_0-qy_0)=\gcd(a,b)$$

其中 $q=\lfloor \frac{a}{b} \rfloor$。

显然我们可以令：

$$x_1=y_0,y_1=x_0-\lfloor \frac{a}{b}\rfloor y_0$$

显然第一步是正确的，后面每一步是由前一步推导出的，所以也都是正确的。

这就是扩展欧几里得算法的流程。

时间复杂度：$O(\log (a+b))$。

扩展中国剩余定理（exCRT）

普通 CRT 感觉不如 exCRT，不讲了。

给定 $a_i,m_i$，解同余方程：

$$\begin{cases} x \equiv a_1 \pmod{m_1} \\ \cdots \\ x \equiv a_n \pmod{m_n} \end{cases}$$

$m_i$ 不一定两两互质。

考虑两两合并，对于 $x \equiv a_1\pmod{m_1}$ 与 $x\equiv a_2\pmod{m_2}$。

这等价于：$x=a_1+y_1m_1$ 与 $x=a_2+y_2m_2$。带入：$a_1+y_1m_1=a_2+y_2m_2$，即 $m_1y_1-m_2y_2=a_2-a_1$，利用扩展欧几里得算法算出一组 $(y_1,y_2)$，算不出来说明整个同余方程无解，算出来了就可以得到 $x\equiv a_1+y_1m_1\pmod{\operatorname{lcm}(m_1,m_2)}$。

时间复杂度：$O(n \log M)$。

大步小步算法（BSGS）

解方程：$a^x \equiv b \pmod p$，$a,p$ 互质。

设 $m=\lceil \sqrt{p}\rceil$，得到 $x=im+j$，其中 $0 \leq i<m,0 \leq j<m$。

带入：

$$a^{im+j}=b \\ a^j=b(a^{-m})^i$$

枚举 $j$，计算 $g^j \bmod p$ 记下来，然后枚举 $i$，计算 $t=b(a^{-m})^i \pmod p$，查找 $t$，找到了就说明此时的 $x$ 是解。如果最终都没找到，说明无解。

时间复杂度：$O(\sqrt{p})$。

扩展大步小步算法（exBSGS）

解方程：$a^x \equiv b\pmod{p}$，$a,p$ 不一定互质。

首先把 $a,b$ 取模，特判 $p=1$ 或 $b=1$ 的情况。

我们计算 $d=\gcd(a,p)$，如果 $b \bmod d \neq 0$，则无解。否则我们可以得到：$\frac{a}{d}a^{x-1}\equiv\frac{b}{d}\pmod{\frac{p}{d}}$。

令：

$$b \leftarrow \frac{b}{d} \\ p \leftarrow\frac{p}{d}$$

不断进行上面的操作，直到 $a,p$ 互质。设我们进行了 $k$ 次操作，$m$ 表示每一轮的 $\frac{a}{d}$ 的乘积，则我们得到了：

$$ma^{x-k}\equiv b\pmod p \\ a^{x-k}\equiv bm^{-1}\pmod p$$

此时 $a,p$ 互质，跑 BSGS 即可。

但是还有一种可能，就是答案 $x<k$，所以在每次计算 $d=\gcd(a,p)$ 之前枚举 $x=0,1,\dots,k-1$ 就好。

时间复杂度：$O(\sqrt{p})$。

Lucas 定理

设 $p$ 为质数，则：

$$\binom{n}{m}\equiv\binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{m}{p}\rfloor}\binom{n\bmod p}{m\bmod p}\pmod p$$

考虑组合意义，设 $n=qp+r,m=sp+t$，$0\leq r<p,0 \leq t<p$，将 $n$ 个球排成一圈，每 $q$ 个分成一组，每组 $p$ 个球，最终剩下 $r$ 个球没有组。设这 $q$ 组为 $G_1,G_2,\dots,G_q$，没有组的为 $G_0$。

考虑从 $n$ 个球中选出 $m$ 个，设在 $G_i$ 中选了 $k_i$ 个球（$1 \leq i \leq q,0 \leq k_i \leq p,0 \leq k_0 \leq r$），则 $\sum k_i+k_0=m$，总方案数为：

$$\binom{n}{m}\equiv\binom{r}{k_0}\prod_{i=1}^q\binom{p}{k_i}\pmod p$$

假设 $\exist 1 \leq i \leq q,1 \leq k_i \leq p-1$，那么有 $\binom{p}{k_i}\equiv0\pmod p$，因此在总方案数中，它们没有贡献，我们只需考虑 $k_i=0$ 或 $k_i=p$ 的部分产生的贡献。又因为 $m=sp+t$，所以一共有 $s$ 个 $k_i=p$，而 $k_0=t$。所以总方案数为：

$$\binom{q}{s}\binom{r}{t}$$

因此：

$$\binom{n}{m}\equiv\binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{m}{p}\rfloor}\binom{n\bmod p}{m\bmod p}\pmod p$$

扩展 Lucas 定理（exLucas）

计算

$$\binom{n}{m} \bmod p$$

其中 $p$ 不一定是质数。

我们考虑把 $p$ 质因数分解，对于每一个 $p^k$ 分别求解，最后用 CRT 或 exCRT 合并。

设 $f(n)$ 表示在 $n!$ 中，含有多少个质因子 $p$。容易知道：

$$f(n)=\sum_{i=1}^\infin \lfloor\frac{n}{p^i}\rfloor=\lfloor\frac{n}{p}\rfloor+f(\lfloor \frac{n}{p}\rfloor)$$

那么 $\binom{n}{m}$ 中包含的质因子 $p$ 的数量 $c=f(n)-f(m)-f(n-m)$。

设 $g(n)$ 表示在 $n!$ 中，除以所有的 $p$ 之后剩下的数对 $p^k$ 取模结果。即：

$$g(n)=\frac{n!}{p^{f(n)}} \bmod p^k$$

推一些式子：

$$n!=\left(\prod_{1 \leq i \leq n,p \nmid i} i\right)\left(\prod_{i=1}^{\lfloor\frac{n}{p}\rfloor}i\cdot p\right) \\ =\left(\prod_{1\leq i\leq n,p \nmid i}i\right)p^{\lfloor\frac{n}{p}\rfloor}\left(\lfloor\frac{n}{p}\rfloor!\right)$$

两边同时除以 $p^{f(n)}$：

$$g(n)=\left(\prod_{1\leq i\leq n,p \nmid i}i\right)g\left(\lfloor\frac{n}{p}\rfloor\right)$$

其中 $g(0)=1$。

考虑快速计算 $\prod_{1\leq i\leq n,p \nmid i}i$，我们知道，$x\equiv x+p^k \pmod{p^k}$，因此我们可以仅仅计算从 $1$ 到 $\frac{n}{p^k}$ 部分中，不被 $p$ 整除的数的乘积。对其求前缀积 $h(n)$ 可得：

$$\left(\prod_{1\leq i\leq n,p \nmid i}i\right)\equiv h\left(p^k\right)^{\lfloor\frac{n}{p^k}\rfloor}h\left(n \bmod \left(p^k\right)\right) \pmod{p^k}$$

综合起来：

$$\binom{n}{m}\equiv p^cg(n)g(m)^{-1}g(n-m)^{-1} \pmod{p^k}$$

然后 CRT 或 exCRT 合并，于是就做完了。时间复杂度：$O(p)$。

整除分块

求：

$$\sum_{i=1}^n \lfloor\frac{n}{i}\rfloor$$

性质：$\lfloor\frac{n}{i}\rfloor$ 取值只有 $O(\sqrt{n})$ 种。

当 $i \leq \sqrt{n}$，$i$ 的取值最多 $\sqrt{n}$ 个，因此原式取值最多 $\sqrt{n}$ 个。

当 $i>\sqrt{n}$，$\lfloor\frac{n}{i}\rfloor<\sqrt{n}$，取值最多 $\sqrt{n}$ 个。

显然每种取值也是连续的。所以我们考虑对于每一段计算。

设当前段左端点为 $l$，右端点为 $r$，那么这一段的值为 $v=\lfloor\frac{n}{l}\rfloor$，$r$ 为满足 $\lfloor\frac{n}{r}\rfloor=v$ 的最大整数。

有如下式子：

$$v\leq\frac{n}{r}<v+1 \\ r \leq \lfloor\frac{n}{l}\rfloor$$

因此 $r=\lfloor\frac{n}{\lfloor\frac{n}{l}\rfloor}\rfloor$。

时间复杂度：$O(\sqrt{n})$。

如果所求式子带权：

$$\sum_{i=1}^n f(i)\lfloor\frac{n}{i}\rfloor$$

处理 $f$ 的前缀和，每一段乘上 $f$ 的区间和即可。

如果是二维整除分块：

$$\sum_{i=1}^{\min(n,m)}\lfloor\frac{n}{i}\rfloor\lfloor\frac{m}{i}\rfloor$$

其他部分不变，只考虑 $r$ 怎么求。类似上面的式子推一下就能得到：$r=\min\left(\lfloor\frac{n}{\lfloor\frac{n}{l}\rfloor}\rfloor,\lfloor\frac{m}{\lfloor\frac{m}{l}\rfloor}\rfloor\right)$。

$k$ 维类似。时间复杂度 $O(k\sqrt{n})$。

莫比乌斯反演

$$\mu(n)=\begin{cases}1&n=1\\(-1)^k&n=p_1p_2p_3\dots p_k,\forall 1 \leq i \leq k,p_i\in\mathbb{P},\forall 1 \leq i<j \leq k,p_i \neq p_j \\0 & \exist p \in \mathbb{P},p^2 \mid n\end{cases}$$

意思就是，如果 $n=1$，那么 $\mu(n)=1$；如果 $n$ 的质因数分解中有 $k$ 个质因子，并且指数都为 $1$，那么 $\mu(n)=(-1)^k$；否则 $\mu(n)=0$。

性质：

$$\sum_{d \mid n}\mu(d)=\begin{cases}1&n=1\\0&n\neq1\end{cases}$$

首先第一项是显然的，我们证明第二项。

容易发现形如 $\prod_{i=1}^k p_i$ 形式的 $d$ 才是有贡献的。设 $n$ 有 $k$ 个不同的质因子，我们从中选取若干个组成了 $d$。如果选取了偶数个，则 $\mu(d)=1$，如果选取了奇数个，则 $\mu(d)=-1$。

设这 $k$ 个不同的质因子中有一个特殊质因子 $p$。我们先在其他 $k-1$ 个质因子中选取，最后决定选不选 $p$。显然，最后选 $p$ 的方案数与不选 $p$ 的方案数相同，但是对于同一种选取其他质因子的方案，选不选 $p$ 之后的集合大小奇偶性不同。因此，选取偶数个的方案数等于选取奇数个的方案数。也就是 $\sum_{d \mid n}\mu(d)=0\ (n\neq1)$。

假设我们已知：

$$f(n)=\sum_{n \mid d} g(d)$$

那么根据反演有：

$$g(n)=\sum_{n\mid d}\mu\left(\frac{d}{n}\right)f(d)$$

$$\sum_{n\mid d}\mu\left(\frac{d}{n}\right)f(d)=\sum_{n \mid d}\mu\left(\frac{d}{n}\right)\sum_{d \mid k}g(k)$$

令 $t=\frac{d}{n}$，有：

$$\sum_{n \mid d}\mu\left(\frac{d}{n}\right)\sum_{d \mid k}g(k)=\sum_{n \mid k}g(k)\sum_{t \mid \frac{k}{n}} \mu(t)=\sum_{n\mid k}g(k)\left[\frac{k}{n}=1\right]=g(n)$$

得证。

$\mu(n)$ 一般使用线性筛求出，代码如下：

mu[1]=1;
for(int i=2;i<=n;i++){
    if(!isp[i]){//是质数
        prm[++pcnt]=i;//记录质数
        mu[i]=-1;
    }
    for(int j=1;j<=pcnt&&i*prm[j]<=n;j++){
        isp[i*prm[j]]=true;//标记合数
        if(i%prm[j]==0){//prm[j]是平方因子
            mu[i*prm[j]]=0;
            break;
        }
        else mu[i*prm[j]]=-mu[i];
    }
}

例题

P5091 【模板】扩展欧拉定理

求 $a^b \bmod m$，$1 \leq a\leq 10^9,1 \leq b \leq 10^{20000000},1 \leq m \leq 10^8$。

模板不讲了。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
ll a,m,mod,x,p;
bool flag;
string b;
ll qpow(ll a,ll b){
    ll res=1%mod;
    while(b){
        if(b&1) res=res*a%mod;
        a=a*a%mod;
        b>>=1;
    }
    return res;
}
int main(){
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    cin>>a>>m>>b;
    a%=m;
    p=mod=m;
    for(ll i=2;i*i<=m;i++){
		if(m%i==0){
			p=p/i*(i-1);
			while(m%i==0) m/=i;
		}
	}
	if(m>1) p=p/m*(m-1);
    for(auto &i:b){
        x=x*10+i-'0';
        if(x>=p){
            flag=true;
            x%=p;
        }
    }
    if(flag) x+=p;
    cout<<qpow(a,x);
    return 0;
}

P5656 【模板】二元一次不定方程 (exgcd)

给定不定方程 $ax+by=c$，无解输出 $-1$；有正整数解输出正整数解数量，正整数解中 $x$ 的最小值、最大值，正整数解中 $y$ 的最小值、最大值；无正整数解有整数解输出整数解中 $x$ 的最小正整数值，$y$ 的最小正整数值。

判无解用裴蜀定理，求特解 $x_0,y_0$ 使用 exgcd，考虑推出通解形式。有：

$$a(x_0+db)+b(y_0-da)=c$$

其中 $db,da$ 为整数。所以，最小的 $d=\frac{1}{\gcd(a,b)}$。令 $d_1=\frac{b}{\gcd(a,b)},d_2=\frac{a}{\gcd(a,b)}$，$s$ 为整数，有通解：

$$x=x_0+sd_1, y=y_0-sd_2$$

当 $x$ 为正整数时，有 $x>0$，即：

$$x_0+sd_1>0 \rightarrow s>-\frac{x_0}{d_1}$$

当 $y$ 为正整数时，有 $y>0$，即：

$$y_0-sd_2>0 \rightarrow s<\frac{y_0}{d_2}$$

分类讨论即可。

代码：

#include<bits/stdc++.h>
using namespace std;
#define LL long long
#define getx(n) (x+(n)*d1)
#define gety(n) (y-(n)*d2)
LL _,a,b,c,x,y,g,d1,d2;
void read(LL &x){
	x=0;
	char ch=getchar();
	int f=1;
	while((ch<'0'||ch>'9')&&ch!='-') ch=getchar();
	if(ch=='-'){
		f=-1;
		ch=getchar();
	}
	while(ch>='0'&&ch<='9'){
		x=x*10+(ch^48);
		ch=getchar();
	}
	x*=f;
	return;
}
void write(const LL &x){
	if(x<0){
		putchar('-');
		write(-x);
		return;
	}
	if(x<=9){
		putchar(x^48);
		return;
	}
	write(x/10);
	putchar(x%10^48);
	return;
}
LL div1(const LL &a,const LL &b){
    if(a%b&&(a<0&&b>0||a>0&&b<0)) return a/b-1;
    return a/b;
}
LL div2(const LL &a,const LL &b){return -div1(-a,b);}
LL exgcd(const LL &a,const LL &b,LL &x,LL &y){
    if(b==0){
        x=1,y=0;
        return a;
    }
    LL g=exgcd(b,a%b,y,x);
    y-=a/b*x;
    return g;
}
void solve(){
    read(a);
	read(b);
	read(c);
    g=exgcd(a,b,x,y);
    if(c%g){
    	printf("-1\n");
        return;
    }
    x*=c/g,y*=c/g;
    d1=b/g,d2=a/g;
    if(div2(y,d2)-1-div1(-x,d1)<=0){
    	write(getx(div1(-x,d1)+1));
    	putchar(' ');
    	write(gety(div2(y,d2)-1));
    	putchar('\n');;
	}
    else{
		write(div2(y,d2)-1-div1(-x,d1));
		putchar(' ');
		write(getx(div1(-x,d1)+1));
		putchar(' ');
		write(gety(div2(y,d2)-1));
		putchar(' ');
		write(getx(div2(y,d2)-1));
		putchar(' ');
		write(gety(div1(-x,d1)+1));
		putchar('\n');
	}
    return;
}
int main(){
    read(_);
    while(_--) solve();
    return 0;
}

P4774 [NOI2018] 屠龙勇士

有 $n$ 条龙，$m$ 把剑，你需要按顺序屠龙。每次屠龙，你需要选择攻击力不高于龙血量 $a_i$ 的一把攻击力最高的剑，如果没有选取攻击力最低的，并攻击龙 $x$ 次。设本次攻击力为 $b_i$，则你需要满足 $b_ix \equiv a_i\pmod{p_i}$，并且 $b_ix\geq a_i$。击败这只龙之后你会获得一把新的剑，旧的剑消失。求最小的正整数 $x$ 或者报告无解。

我们先不考虑这个同余式和不等式，思考应该怎么选剑。显然，我们可以用一个 multiset 和二分查找完成。这样我们就确定了每一个 $b_i$。然后我们跑一个 exCRT 求出最小的 $x$，再利用通解公式求出满足不等式的最小解即可。

简单说一下怎么求通解。设 $P=\operatorname{lcm}_{i=1}^n(p_i)$，我们求出来的解为 $x_0$，则有：

$$x_0\equiv C\pmod P \\ x\equiv C\pmod P$$

求出 $C$，可以得到：$x=sP+C$，其中 $s$ 为整数。假设不等式的解集为 $x \geq T$，则 $sP+C \geq T$，$s \geq \lceil\frac{T-C}{P}\rceil$，答案即为 $\lceil\frac{T-C}{P}\rceil P+C$。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const ll MAXN=1e5+15;
ll _,n,m,atk,lcm,low,m1,m2,g,c,x,y,ans;
ll a[MAXN],p[MAXN],q[MAXN],b[MAXN];
multiset<ll> mset;
void exgcd(const ll &a,const ll &b,ll &x,ll &y){
    if(b==0){
        x=1,y=0;
        return;
    }
    exgcd(b,a%b,y,x);
    y-=a/b*x;
    return;
}
ll inv(const ll &a,const ll &mod){
    ll x,y;
    exgcd(a,mod,x,y);
    return (x%mod+mod)%mod;
}
void solve(){
    mset.clear();
    cin>>n>>m;
    for(ll i=1;i<=n;i++) cin>>a[i];
    for(ll i=1;i<=n;i++) cin>>p[i];
    for(ll i=1;i<=n;i++) cin>>q[i];
    for(ll i=1;i<=m;i++){
        cin>>atk;
        mset.insert(atk);
    }
    for(ll i=1;i<=n;i++){
        auto it=mset.upper_bound(a[i]);
        if(it!=mset.begin()) it--;
        b[i]=*it;
        mset.erase(it);
        mset.insert(q[i]);
    }
    lcm=1,low=ans=0;
    for(ll i=1;i<=n;i++){
        low=max((a[i]+b[i]-1)/b[i],low);
        g=__gcd(b[i],p[i]);
        if(a[i]%g){
            cout<<"-1\n";
            return;
        }
        a[i]/=g,b[i]/=g,p[i]/=g;
        a[i]=(__int128)a[i]%p[i]*inv(b[i]%p[i],p[i])%p[i];
        m1=lcm,m2=p[i];
        g=__gcd(m1,m2);
        c=a[i]-ans;
        if(c%g){
            cout<<"-1\n";
            return;
        }
        m1/=g,m2/=g,c/=g;
        exgcd(m1,m2,x,y);
        lcm*=m2;
        x=((__int128)x*c%m2+m2)%m2;
        ans=((__int128)x*g%lcm*m1%lcm+ans)%lcm;
    }
    if(ans<low){
        low=(low-ans+lcm-1)/lcm;
        ans+=low*lcm;
    }
    cout<<ans<<'\n';
    return;
}
int main(){
    ios::sync_with_stdio(false);
    cin.tie(0);
    cout.tie(0);
    cin>>_;
    while(_--) solve();
    return 0;
}

P4345 [SHOI2015] 超能粒子炮·改

求：

$$\sum_{i=0}^k \binom{n}{i} \bmod 2333$$

设 $p=2333$，再设函数 $f$：

$$f(n,k)=\sum_{i=0}^k \binom{n}{i} \bmod p \\ =\sum_{i=0}^k \binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{i}{p}\rfloor}\binom{n \bmod p}{i \bmod p}\bmod p \\ =\left(\sum_{i=0}^{\left\lfloor\frac{k}{p}\right\rfloor-1}\binom{\lfloor\frac{n}{p}\rfloor}{i}\left(\sum_{j=0}^{p-1}\binom{n \bmod p}{j}\right)+\binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{k}{p}\rfloor}\left(\sum_{j=0}^{k\ \bmod\ p}\binom{n \bmod p}{j}\right)\right) \bmod p \\ =\left(2^{n\ \bmod\ p}\sum_{i=0}^{\left\lfloor\frac{k}{p}\right\rfloor-1}\binom{\lfloor\frac{n}{p}\rfloor}{i}+\binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{k}{p}\rfloor}\left(\sum_{j=0}^{k\ \bmod\ p}\binom{n \bmod p}{j}\right)\right) \bmod p \\ =\left(2^{n\ \bmod\ p}f\left(\left\lfloor\frac{n}{p}\right\rfloor,\left\lfloor\frac{k}{p}\right\rfloor-1\right)+\binom{\lfloor\frac{n}{p}\rfloor}{\lfloor\frac{k}{p}\rfloor}\left(\sum_{j=0}^{k\ \bmod\ p}\binom{n \bmod p}{j}\right)\right) \bmod p$$

加号后面的可以预处理前缀和，递归计算即可。

时间复杂度：$O(p^2+T\log_p^2n)$。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const ll MAXN=3005,MOD=2333;
ll _,n,k;
ll fact[MAXN],invfac[MAXN],pows[MAXN],pre[MAXN][MAXN];
ll qpow(ll a,ll b){
    ll res=1;
    while(b){
        if(b&1) res=res*a%MOD;
        a=a*a%MOD;
        b>>=1;
    }
    return res;
}
ll C(const ll &n,const ll &m){
    if(m<0||n<m) return 0;
    return fact[n]*invfac[m]%MOD*invfac[n-m]%MOD;
}
void init(){
    fact[0]=1;
    for(ll i=1;i<MOD;i++) fact[i]=i*fact[i-1]%MOD;
    invfac[MOD-1]=qpow(fact[MOD-1],MOD-2);
    for(ll i=MOD-1;i;i--) invfac[i-1]=i*invfac[i]%MOD;
    pows[0]=1;
    for(ll i=1;i<MOD;i++) pows[i]=pows[i-1]*2%MOD;
    for(int i=0;i<MOD;i++){
        pre[i][0]=1;
        for(int j=1;j<MOD;j++) pre[i][j]=(pre[i][j-1]+C(i,j))%MOD;
    }
    return;
}
ll lucas(const ll &n,const ll &m){
    if(m<0||n<m) return 0;
    if(n<MOD) return C(n,m);
    return lucas(n/MOD,m/MOD)*lucas(n%MOD,m%MOD)%MOD;
}
ll f(const ll &n,const ll &k){
    ll res=0;
    if(n<MOD){
        for(ll i=0;i<=min(k,MOD-1);i++) res=(res+C(n,i))%MOD;
        return res;
    }
    res=pows[n%MOD]*f(n/MOD,k/MOD-1);
    res=(res+lucas(n/MOD,k/MOD)*pre[n%MOD][k%MOD]%MOD)%MOD;
    return res;
}
void solve(){
    cin>>n>>k;
    k=min(n,k);
    cout<<f(n,k)<<'\n';
    return;
}
int main(){
    cin>>_;
    init();
    while(_--) solve();
    return 0;
}

P3773 [CTSC2017] 吉夫特

给定长度为 $n$ 的序列 $a$，求有多少个长度至少为 $2$ 的不升子序列 $a'$ 满足：

$$\prod_{i=2}^k \binom{a'_{i-1}}{a'_{i}}\equiv1\pmod 2$$

求出结果对 $10^9+7$ 取模。

首先根据 Lucas 定理有：

$$\prod_{i=2}^k \binom{a'_{i-1}}{a'_{i}}\bmod 2 \\ =\prod_{i=2}^k\binom{\left\lfloor\frac{a'_{i-1}}{2}\right\rfloor}{\left\lfloor\frac{a'_{i}}{2}\right\rfloor}\binom{a'_{i-1}\bmod2}{a'_i\bmod2}\bmod2>0$$

说明 $\binom{a'_{i-1}\ \bmod\ 2}{a'_i\ \bmod\ 2}=1$。

继续展开左边的这一坨分式，我们发现整体式子非 $0$ 当且仅当 $a'_{i-1}$ 与 $a'_i$ 在二进制每一位上 $a'_{i-1,j}\geq a'_{i,j}$ 成立。也就是定义 $\&$ 表示按位与，$a'_{i-1}\&a'_{i}=a'_i$。这也可以看作 $a'_i$ 是 $a'_{i-1}$ 的子集。

注意到 $a_i$ 互不相同，设 $dp_i$ 表示结尾是 $i$ 的序列数量，枚举子集 $j$ 刷表法转移即可。

时间复杂度：$O(3^{\log_2\max a_i})$。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const ll MAXN=1e6+16,MOD=1e9+7;
ll n,ans;
ll a[MAXN],dp[MAXN];
int main(){
    cin>>n;
    for(ll i=1;i<=n;i++){
        cin>>a[i];
        for(ll j=a[i]&a[i]-1;j;j=j-1&a[i]) dp[j]=(dp[j]+dp[a[i]]+1)%MOD;
        ans=(ans+dp[a[i]])%MOD;
    }
    cout<<ans;
    return 0;
}

P2260 [清华集训 2012] 模积和

求：

$$\sum_{i=1}^n\sum_{j=1}^m(n\bmod i)(m\bmod j),i\neq j$$

规定 $n\leq m$，容易想到把取模变成整除形式：

$$\sum_{i=1}^n\sum_{j=1}^m(n\bmod i)(m\bmod j),i\neq j \\ =\sum_{i=1}^n\sum_{j=1}^m\left(n-i\left\lfloor\frac{n}{i}\right\rfloor\right)\left(m-j\left\lfloor\frac{m}{j}\right\rfloor\right)-\sum_{i=1}^{n}\left(n-i\left\lfloor\frac{n}{i}\right\rfloor\right)\left(m-i\left\lfloor\frac{m}{i}\right\rfloor\right) \\ =\sum_{i=1}^n\sum_{j=1}^m nm-mi\left\lfloor\frac{n}{i}\right\rfloor-nj\left\lfloor\frac{m}{j}\right\rfloor+ij\left\lfloor\frac{n}{i}\right\rfloor \left\lfloor\frac{m}{j}\right\rfloor-\sum_{i=1}^n nm-mi\left\lfloor\frac{n}{i}\right\rfloor-ni\left\lfloor\frac{m}{i}\right\rfloor+i^2\left\lfloor\frac{n}{i}\right\rfloor \left\lfloor\frac{m}{i}\right\rfloor \\ =n^2m^2-n^2m-m^2\sum_{i=1}^n i\left\lfloor\frac{n}{i}\right\rfloor-n^2\sum_{i=1}^m i\left\lfloor\frac{m}{i}\right\rfloor+\left(\sum_{i=1}^n i\left\lfloor\frac{n}{i}\right\rfloor\right)\left(\sum_{i=1}^mi\left\lfloor\frac{m}{i}\right\rfloor\right)+m\sum_{i=1}^ni\left\lfloor\frac{n}{i}\right\rfloor+n\sum_{i=1}^ni\left\lfloor\frac{m}{i}\right\rfloor-\sum_{i=1}^ni^2\left\lfloor\frac{n}{i}\right\rfloor \left\lfloor\frac{m}{i}\right\rfloor$$

我们知道：

$$\sum_{i=1}^n i^2=\frac{n(n+1)(2n+1)}{6}$$

因此整除分块即可。注意一些循环边界。

时间复杂度：$O(\sqrt{n})$。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const ll MOD=19940417;
ll n,m,ans,x,y;
int main(){
    cin>>n>>m;
    if(n>m) swap(n,m);
    ans=n*n%MOD*m%MOD*m%MOD;
    ans-=n*n%MOD*m%MOD;
    ans=(ans%MOD+MOD)%MOD;
    for(ll i=1,j=1;j<=n;i=j+1){
        if(n/i==0) break;
        j=n/(n/i);
        x=(x+n/i*((i+j)*(j-i+1)/2%MOD)%MOD)%MOD;
    }
    ans-=x*m%MOD*m%MOD;
    ans=(ans%MOD+MOD)%MOD;
    for(ll i=1,j=1;j<=m;i=j+1){
        if(m/i==0) break;
        j=m/(m/i);
        y=(y+m/i*((i+j)*(j-i+1)/2%MOD)%MOD)%MOD;
    }
    ans-=y*n%MOD*n%MOD;
    ans=(ans%MOD+MOD)%MOD;
    ans=(ans+x*y%MOD)%MOD;
    ans=(ans+x*m%MOD)%MOD;
    x=0;
    for(ll i=1,j=1;i<=n;i=j+1){
        if(m/i==0) break;
        j=min(n,m/(m/i));
        x=(x+m/i*((i+j)*(j-i+1)/2%MOD)%MOD)%MOD;
    }
    ans=(ans+x*n%MOD)%MOD;
    x=y=0;
    for(ll i=1,j=1;j<=n;i=j+1){
        if(n/i==0||m/i==0) break;
        j=min(n/(n/i),m/(m/i));
        x=(__int128)j*(j+1)*(j*2+1)/6%MOD;
        x-=(__int128)(i-1)*i*((i-1)*2+1)/6%MOD;
        x=(x%MOD+MOD)%MOD;
        y=(y+(n/i)*(m/i)%MOD*x%MOD)%MOD;
    }
    ans=(ans-y+MOD)%MOD;
    cout<<ans;
    return 0;
}

P2522 [HAOI2011] Problem b

给出 $n$ 次询问，每次给定 $a,b,c,d,k$，求：

$$\sum_{i=a}^b\sum_{j=c}^d \left[\gcd(i,j)=k\right]$$

设：

$$f(a,b,c,d)=\sum_{i=a}^b\sum_{j=c}^d \left[\gcd(i,j)=k\right]$$

容斥则有：

$$ans=f(1,b,1,d)-f(1,a-1,1,d)-f(1,b,1,c-1)+f(1,a-1,1,c-1)$$

钦定 $n \leq m$，所以只需要求：

$$T(n,m)=\sum_{i=1}^n\sum_{j=1}^m\left[\gcd(i,j)=k\right] \\$$

设 $F(d)$ 表示最大公约数为 $d$ 的数对数量，$G(d)$ 表示满足 $d\mid \gcd(i,j)$ 的数对数量。得到：

$$G(d)=\sum_{d\mid k}F(k) \\ G(d)=\left\lfloor\frac{n}{d}\right\rfloor\left\lfloor\frac{m}{d}\right\rfloor$$

明显的莫比乌斯反演：

$$F(k)=\sum_{k\mid d}\mu\left(\frac{d}{k}\right)G(d) \\ =\sum_{i=1}^{\left\lfloor\frac{n}{k}\right\rfloor}\mu(i)\left\lfloor\frac{n}{ki}\right\rfloor\left\lfloor\frac{m}{ki}\right\rfloor$$

令 $N=\left\lfloor\frac{n}{k}\right\rfloor,M=\left\lfloor\frac{m}{k}\right\rfloor$：

$$F(k)=\sum_{i=1}^{N}\mu(i)\left\lfloor\frac{N}{i}\right\rfloor\left\lfloor\frac{M}{i}\right\rfloor$$

整除分块即可。

时间复杂度：$O(N+n\sqrt{N})$。

代码：

#include<bits/stdc++.h>
using namespace std;
#define ll long long
const ll MAXN=5e4+54;
ll _,a,b,c,d,k,pcnt;
ll mu[MAXN],pre[MAXN],prm[MAXN];
bool isp[MAXN];
void init(){
    mu[1]=1;
    for(ll i=2;i<=50000;i++){
        if(!isp[i]){
            prm[++pcnt]=i;
            mu[i]=-1;
        }
        for(ll j=1;j<=pcnt&&i*prm[j]<=50000;j++){
            isp[i*prm[j]]=true;
            if(i%prm[j]==0){
                mu[i*prm[j]]=0;
                break;
            }
            else mu[i*prm[j]]=-mu[i];
        }
    }
    for(ll i=1;i<=50000;i++) pre[i]=mu[i]+pre[i-1];
    return;
}
ll f(const ll &n,const ll &m){
    ll tn=n/k,tm=m/k,res=0;
    if(tn>tm) swap(tn,tm);
    for(ll i=1,j=1;i<=tn;i=j+1){
        if(tn/i==0||tm/i==0) break;
        j=min(tn/(tn/i),tm/(tm/i));
        res+=(pre[j]-pre[i-1])*(tn/i)*(tm/i);
    }
    return res;
}
void solve(){
    cin>>a>>b>>c>>d>>k;
    cout<<f(b,d)-f(a-1,d)-f(b,c-1)+f(a-1,c-1)<<'\n';
    return;
}
int main(){
    init();
    cin>>_;
    while(_--) solve();
    return 0;
}

P2257 YY的GCD

求：

$$\sum_{i=1}^n\sum_{j=1}^m\operatorname{isp}(\gcd(i,j))$$

其中 $\operatorname{isp}(x)$ 表示 $x$ 是否为质数。

依旧令 $n\leq m$，推式子：

$$\sum_{i=1}^n\sum_{j=1}^m\operatorname{isp}(\gcd(i,j))\\ =\sum_{p\in\mathbb{P}}\sum_{i=1}^n\sum_{j=1}^m\left[\gcd(i,j)=p\right] \\$$

显然和上一题一模一样。令 $N=\left\lfloor\frac{n}{p}\right\rfloor,M=\left\lfloor\frac{m}{p}\right\rfloor$：

$$ans=\sum_{p\in\mathbb{P}}\sum_{i=1}^{N}\mu(i)\left\lfloor\frac{N}{i}\right\rfloor\left\lfloor\frac{M}{i}\right\rfloor\\ =\sum_{i=1}^n\left\lfloor\frac{n}{i}\right\rfloor\left\lfloor\frac{m}{i}\right\rfloor\left(\sum_{p\mid i,p\in\mathbb{P}}\mu\left(\frac{i}{p}\right)\right)\\$$

右边的一坨令其为 $f(i)$，预处理即可。

时间复杂度：$O(n\log\log n+T\sqrt{n})$。

代码：

#include<bits/stdc++.h>
using namespace std;
const int MAXN=1e7+17;
int _,n,m,pcnt;
long long ans;
long long mu[MAXN],pre[MAXN],f[MAXN],prm[1000005];
bool isp[MAXN];
void init(){
    mu[1]=1;
    for(int i=2;i<MAXN;i++){
        if(!isp[i]){
            prm[++pcnt]=i;
            mu[i]=-1;
        }
        for(int j=1;j<=pcnt&&i*prm[j]<MAXN;j++){
            isp[i*prm[j]]=true;
            if(i%prm[j]==0){
                mu[i*prm[j]]=0;
                break;
            }
            else mu[i*prm[j]]=-mu[i];
        }
    }
    for(int i=1;i<=pcnt;i++) for(int j=prm[i];j<MAXN;j+=prm[i]) f[j]+=mu[j/prm[i]];
    for(int i=1;i<MAXN;i++) pre[i]=f[i]+pre[i-1];
    return;
}
void solve(){
    ans=0;
    cin>>n>>m;
    if(n>m) swap(n,m);
    for(int i=1,j=1;i<=n;i=j+1){
        if(n/i==0||m/i==0) break;
        j=min(n/(n/i),m/(m/i));
        ans+=(pre[j]-pre[i-1])*(n/i)*(m/i);
    }
    cout<<ans<<'\n';
    return;
}
int main(){
    init();
    cin>>_;
    while(_--) solve();
    return 0;
}