ORZ foreverlasting聚聚，QQ上问了他好久才会了这题（所以我们又聊了好久的Gal）
我们先来尝试推导一下(S)的性质，我们利用狄利克雷卷积来推：
[2^omega=Iast|mu|
]
这个很好理解吧，考虑一下它的组合意义即可
然后两边同卷上(I)有：
[2^omega ast I=Iast Iast |mu|=dast |mu|
]
后面还是同样，考虑(dast |mu|)的组合意义，一正一反的情况下其实就是(d(n^2))
因此我们有了(2^omegaast I=d(n^2))，即(S(n)=d(n^2))
那么显然(S)现在是个积性函数了，答案又是阶乘的形式，因此可以从(n-1)的答案推到(n)来
考虑一个非常暴力的过程，每次暴力分解质因数，复杂度大概是(O(nsqrt n))的
然后你只需要一台好一点的电脑我仿佛已经闻到了CPU的香气
然后考虑怎么优化这个过程，我们发现类似于某个套路，这种方法之所以慢是因为会出现不必要的枚举，因此我们只需要记录一下每个数的最小质因数，然后每次直接除去即可，顺带把贡献算一下
这样的复杂度很迷啊，加藤聚聚说是一个(log)的，我感觉还要再少点，毕竟向下除至少去掉一个(2)
那么我们就可以很快的做掉这道题了（用自己的笔记本跑了2s就出来了）
#include<cstdio>
#define RI register int
#define CI const int&
using namespace std;
const int N=10000000,mod=1000000087;
int prime[N+5],cnt,mnp[N+5],bkt[N+5],inv[(N<<1)+5],ret=1,ans;
#define Pi prime[j]
inline void init(void)
{
RI i,j; for (mnp[1]=1,i=2;i<=N;++i)
{
if (!mnp[i]) mnp[i]=i,prime[++cnt]=i;
for (j=1;j<=cnt&&1LL*i*Pi<=N;++j)
{
mnp[i*Pi]=Pi; if (i%Pi==0) break;
}
}
for (inv[0]=inv[1]=1,i=2;i<=(N<<1)+1;++i)
inv[i]=1LL*inv[mod%i]*(mod-mod/i)%mod;
}
#undef Pi
inline void inc(int& x,CI y)
{
if ((x+=y)>=mod) x-=mod;
}
inline int sum(CI x,CI y)
{
int t=x+y; return t>=mod?t-mod:t;
}
int main()
{
freopen("ans.txt","w",stdout);
init(); for (RI i=2;i<=N;++i)
{
ret=1LL*ret*inv[bkt[i]+1]%mod; inc(bkt[i],2); ret=1LL*ret*(bkt[i]+1)%mod;
for (int x=i;x!=mnp[x];x/=mnp[x])
{
if (x/mnp[x]==mnp[x])
{
ret=1LL*ret*inv[bkt[x]+1]%mod*inv[bkt[mnp[x]]+1]%mod;
inc(bkt[mnp[x]],sum(bkt[x],bkt[x]));
ret=1LL*ret*(bkt[mnp[x]]+1)%mod; bkt[x]=0;
} else
{
ret=1LL*ret*inv[bkt[x]+1]%mod*inv[bkt[mnp[x]]+1]%mod*inv[bkt[x/mnp[x]]+1]%mod;
inc(bkt[mnp[x]],bkt[x]); inc(bkt[x/mnp[x]],bkt[x]);
ret=1LL*ret*(bkt[mnp[x]]+1)%mod*(bkt[x/mnp[x]]+1)%mod; bkt[x]=0;
}
}
inc(ans,ret);
}
return printf("%d",ans),0;
}