\[枉过路，经行处，寒水凄凉漱黄土 \]\[相逢途，舟难住，雾凇封树冰刺骨 \]\[虽未闻望地如何，亦誓要从渃河渡 \]
首先，考虑可能的 \(b\) 序列，对于 \(b_i\)，设 \(p=\prod_{j\neq i}b_j\)，那么如果 \((b_i-1)p\ge k\)，则用 \(b_i-1\) 替换 \(b_i\) 一定更优。所以，我们使用 \(b_i=x\)，一定满足 \(x\) 是某个最小的 \(i\) 使得对于给定的 \(y\) 有 \(\lceil \dfrac{k}{i}\rceil=y\)。
并且，同时我们可以得出，从 \(b\) 中任意选出若干个 \(x\)，一定存在 \(i\) 使得 \(\prod x=\lceil \dfrac{k}{i}\rceil\)。
那么，我们就求出所有可能的 \(\lceil \dfrac{k}{i}\rceil\)。
for(int i=1,j;i<=k;i=j+1){
int h=(k+i-1)/i;
if(h==1)j=k;
else j=k/(h-1)-(!(k%(h-1)));
id[h]=++cnt,rs[cnt]=h;
}
如上代码即可以求出，我们每次求出的 \(i,j\) 即满足 \([i,j]\) 内所有整数就是满足 \(\lceil \dfrac{k}{x}\rceil=\lceil \dfrac{k}{i}\rceil\) 的整数集。
然后，我们设 \(dp_{i,j}\) 表示在前 \(i\) 个维度已经切完，在剩下的维度至少要切成 \(j\) 份的最大的 \(k\prod_{x\le i}\lceil\dfrac{a_i}{b_i}\rceil\dfrac{1}{a_i}\)。
而我们枚举从当前位置如何分 \(j\) 的时候，也使用和分 \(k\) 类似的方法——对于所有的 \(\lceil \dfrac{k}{i}\rceil\) 进行计算。因为 \(\lceil \dfrac{\lceil \dfrac{k}{i}\rceil}{j}\rceil\) 显然一定是 \(\lceil \dfrac{k}{i}\rceil\) 的形式，所以就可以直接枚举当前 \(j\) 的分因数进行转移。
所有可能的 \(\lceil \dfrac{k}{i}\rceil\) 一半分布在 \([1,\sqrt{k}]\)，一半分布在 \((\sqrt k,k]\)，所以可以证明总数是 \(O(\sqrt k)\) 的，也就是状态数量 \(O(n\sqrt k)\)。
然后，我们考虑 \(\lceil \dfrac{k}{i}\rceil\) 的分布，如果我们将其看作连续的则分开处理 \(\sqrt{k}\) 两边的两块，那么复杂度就是 \(\sum_{i\le \sqrt{k}}(\sqrt{\dfrac{k}{i}}+\sqrt{i})\)。考虑分开计算 \(\sum_{i\le \sqrt k}{\sqrt{i}}\) 和 \(\sum_{i\le \sqrt{k}}\sqrt{\dfrac{k}{i}}\)
对于 \(\sum_{i\ge \sqrt k}{\sqrt{i}}\)，考虑序列 \(a_i=i^{1/2}\)，其有限积分结果为 \(b_i=\dfrac{2}{3}i^{3/2}\)，则其前 \(\sqrt{k}\) 项和为 \(b_{\sqrt k}-b_0=\dfrac{2}{3}k^{3/4}\)
对于 \(\sum_{i\le \sqrt{k}}\sqrt{\dfrac{k}{i}}\) 很明显可以提出 \(\sqrt{k}\)，得到 \(\sqrt{k}\sum{\dfrac{1}{\sqrt{i}}}\)。求数列 \(a_i=i^{-1/2}\) 的前 \(k\) 项和，进行积分，得到 \(b_i=2i^{1/2}\)，在 \(\sqrt{k}\) 的取值为 \(2k^{1/4}\)，乘上外面的的 \(k^{1/2}\)，得到 \(2k^{3/4}\)。
所以，总的复杂度是 \(O(nk^{3/4})\) 的。
int n,k,a[105],id[20000005],rs[10005],cnt;
ld dp[105][10005];
signed main(){
ios::sync_with_stdio(false);
cin.tie(0);
cin>>n>>k;
rp(i,n)cin>>a[i];
for(int i=1,j;i<=k;i=j+1){
int h=(k+i-1)/i;
if(h==1)j=k;
else j=k/(h-1)-(!(k%(h-1)));
id[h]=++cnt,rs[cnt]=h;
}dp[0][id[k]]=k;ld ans=0;
for(int x=1;x<=n;x++){
for(int y=1;y<=cnt;y++)if(dp[x-1][y]!=0.0){
for(int i=1,j;i<=rs[y];i=j+1){
int h=(rs[y]+i-1)/i;
if(h==1)j=rs[y];
else j=rs[y]/(h-1)-(!(rs[y]%(h-1)));
if(i<=a[x])dp[x][id[h]]=max(dp[x][id[h]],dp[x-1][y]*((a[x]/i)*1.0/a[x]));
}
if(a[x]>=rs[y])ans=max(ans,dp[x-1][y]*(a[x]/rs[y]*1.0/a[x]));
}
}cout<<fixed<<setprecision(15)<<ans<<endl;
return 0;
}
//Crayan_r
还有别的不同实现，我们先将 \(k=k-1\)，考虑别的所有 \(b_i\) 是 \(\lfloor\dfrac{k}{i}\rfloor\)，然后只有一个 \(b_i\) 是 \(\lfloor\dfrac{k}{i}\rfloor+1\)，这样的结果可以满足 \(\prod b_i\ge k+1\)，也就是大于等于原来的 \(k\)，好处是有更方便的枚举实现和更小的常数，但是实际上和我们原来的做法只有这么一处区别而已。
\(\text{Code by sjcsjcsjc}\)
int n,k,a[110],rev[10010],id[10000010];
double dp[10010];
signed main()
{
ios::sync_with_stdio(false);cin.tie(0);
cin>>n>>k;k--;
int cnt=0;
for(int i=1,j;i<=k;i=j+1){
j=k/(k/i);
id[k/i]=++cnt;rev[cnt]=(k/i);
}
for(int i=1;i<=n;i++) cin>>a[i];
dp[id[k]]=(k+1)*1.0;
double ans=0.0;
for(int i=1;i<=n;i++){
int x=a[i];
for(int j=cnt;j>=1;j--) if(dp[j]!=0.0){
int y=rev[j];
for(int I=2,J;I<=y;I=J+1){
int u=y/I;
J=y/u;
int nw=id[u];
if(I<=x) dp[nw]=max(dp[nw],dp[j]*(x/I)*1.0/x);
}
if(y+1<=x) ans=max(ans,dp[j]*(x/(y+1))*1.0/x);
}
}
cout<<fixed<<setprecision(15)<<ans<<'\n';
return 0;
}
但是还有不同的做法。我们发现，在我们暴力计算的时候，一定会有一个位置满足前面所有项的积小于等于 \(\sqrt{k}\)，后面所有项的积也小于等于 \(\sqrt{k}\)。那么我们可以 \(\text{Meet in the middle}\)，通过 \(dp\) 暴力做出前缀和后缀不超过 \(\sqrt{k}\) 的结果。然后暴力枚举我们的中间项是哪一项，将左右合并得到答案。我们的 \(dp\) 第二维就只需要 \(\sqrt{k}\) 个状态，然后就可以对于所有 \(j\) 暴力枚举 \(1\sim \sqrt j\) 中的所有数进行转移。这种暴力几乎是最基础的，但是我们前面也证明了它的复杂度是 \(O(nk^{3/4})\)。
\(\text{Code by Um_nik}\)
const int N = 102;
const int M = 3233;
const int B = 301;
const int C = (int)1e7 / B + 3;
int n, k;
int a[N];
ld pref[N][M], suf[N][M];
int main()
{
startTime = clock();
scanf("%d%d", &n, &k);
int kk = k;
for (int i = 0; i < n; i++) {
scanf("%d", &a[i]);
kk = (kk + a[i] - 1) / a[i];
}
if (kk > 1) {
printf("0\n");
return 0;
}
if (k == 1) {
printf("1\n");
return 0;
}
pref[0][1] = 1;
suf[n][1] = 1;
for (int i = 0; i < n; i++)
for (int x = 1; x < M; x++)
for (int y = 1; y <= a[i] && x * y < 2 * M; y++)
pref[i + 1][min(M - 1, x * y)] = max(pref[i + 1][min(M - 1, x * y)], pref[i][x] * (ld)(a[i] / y) / (ld)a[i]);
for (int i = n - 1; i >= 0; i--)
for (int x = 1; x < M; x++)
for (int y = 1; y <= a[i] && x * y < 2 * M; y++)
suf[i][min(M - 1, x * y)] = max(suf[i][min(M - 1, x * y)], suf[i + 1][x] * (ld)(a[i] / y) / (ld)a[i]);
for (int i = 0; i <= n; i++)
for (int x = M - 2; x > 0; x--) {
pref[i][x] = max(pref[i][x], pref[i][x + 1]);
suf[i][x] = max(suf[i][x], suf[i][x + 1]);
}
ld ans = 0;
for (int i = 0; i < n; i++) {
for (int x = 1; x < M; x++)
for (int y = 1; y < M && x * y < C; y++) {
int z = (k + x * y - 1) / (x * y);
if (z > a[i]) continue;
ld cur = pref[i][x] * suf[i + 1][y];
cur *= (ld)(a[i] / z) / (ld)a[i];
ans = max(ans, cur);
}
for (int z = 1; z <= a[i] && z <= B; z++)
for (int x = 1; x < M; x++) {
int y = (k + x * z - 1) / (x * z);
if (y >= M) continue;
ld cur = pref[i][x] * suf[i + 1][y];
cur *= (ld)(a[i] / z) / (ld)a[i];
ans = max(ans, cur);
}
}
ans *= k;
printf("%.14lf\n", (double)ans);
return 0;
}