loj6077

题解：

网上的做法好像都是容斥

那就先说一下容斥

首先问题等价于求下面这个式子的方案数

$$\sum_{i=1}^{n} ai (0<ai<i) =k$$

直接$dp$复杂度是$nk$的,无法接受

我们考虑容斥

答案=所有-至少1个不满足+至少2个不满足-....

而不满足等价于$ai>=i$ 那么我们只需要令$b[i]=a[i]-i$就可以与其他等价了

实际上这样子看上去这个容斥并没有用处

因为要$dp$出容斥系数看上去至少是$O(n^2)$

现在我们需要算出，用i个不相同的数相加=k的方案数

我们注意一件事情，i的取值只有$\sqrt{k}$，所以如果状态设计得好，是可以的

要求数不同有一个经典套路就是 把操作转化成

1.给当前数组内全体数+1

2.给当前数组内全体数+1再在末尾放一个1

于是$f[i][j]$表示数组内有i个数，和为j的方案数 $$f[i][j]=f[i-1][j-i]+f[i][j-i]$$复杂度$n\sqrt{n}$

然后你写一下，肯定会发现出问题了（我没写这种做法但下面这种的时候遇到了同样的问题）

因为这么做可能有数>n了

怎么保证没有数>n呢，我们只需要令每个时刻的$f[i][j]$都是<=n的

那么出现大于n的只可能是n+1，$f[i][j]-=f[i-1][j-n-1]$就可以了

还是比较难想的。。

这题还有个做法是生成函数，在思维上就简单很多

$$f(x)=\frac{\prod_{i=1}^{n} (1-x^i)}{(1-x)^n}$$ 

然后要求这个东西，暴力是n^2logn的

然后乘法比较显然的是两边取ln，得到

$$ln(f(x))=\sum_{i=1}^{n} {ln(1-x^i)} - n*ln(1-x)$$

这样子我们发现要解决的就是$ln(1-x^i)$

这个东西就两项但我们要用$nlogn$的多项式求逆是不是太浪费了

于是打表找规律，可以发现（不过要是不知道结论考场谁有时间去打个多项式求ln找规律啊。。）

$ln(1-x)=\frac{1}{1}x+\frac{1}{2}x^2+\frac{1}{3}x^3+...$

$ln(1-x^2)=\frac{1}{1}x^2+\frac{1}{2}x^4+\frac{1}{3}x^6+...$

下面的规律同理

于是我们可以利用筛法在$nlogn$时间内求出

当然由于$$f(x)=\sum {d|x} {}{ \ \ \ g1(d)g2(\frac{x}{d}) }$$ 其中$g1(i)=\frac{1}{i},g2(i)=1$

