\(3\leq n\leq 8,2\leq L \leq 10^9,1\leq X<Y\leq L,1\leq K \leq n\)
显然最大的数和其他的数字可以分开处理 , 即先计算出其他数字的和再计算最大值的取值方案 .
考虑用容斥解决相同最多 \(K\) 个的限制 .
考虑把 \(n-1\) 各元素分成若干个子集 , 每个子集内部元素强制相等 , 子集之间不做限制 .
首先计算大小为 \(i\) 的集合的容斥系数 .
实际大小为 \(i\) 的集合应该贡献的系数是 \(f_i=[i\leq K]i!\) 即先判断这个集合是否合法 , 然后因为计算的是无序排列数 , 所以需要先乘上每个连通块大小的阶乘 , 最后再除以 \((n-1)\) 的阶乘 .
那么设容斥系数是 \(g_i\) , \(G,F\) 分别是两者的 EGF , 那么有 \(F=\exp(G),G=\ln(F)\)
那么容斥系数就求完了 .
然后需要处理的问题是给定若干个集合大小 , 每个集合内部元素相等且取值在 \([0,X]\) 之间 , 设这些集合所有元素的和为 \(sum\), 那么一种方案的贡献为 \(\max(0,\min(L+1,sum)-Y)\), 求所有方案贡献之和 .
再次容斥哪些集合超出了 \(X\) 的限制 , 现在每个集合的取值没有限制 .
设我们钦定超过的集合的大小之和为 \(over\)
那么现在方案的 OGF 就是
\(\displaystyle x^{over*(X+1)}\prod_{siz} \frac{1}{1-x^{siz}}\) 其中 \(siz\) 枚举了该划分方案中的全部大小 .
考虑分段计算贡献
设方案的 OGF 为 \(F\)
那么上述贡献即为
\(\displaystyle\sum\limits_{i=Y+1}^{L+1}[x^i]F(x)(i-Y)\)
\(\displaystyle\sum\limits_{i=L+2}^{+\infty}[x^i]F(x)(L+1-Y)\)
我们只需要求解以下形式的问题即可
求
\(\displaystyle\sum\limits_{i=l}^{r}[x^i]F(x)\) 和 \(\sum\limits_{i=l}^{r}[x^i]F(x)i\)
对于后者只需求导后乘 \(x\) 即可转化为前者 .
前者可以将 \(F\) 乘上 \(\frac{1}{1-x}\) 将系数转化为前缀和就变成了求
\([x^n]\frac{F(x)}{1-x}\)
使用分式求第 \(n\) 项的算法即可 \(O(deg^2\log n)\) 求出 , 其中 \(deg^2\) 是暴力多项式乘法的复杂度 .
注意到求和上界为 \(+\infty\) 不好处理 , 实际上只需统一取一个足够大的数字即可 , 比如 \(nL\) .
#include<bits/stdc++.h> typedef long long ll; using namespace std; int read() { int ret=0;bool f=0;char c=getchar(); while(c>'9'||c<'0')f|=(c=='-'),c=getchar(); while(c>='0'&&c<='9')ret=(ret<<3)+(ret<<1)+(c^48),c=getchar(); return f?-ret:ret; } const int maxn=20,mod=1e9+7; int qpow(int a,int b) { int ret=1; for(;b;b>>=1,a=(ll)a*a%mod)if(b&1)ret=(ll)ret*a%mod; return ret; } int n,K,L,Y,X; struct poly { vector<int>v; poly()=default; poly(initializer_list<int>_v){v=_v;} void pb(int x){v.push_back(x);} int& operator[](const int &i){return v[i];} void set(int l){v.resize(l);} int len(){return v.size();} poly operator *(poly x) { poly ret; ret.set(len()+x.len()-1); for(int i=0;i<len();i++) for(int j=0;j<x.len();j++)(ret[i+j]+=(ll)v[i]*x[j]%mod)%=mod; return ret; } void operator *=(poly x){poly ret=*this*x;*this=ret;} poly getd() { poly ret; for(int i=1;i<len();i++)ret.pb((ll)i*v[i]%mod); if(ret.v.empty())ret.pb(0); return ret; } void print() { for(int &i:v)printf("%d ",i); putchar('\n'); } }; int coef[maxn],ifac[maxn],fac[maxn]; void init() { poly bas={0},now={1}; fac[0]=1;for(int i=1;i<=n;i++)fac[i]=(ll)fac[i-1]*i%mod; ifac[n]=qpow(fac[n],mod-2);for(int i=n-1;i>=0;i--)ifac[i]=(ll)ifac[i+1]*(i+1)%mod; for(int i=1;i<=K;i++)bas.pb(1); for(int i=1;i<=n-1;i++) { now*=bas; int inv=qpow(i,mod-2); if(now.len()>n)now.set(n); if(i&1)for(int j=1;j<now.len();j++)(coef[j]+=(ll)now[j]*inv%mod)%=mod; else for(int j=1;j<now.len();j++)(coef[j]+=mod-(ll)now[j]*inv%mod)%=mod; } for(int i=1;i<=n-1;i++)coef[i]=(ll)coef[i]*fac[i]%mod; } int calc(poly f,poly g,ll n)// [x^n] f/g { for(;n;n>>=1) { poly tmp=g; for(int i=1;i<tmp.len();i+=2)tmp[i]=(mod-tmp[i])%mod; f*=tmp;g*=tmp; int to=-1;for(int i=n&1;i<f.len();i+=2)f[i/2]=f[i],to=i/2;f.set(to+1); to=-1;for(int i=0;i<g.len();i+=2)g[i/2]=g[i],to=i/2;g.set(to+1); } return f.len()?(ll)f[0]*qpow(g[0],mod-2)%mod:0; } int ans; void calc(vector<int>B) { int coe=1; static int cnt[maxn]; for(int &i:B)coe=(ll)coe*coef[i]%mod,cnt[i]++; for(int i=1;i<=n-1;i++){coe=(ll)coe*ifac[cnt[i]]%mod;cnt[i]=0;} for(int &i:B)coe=(ll)coe*ifac[i]%mod; poly down={1}; for(int &i:B) { poly now;now.pb(1);for(int j=1;j<=i-1;j++)now.pb(0);now.pb(mod-1); down*=now; } static int f[2][maxn],g[2][maxn]; memset(f,0,sizeof f);f[0][0]=1; for(int &i:B) { memset(g,0,sizeof g); for(int j=0;j<=n-1;j++) { if(f[0][j])g[1][j+i]+=f[0][j],g[0][j]+=f[0][j]; if(f[1][j])g[0][j+i]+=f[1][j],g[1][j]+=f[1][j]; } memcpy(f,g,sizeof g); } int ret=0; poly pre=down*(poly){1,mod-1},pre0up=down.getd()*(poly){0,1},pre0down=down*down*(poly){1,mod-1}; for(int &i:pre0down.v)i=(mod-i)%mod; for(int i=0;i<=n-1;i++) { if(!(f[0][i]-f[1][i]))continue; ll be=max(L-(ll)i*(X+1),-1ll),ed=max(-1ll,(ll)n*L-(ll)i*(X+1)); int tmp=0; if(ed>=0)(tmp+=(ll)calc((poly){1},pre,ed)*(L+1-Y)%mod)%=mod; if(be>=0)(tmp+=mod-(ll)calc((poly){1},pre,be)*(L+1-Y)%mod)%=mod; be=max(Y-(ll)i*(X+1),-1ll);ed=max(-1ll,L-(ll)i*(X+1)); if(ed>=0) { (tmp+=mod-(ll)calc((poly){1},pre,ed)*(Y-(ll)i*(X+1)%mod+mod)%mod)%=mod; (tmp+=(ll)calc(pre0up,pre0down,ed))%=mod; } if(be>=0) { (tmp+=(ll)calc((poly){1},pre,be)*(Y-(ll)i*(X+1)%mod+mod)%mod)%=mod; (tmp+=mod-(ll)calc(pre0up,pre0down,be))%=mod; } (ret+=(ll)tmp*(f[0][i]-f[1][i]+mod)%mod)%=mod; } ret=(ll)ret*coe%mod; (ans+=ret)%=mod; } void dfs(int res,int lim,vector<int>now) { if(!res){calc(now);return;} if(lim<res)dfs(res,lim+1,now); if(lim<=res){now.push_back(lim);dfs(res-lim,lim,now);} } int main() { freopen("dish.in","r",stdin); freopen("dish.out","w",stdout); n=read();L=read();Y=read();X=read();K=read(); init(); dfs(n-1,1,{}); printf("%d\n",ans); return 0; }