\(\mathcal{Description}\)

  Link

  （稍作简化：）对于变量 \(p_{1..n}\)，满足 \(p_i\in[0,1],~\sum p_i=1\) 时，求 \(\max \sum_{i=1}^n(p_i-p_i^2)i\)。

  数据组数 \(T\le10^5\)，\(n\le10^6\)。

\(\mathcal{Solution}\)

  Lagrange 乘子法的板题，可惜我不会。（

  先忽略 \(p_i\in[0,1]\) 的限制，发现这是一个带约数的最优化问题：

\[\max~~~~f=\sum_{i=1}^n(p_i-p_i^2)i,\\ \operatorname{s.t.}~~~~g=\sum_{i=1}^np_i-1=0. \]

考虑在无约束时，\(\forall i,~\frac{\partial f}{\partial p_i}=0\) 时能取到 \(f\) 的极值，本题由于偏导是一次式，仅有一个解，所以一定是最值。我们尝试将约束 \(g=0\) 变成 \(f\) 的一个维度，使得当 \(f'\) 关于 \(g\) 的偏导为 \(0\) 时恰有 \(g=0\)，就能化归为无约束的情况了。具体地，构造

\[L=f+\lambda g \]

此时 \(g\) 被作为引入变量 \(\lambda\) 的系数，所以当 \(\frac{\partial L}{\partial \lambda}=0\) 时，自然有 \(g=0\)。这就是 Lagrange 乘子法。

  回到本题，尝试直接解出 \(\lambda\)，先表示所有 \(p\)：

\[\frac{\partial L}{\partial p_i}=0=\lambda+2i-4ip_i\\ \Rightarrow p_i=\frac{\lambda+2i}{4i}. \]

代入约束 \(g=0\) 中：

\[\sum_{i=1}^n\frac{\lambda+2i}{4i}-1=0\\ \Rightarrow \lambda=\frac{4-2n}{h_n}. \]

其中 \(h_n\) 为调和级数前缀和。注意到直接代入 \(\lambda\) 可能时一段 \(p\) 的前缀 \(<0\)，所以只好二分钦定一个前缀为 \(0\)。假设钦定 \(p_{1..t-1}=0\)，类似地有

\[\lambda_t=\frac{4-2(n-t+1)}{h_n-h_{t-1}}. \]

代入 \(f_t\) 后简单化简有

\[f_t=\frac{(n+t)(n-t+1)}{4}-\frac{\lambda^2}{8(h_n-h_{t-1})}. \]

就能直接计算了。单次复杂度是二分的 \(\mathcal O(\log n)\)。可以通过单调滑动 \(t\) 值预处理做到 \(\mathcal O(n)-\mathcal O(1)\)。

\(\mathcal{Code}\)

/*~Rainybunny~*/

#include <cstdio>

#define rep( i, l, r ) for ( int i = l, rep##i = r; i <= rep##i; ++i )
#define per( i, r, l ) for ( int i = r, per##i = l; i >= per##i; --i )

const int MAXN = 1e6;
double h[MAXN + 5];

inline double calc( const int n, const int t ) {
    double c = n - t + 1, s = 0.5 * ( n + t ) * c,
      lam = ( 4 - 2 * c ) / ( h[n] - h[t - 1] );
    if ( ( lam + 2 * t ) / ( 4 * t ) < 0 ) return -1.;
    return 0.5 * s - 0.125 * lam * lam * ( h[n] - h[t - 1] );
}

int main() {
    freopen( "sample.in", "r", stdin );
    freopen( "sample.out", "w", stdout );

    rep ( i, 1, MAXN ) h[i] = h[i - 1] + 1. / i;

    int T, n; scanf( "%d", &T );
    while ( T-- ) {
        scanf( "%d", &n );
        
        int l = 1, r = n;
        while ( l < r ) {
            int mid = l + r >> 1;
            if ( calc( n, mid ) != -1. ) r = mid;
            else l = mid + 1;
        }
        printf( "%.12f\n", calc( n, l ) );
    }
    return 0;
}