因为$g1(i),g2(i)$是完全积性函数(只要是积性就可以了），所以他们的卷积也是积性函数

所以可以用线性筛做到$O(n)$ 我们类似于做约数和再维护一个h(i)表示不包括i的最小素因子的g1的和

于是就是多项式exp模板题了

#pragma GCC optimize(2)
#include <bits/stdc++.h>
using namespace std;
#define rint register int
#define IL inline
#define rep(i,h,t) for(int i=h;i<=t;i++)
#define dep(i,t,h) for(int i=t;i>=h;i--)
#define ll long long
#define me(x) memset(x,0,sizeof(x))
#define mep(x,y) memcpy(x,y,sizeof(y))
#define mid (t<=0?(h+t-1)/2:(h+t)/2)
namespace IO{
    char ss[1<<24],*A=ss,*B=ss;
    IL char gc()
    {
        return A==B&&(B=(A=ss)+fread(ss,1,1<<24,stdin),A==B)?EOF:*A++;
    }
    template<class T> void read(T &x)
    {
        rint f=1,c; while (c=gc(),c<48||c>57) if (c=='-') f=-1; x=(c^48);
        while (c=gc(),c>47&&c<58) x=(x<<3)+(x<<1)+(c^48); x*=f; 
    }
    char sr[1<<24],z[20]; ll Z,C1=-1;
    template<class T>void wer(T x)
    {
        if (x<0) sr[++C1]='-',x=-x;
        while (z[++Z]=x%10+48,x/=10);
        while (sr[++C1]=z[Z],--Z);
    }
    IL void wer1()
    {
        sr[++C1]=' ';
    }
    IL void wer2()
    {
        sr[++C1]='\n';
    }
    template<class T>IL void maxa(T &x,T y) {if (x<y) x=y;}
    template<class T>IL void mina(T &x,T y) {if (x>y) x=y;} 
    template<class T>IL T MAX(T x,T y){return x>y?x:y;}
    template<class T>IL T MIN(T x,T y){return x<y?x:y;}
};
using namespace IO;
const double ee=1.00000000000000;
const double pi=acos(-1.0);
const int N=4e5+10;
const int mo=1e9+7;
int r[N],n,m,l,inv[N],n1[N],n2[N];
struct cp{
    double a,b;
    IL cp operator +(const cp &o) const
    {
        return (cp){a+o.a,b+o.b};
    }
    IL cp operator -(const cp &o) const
    {
        return (cp){a-o.a,b-o.b};
    }
    IL cp operator *(register const cp &o) const
    {
        return (cp){a*o.a-b*o.b,o.a*b+o.b*a};
    }
}a[N],b[N],c[N],d[N];
IL int fsp(int x,int y)
{
    ll now=1;
    while (y)
    {
        if (y&1) now=now*x%mo;
        x=1ll*x*x%mo;
        y>>=1;
    }
    return now;
}
IL void clear()
{
    for (int i=0;i<=n;i++) a[i].a=a[i].b=b[i].a=b[i].b=c[i].a=c[i].b=d[i].a=d[i].b=0;
}
cp *w[N],tmp[N*2];
int p;
IL void init()
{
    cp *now=tmp;
    for (int i=1;i<=p;i<<=1)
    {
        w[i]=now;
        for (int j=0;j<i;j++) w[i][j]=(cp){cos(pi*j/i),sin(pi*j/i)};
        now+=i;
    }
}
IL void fft_init()
{
    l=0; for (n=1;n<=m;n<<=1) l++;
    for (int i=0;i<n;i++) r[i]=(r[i/2]/2)|((i&1)<<(l-1));
}
void fft(cp *a,int o)
{
    for (int i=0;i<n;i++) if (i>r[i]) swap(a[i],a[r[i]]);
    for (int i=1;i<n;i<<=1)
        for (int j=0;j<n;j+=(i*2))
        {
          cp *x1=a+j,*x2=a+i+j,*W=w[i];
          for (int k=0;k<i;k++,x1++,x2++,W++)
          {
              cp x=*x1,y=(cp){(*W).a,(*W).b*o}*(*x2); 
              *x1=x+y,*x2=x-y;
          }
        }
    if (o==-1) for(int i=0;i<n;i++) a[i].a/=n;
}
IL void getcj(int *A,int *B,int len)
{
    rep(i,0,len)
    {
        A[i]=(A[i]+mo)%mo,B[i]=(B[i]+mo)%mo;
    }
    for (int i=0;i<len;i++)
    {
       a[i]=(cp){A[i]&32767,A[i]>>15};
       b[i]=(cp){B[i]&32767,B[i]>>15};
    }
    m=len*2; fft_init();
    fft(a,1); fft(b,1);
    for (int i=0;i<n;i++)
    {
        int j=(n-1)&(n-i);
        c[j]=(cp){0.5*(a[i].a+a[j].a),0.5*(a[i].b-a[j].b)}*b[i];
        d[j]=(cp){0.5*(a[i].b+a[j].b),0.5*(a[j].a-a[i].a)}*b[i];
    }
    fft(c,1); fft(d,1);
    double inv=ee/n;
    rep(i,0,n) c[i].a*=inv,c[i].b*=inv;
    rep(i,0,n) d[i].a*=inv,d[i].b*=inv;
    rep(i,0,len)
    {
        ll a1=c[i].a+0.5,a2=c[i].b+0.5;
        ll a3=d[i].a+0.5,a4=d[i].b+0.5;
        B[i]=(a1+((a2+a3)%mo<<15)+((a4%mo)<<30))%mo;
    }
    clear();
}
void getinv(int *A,int *B,int len)
{
    if (len==1) { B[0]=fsp(A[0],mo-2); return;};
    getinv(A,B,(len+1)/2);
    int C[N]={};
    rep(i,0,len-1) C[i]=A[i];
    getcj(B,C,len);
    getcj(B,C,len);
    for (int i=0;i<len;i++) B[i]=((2ll*B[i]-C[i])%mo+mo)%mo;
}
IL void getDao(int *a,int *b,int len)
{
    for (int i=1;i<len;i++) b[i-1]=1ll*i*a[i]%mo;
    b[len-1]=0;
}
IL void getjf(int *a,int *b,int len)
{
    for (int i=0;i<len;i++) b[i+1]=1ll*a[i]*inv[i+1]%mo;
    b[0]=0;
}
IL void getln(int *A,int *B,int len)
{
    int C[N]={},D[N]={};
    getDao(A,C,len);
    getinv(A,D,len);
    getcj(C,D,len);
    getjf(D,B,len);
}
IL void getexp(int *A,int *B,int len)
{
    if (len==1) {B[0]=1; return;}
    getexp(A,B,(len+1)>>1);
    int C[N]={};
    getln(B,C,len);
    for(int i=0;i<len;i++) C[i]=(-C[i]+A[i])%mo;
    C[0]=(C[0]+1)%mo;
    getcj(C,B,len);
}
int main()
{
    freopen("1.in","r",stdin);
    freopen("1.out","w",stdout);
    inv[1]=1;
    rep(i,2,1e5+20) inv[i]=(1ll*inv[mo%i]*(mo-(mo/i)))%mo;
    int n,k;
    read(n); read(k); k++;
    n1[0]=1; n1[1]=-1;
    p=k<<1;
    init();
    getln(n1,n2,k);
    rep(i,0,k) n1[i]=-1ll*n*n2[i]%mo;
    rep(i,1,n)  
    {
      for(int j=1;j*i<=k;j++)
        (n1[i*j]-=inv[j])%=mo;
    }
    me(n2); 
    getexp(n1,n2,k);
    cout<<(n2[k-1]+mo)%mo<<endl;
    return 0;
}