[bzoj5472] 数列

Description

输入一个长度为n的数组{ai}(1 <= i <= n)

问有多少个长度为n的数组{xi}(1 <= i <= n)，满足1 <= xi <= ai。

并且相邻两项的最大公约数小于等于k。

换句话说，对于1<i <= n，满足gcd(xi−1,xi) <= k。

问这样的数组{xi}有多少个，答案对1000000007取模。

Input

第一行两个整数n,k。接下来一行n个整数，表示数组{ai}。

1 <= ai,k <= 1000000。Sigma(ai) <= 1000000。

Output

一行一个整数表示这样的数组的个数，对1000000007取模。

Sample Input

9 2
1 2 3 4 5 6 7 8 9

Sample Output

168852

Solution

考虑\(dp\)，

设\(f_{i,j}\)表示当前选到第\(i\)个数了，这个数选的是\(j\)。

显然可以得到转移：

\[f_{i,j}=\sum_{x=1}^{lim_{i-1}}f_{i-1,x}[\gcd(j,x)\leqslant k] \]

这个是\(O(n^2)\)的，不能接受，考虑针对后面的\(\gcd\)化简式子，莫比乌斯反演一波可以得到：

\[f_{i,j}=\sum_{T|j}\sum_{t|T}\mu(t)[\frac{T}{t}\leqslant k]\sum_{x=1}^{\frac{lim_{i-1}}{T}}f_{i-1,xT} \]

中间都是一些套路的化简，这里不赘述了。

然后设：

\[g(n)=\sum_{t|n}\mu(t)[\frac{n}{t}\leqslant k] \]

这个可以\(O(n\log n)\)预处理出来。

设：

\[h(x)=\sum_{i=1}^{\frac{lim_{i-1}}{x}}f_{i-1,ix} \]

这个每次转移都是\(O(lim\cdot \log lim)\)，总复杂度\(O(n\log n)\)。

所以转移可以写成这样：

\[f_{i,j}=\sum_{t|j}g(t)h(t) \]

这个可以枚举约数\(O(n\log n)\)。

所以，总复杂度为\(O(n\log n)\)。

#pragma GCC optimize(3)
#include<bits/stdc++.h>
using namespace std;

const int mod = 1000000007;

void read(int &x) {
    x=0;int f=1;char ch=getchar();
    for(;!isdigit(ch);ch=getchar()) if(ch=='-') f=-f;
    for(;isdigit(ch);ch=getchar()) x=x*10+ch-'0';x*=f;
}
 
void print(int x) {
    if(x<0) putchar('-'),x=-x;
    if(!x) return ;print(x/10),putchar(x%10+48);
}
void write(int x) {if(!x) putchar('0');else print(x);putchar('\n');}

const int maxn = 2e6+10;

vector <int > f[maxn];

int g[maxn],mu[maxn],pri[maxn],tot,vis[maxn],k,n,lim[maxn],h[maxn],G[maxn];

void sieve() {
	mu[1]=1;
	for(int i=2;i<maxn;i++) {
		if(!vis[i]) pri[++tot]=i,mu[i]=-1;
		for(int j=1;j<=tot&&i*pri[j]<maxn;j++) {
			vis[i*pri[j]]=1;
			if(i%pri[j]==0) break;
			mu[i*pri[j]]=-mu[i];
		}
	}
	for(int t=1;t<maxn;t++)
		for(int d=t;d<maxn&&d<=k*t;d+=t) g[d]+=mu[t];
}

int main() {
	read(n),read(k);
	sieve();
	for(int i=1;i<=n;i++) {
		read(lim[i]);
		for(int j=0;j<=lim[i]+1;j++) f[i].push_back(0);
	}
	for(int i=1;i<=lim[1];i++) f[1][i]=1;
	for(int i=2;i<=n;i++) {
		for(int j=1;j<=lim[i];j++) h[j]=0;   // ATTENTION lim[i] !
		for(int j=1;j<=lim[i-1];j++)
			for(int u=j;u<=lim[i-1];u+=j) h[j]=(h[j]+f[i-1][u])%mod;
		for(int t=1;t<=lim[i];t++)
			for(int d=t;d<=lim[i];d+=t) 
				f[i][d]=(1ll*f[i][d]+1ll*h[t]*g[t]%mod)%mod;
	}
	int ans=0;
	for(int i=1;i<=lim[n];i++) ans=(ans+f[n][i])%mod;
	write((ans%mod+mod)%mod);
	return 0;
}