题目描述
定义：
\[F(0)=1,F(n)=n \\ A_{i}^j=\frac{F(i)}{F(j)} \]求
\[\sum_{i=1}^n\sum_{j=1}^i\sum_{k=1}^{j}\gcd(A_{i-j}^j\times F(j),A_{j-k}^k\times F(k)) \]\(\mathcal{Solution}\)
式子看起来挺花里胡哨的，把里面的东西简单展开一下其实就是求
\[\sum_{i=1}^n\sum_{j=1}^i\sum_{k=1}^{j}\gcd((i-j)!,(j-k)!) \\ \sum_{i=1}^n\sum_{j=1}^i\sum_{k=1}^{j}\min(i-j,j-k)! \]遇到多个求和号在一起有很多突破的方向，这里选择考虑贡献。
这个 \(\min\) 似乎很难搞，但是终究是对若干个阶乘求和，则考虑对于一个 \(d\) ，计算 \(d!\) 的贡献，容易发现 \(d\) 在 \(\left[1,\left \lfloor \frac{n-1}{2} \right \rfloor\right ]\)内。
上面那个式子放在数轴上考虑，其实就是取出两个连在一起的区间 \([k,j],[j,i]\)，答案加上这两个区间较短区间的长度的阶乘。（这里的区间长度为右端点减左端点，即所求式子中的形式，后同）
把这两个区间合并成一个区间 \([k,i]\)，因为较短的区间长度为 \(d\)，则 \([k,i]\) 长度至少为 \(2d\)。
当 \([k,i]\) 的长度 \(len=2d\) 时没有区别，当 \(len>2d\) 时较短区间在左侧或右侧都可以。
那么就可以写出答案的式子：
\[\sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!\left(\sum_{len=2d+1}^{n}\sum_{i=1}^{n-len}2+\sum_{i=1}^{n-2d}1\right) \]\(d,len\) 均为上文定义，\(i\) 是枚举的大区间的左端点，也就是上文中 \([k,i]\) 的 \(k\)。
后面的括号的左半部分是 \(len>2d\) 时的贡献次数，右半部分是 \(len=2d\) 时的贡献次数。
化简一下：
\[\sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!\left(2\sum_{len=2d+1}^{n}(n-len)+(n-2d)\right) \\ \sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!\left(\sum_{len=0}^{n-2d}len+(n-2d)\right) \\ \sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!\left(2\times\frac{(n-2d)(n-2d-1)}{2}+(n-2d)\right) \\ \sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!\left(n-2d\right)^2 \]最后这个东西长的十分好看，事实上也很好求，把里面的完全平方式展开一下：
\[n^2\sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!-n\sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!d+\sum_{d=1}^{\left \lfloor \frac{n-1}{2} \right \rfloor}d!4d^2 \]发现里面三个求和号都能在 \(\mathcal{O}(n)\) 的复杂度内预处理，预处理后直接计算即可。
时间复杂度 \(\mathcal{O}(\max{n}+T)\)
#include<iostream>
#include<cstdio>
#define ll long long
const ll mod = 10086001;
template <typename T> T Max(T x, T y) { return x > y ? x : y; }
template <typename T> T Min(T x, T y) { return x < y ? x : y; }
template <typename T> T Add(T x, T y) { return (x + y > mod) ? (x + y - mod) : (x + y); }
template <typename T> T Mod(T x) { return (x >= mod) ? (x - mod) : x; }
template <typename T>
T& read(T& r) {
r = 0; bool w = 0; char ch = getchar();
while(ch < '0' || ch > '9') w = ch == '-' ? 1 : 0, ch = getchar();
while(ch >= '0' && ch <= '9') r = r * 10 + (ch ^ 48), ch = getchar();
return r = w ? -r : r;
}
const int N = 1000100;
ll fac[N], f1[N], f2[N], f3[N];
int T, n[N], mx;
signed main() {
read(T);
for(int i = 1; i <= T; ++i) mx = Max(mx, read(n[i]));
fac[0] = f1[0] = 1;
for(int i = 1; i <= mx; ++i)
fac[i] = fac[i-1] * i % mod,
f1[i] = Add(f1[i-1], fac[i]),
f2[i] = Add(f2[i-1], fac[i] * 4 * i % mod),
f3[i] = Add(f3[i-1], fac[i] * 4 * i % mod * i % mod);
for(int i = 1; i <= T; ++i) {
int m = (n[i]-1)/2;
printf("%lld\n", Add(Mod(1ll * n[i] * n[i] % mod * f1[m] % mod - n[i] * f2[m] % mod + mod), f3[m]));
}
return 0;
}