题解:P3766 核心密码B

__Luna__

2024-11-18 21:05:34

Solution

思路

首先,我们给 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]); } ```