思路
首先,我们给 f(n) 做一个变形:
f(n)=\sum_{i=2}^n\frac{g(i)}{i}=\sum_{k=2}^{+\infin}\sum_{j=2}^{\lfloor\sqrt[k]n\rfloor}\frac{1}{j^k}
根据变形后的式子,我们可以分开处理不同幂次的完全 k 次方数。由于数据范围是 2\leq n\leq10^{18},而经计算有 2^{60}\approx1.1529215\times10^{18}>10^{18},故实际计算时指数 k 处理到 60 即可。
接下来进行对 \displaystyle\sum_{j=2}^{\lfloor\sqrt[k]n\rfloor}\frac{1}{j^k} 的计算,需要依据 k 分情况处理。
对于 k\ge3 的情况
由于底数 j 最大只有 \sqrt[3]{10^{18}}=10^6,可以暴力累加 \displaystyle\frac{1}{j^k}。对于多组询问可以使用一个技巧,即将询问全部读入,从小到大排序后依次处理。这样可以避免重复计算,也不需要记忆化占用大量内存。
对于 k=2 的情况
我们取一个不大的数 $N$,当 $j\le N$ 时,依然可以使用暴力。当 $j>N$ 时,我们可以将 $\displaystyle\frac{1}{j^2}$ 近似看成 $\displaystyle\frac{1}{j^2-0.25}$,裂项后即为 $\displaystyle\frac{1}{j-0.5}-\frac{1}{j+0.5}$。用裂项后的结果代替 $j>N$ 的各项,得到当 $j>N$ 时,有:
$$\sum_{j=2}^{\lfloor\sqrt n\rfloor}\frac{1}{j^2}\approx\sum_{j=2}^{N}\frac{1}{j^2}+\frac{1}{N+0.5}-\frac{1}{\lfloor\sqrt n\rfloor+0.5}$$
下面我们需要找到一个合适的 $N$。分析误差,可知经过近似后,每一项比准确值大了 $\displaystyle\frac{1}{j^2-0.25}-\frac{1}{j^2}=\frac{1}{4j^4-j^2}$。总误差为:
$$
\begin{align*}
\sum_{j=N+1}^{\lfloor\sqrt n\rfloor}\frac{1}{4j^4-j^2}&\le\sum_{j=N}^{+\infin}\frac{1}{4j^4-j^2}\\
&\approx\sum_{j=N}^{+\infin}\frac{1}{4j^4}\\
&\le\int_{N}^{+\infin}{4x^4}dx\\
&=\frac{1}{12N^3}
\end{align*}$$
根据题意,我们要求 $\displaystyle\frac{1}{12N^3}\le2\times10^{-14}$,解得 $N$ 最小为 $16092$。保险起见,我们可以取 $N=20000$。
## 代码
在编码细节上,需要注意精度问题。
[$283~\mathrm{ms}$](https://www.luogu.com.cn/record/189813395) 的 AC 代码:
```cpp
#include<bits/stdc++.h>
using namespace std;
const long double _1=1.0,m=_1/(2e4+0.5);
__int128 pow_(__int128 a,int x)
{
__int128 r=1;
while(x)
{
if(x&1)r*=a;
a*=a;
x>>=1;
}
return r;
}
struct nd
{
long long n;
int r;
nd(long long n=0,int r=0):n(n),r(r){}
bool operator<(nd x){return n<x.n;}
}a[50005];
long double ans[50005];
int main()
{
int T;
scanf("%d",&T);
for(int i=0;i<T;i++)
{
scanf("%lld",&a[i].n);
a[i].r=i;
}
sort(a,a+T);
for(int i=2,j,k;i<=60;i++)
{
k=2;
__int128 l=pow(k,i);
long double s=0;
for(j=0;j<T;j++) //暴力
{
if(i==2&&a[j].n>=4e8) //退出暴力
{
break;
}
for(;l<=a[j].n;k++,l=pow_(k,i))
{
s+=_1/l;
}
ans[a[j].r]+=s;
}
for(;j<T;j++) //近似处理
{
for(;k<=2e4;k++,l=pow(k,i))
{
s+=_1/l;
}
int n=sqrtl(a[j].n);
ans[a[j].r]+=s+(m-_1/(n+0.5));
}
}
for(int i=0;i<T;i++)printf("%.15llf\n",ans[i]);
}
```