神仙题。
首先我们每次肯定会按 \(1\sim k\),\(k+1\sim 2k\),\(2k+1\sim 3k\) 的顺序选择勘测的洞穴,显然如果出现两个洞穴在某一时刻被勘测次数差 \(\ge 2\) 那肯定是不优的,这个感性理解一下即可。
那么我们假设 \(E_i\) 为在 \(i\) 处的宝藏被勘测到的期望天数,那么由于我们按 \(1\sim k,k+1\sim 2k,2k+1\sim 3k\) 的顺序勘测,有 \(E_{i+k}=E_{i}+1(i+k\le n)\),以及 \(E_{i+k-n}=E_{i}·\dfrac{1}{2}+1(i+k>n)\),最终答案就是 \(\dfrac{1}{n}·\sum\limits_{i=1}^nE_i\)。使用手动高斯消元可以做到 \(\mathcal O(n)\) 的复杂度,但还是不足以解决这道题。
我们考虑先手玩一些 \(n,k\) 比较小的情况找找规律,譬如 \(n=11,k=4\),那么我们将 \(\bmod 4\) 相同的归为一组,有 \(E_4+1=E_8\),\(E_1+2=E_5+1=E_9\),\(E_2+2=E_6+1=E_{10}\),\(E_3+2=E_7+1=E_{11}\),我们将 \(ans\) 表示成 \(E_1\sim E_4\) 的线性组合,可以得到 \(ans=\dfrac{1}{11}·(3E_1+3+3E_2+3+3E_3+3+2E_4+1)\),这样我们只要能求出 \(3\sum\limits_{i=1}^3E_i+2E_4\) 即可算出答案,而类似地我们也可以得到恒等式 \(E_1=\dfrac{1}{2}E_8+1=\dfrac{1}{2}(E_4+1)+1\),\(E_2=\dfrac{1}{2}E_9+1=\dfrac{1}{2}(E_1+2)+1\),这样我们建立了恒等式 \(E_1=\dfrac{1}{2}E_4+\dfrac{3}{2},E_2=\dfrac{1}{2}E_1+2,E_3=\dfrac{1}{2}E_2+2\)。问题就转化为 \(n=4,k=3\) 的情况。
受到这个问题的启发,我们考虑递归地解决问题,solve(N, K, A, B, f1, f2)
表示目前问题规模 \(n=N,k=K\),\(A,B,f1,f2\) 为两个一次函数,表示 \(E_{i+k}=A(E_i)(i+k\le n)\),\(E_{i+k-n}=B(E_i)\),\(ans=\sum\limits_{i=1}^kf1(E_i)+\sum\limits_{i=k+1}^nf2(E_i)\)。递归求解即可,具体实现见代码。
时间复杂度 \(T\log^2n\),1 个 log 是欧几里得的复杂度,1 个 log 是复合逆的快速幂。
const int MOD = 1e9 + 7;
const int INV2 = MOD + 1 >> 1;
int qpow(int x, int e) {
int ret = 1;
for (; e; e >>= 1, x = 1ll * x * x % MOD)
if (e & 1) ret = 1ll * ret * x % MOD;
return ret;
}
struct line {
int k, b;
line() {k = b = 0;}
line(int _k, int _b): k(_k), b(_b) {}
int at(int x) {return (1ll * x * k + b) % MOD;}
line operator + (const line &rhs) {return line((k + rhs.k) % MOD, (b + rhs.b) % MOD);}
line operator - (const line &rhs) {return line((k - rhs.k + MOD) % MOD, (b - rhs.b + MOD) % MOD);}
void print() {printf("y = %dx + %d\n", k, b);}
};
line merge(line F, line G) {//G(F(x))
return line(1ll * F.k * G.k % MOD, (1ll * G.k * F.b + G.b) % MOD);
}
line getinv(line F) {
int iv = qpow(F.k, MOD - 2);
return line(iv, (MOD - 1ll * iv * F.b % MOD) % MOD);
}
line qpow_l(line x, int e) {
line res = line(1, 0);
for (; e; e >>= 1, x = merge(x, x)) if (e & 1) res = merge(res, x);
return res;
}
line getsum(line x, int e) {
if (!e) return line(0, 0);
if (x.k == 1) return line(e, 1ll * e * (e + 1) / 2 % MOD * x.b % MOD);
line pw = qpow_l(x, e + 1) - x;
line res; res.k = 1ll * pw.k * qpow(x.k - 1, MOD - 2) % MOD;
res.b = 1ll * x.b * (res.k - e + MOD) % MOD * qpow(x.k - 1, MOD - 2) % MOD;
return res;
}
int solve(int n, int k, line A, line B, line F1, line F2) {
if (!k) {
int E1 = 1ll * A.b * qpow((1 - A.k + MOD) % MOD, MOD - 2) % MOD;
return F2.at(E1);
}
line f1 = F1 + merge(getsum(A, n / k), F2);
line f2 = F1 + merge(getsum(A, n / k - 1), F2);
line a = merge(getinv(B), qpow_l(getinv(A), n / k - 1));
line b = merge(getinv(B), qpow_l(getinv(A), n / k));
return (solve(k, n % k, a, b, line(f1.k, 0), line(f2.k, 0)) + 1ll * (n % k) * f1.b + 1ll * (k - n % k) * f2.b) % MOD;
}
int main() {
freopen("silicon20.in", "r", stdin);
freopen("silicon20.out", "w", stdout);
int qu; scanf("%d", &qu);
while (qu--) {
int n, k; scanf("%d%d", &n, &k); int d = __gcd(n, k); n /= d; k /= d;
printf("%d\n", 1ll * solve(n, k, line(1, 1), line(INV2, 1), line(1, 0), line(1, 0)) * qpow(n, MOD - 2) % MOD);
}
return 0;
}