定义
约定1：以下分数都是最简，且令$frac{1}{0}$有意义，其大于其余分数，并称平行于$y$轴的直线斜率为$-frac{1}{0}$
分数加：对于分数$a=frac{a_{1}}{a_{0}}$、$b=frac{b_{1}}{b_{0}}$，定义$aoplus b=frac{a_{1}+b_{1}}{a_{0}+b_{0}}$
Farey neighbor：对于两个非负最简分数$a=frac{a_{1}}{a_{0}}$和$b=frac{b_{1}}{b_{0}}$，其为Farey neighbor当且仅当$a_{0}b_{1}-a_{1}b_{0}=1$
约定2：二叉树根深度为0
Stern-Brocot Tree：这是一棵满二叉树，其每一个节点的元素是一个分数
其中根的分数为$frac{1}{1}$，第$h$层（$hge 1$）的分数通过如下方式构造：
对于前$h-1$层，其中序遍历依次为$a_{1},a_{2},...,a_{2^{h}-1}$，再令$a_{0}=frac{0}{1}$，$a_{2^{h}}=frac{1}{0}$，则第$h$层的元素从左到右依次为$a_{0}oplus a_{1},a_{1}oplus a_{2},...,a_{2^{h}-1}oplus a_{2^{h}}$
约定3：$x_{A}$和$y_{A}$分别表示点$A$的$x$和$y$坐标，$A$为整点当且仅当$x_{A},y_{A}in Z$
约定4："图形内"不包含边界，严格在某一条线上指不包含该线端点的部分
题解
题目即求$sum_{i=1}^{N}lfloorfrac{N}{i} floor$，这个用数论分块做可以做到$o(Tsqrt{N})$，但还是无法通过
考虑这个式子的几何意义：对于$xy=N$在第一象限内的曲线$L$，其与$x$轴和$y$轴所构成的图形内以及严格在$L$上的整点数
定义$S(P,a,b)$：令$A$和$B$分别表示$L$上导数为$-a$和$-b$的两点，$t_{A}$为$L$过点$A$的切线，则$S(P,a,b)$为$t_{A}$、$t_{B}$与$L$所构成的图形内部以及严格在$L$的$AB$上的整点数，并令$T$为$t_{A}$与$t_{B}$的交点
（虽然与$P$无关，但没有关系）
设置$S(P,a,b)$的定义域为：
1.$a$和$b$为Farey neighbor
2.$P$是整点，过点$P$作两条斜率为$-a$和$-b$的直线，记作$p_{a}$和$p_{b}$，称$p_{a}$对应的直线为$t_{A}$，$p_{b}$对应的直线为$t_{B}$，则满足——
(1)$p_{a},p_{b}$在其对应的直线下方（可以重合）
(2)这四条直线所构成的无限大的图形内没有整点
(3)$t_{A},t_{B}$分别满足：与对应其的直线重合或其上没有整点
初始即求$S(P(0,0),frac{0}{1},frac{1}{0})$，符合定义域
$S(P,a,b)$的计算
考虑求$S(P,a,b)$，设切点别为$A$和$B$，$a=frac{a_{1}}{a_{0}}$，$b=frac{b_{1}}{b_{0}}$，令$c=aoplus b$，选取点$C$使得$x_{B}<x_{C}<x_{A}$且$f'(x_{C})=-c$，由于$f'(x)$（在$x>0$）是一个单调递增且连续的函数，因此显然存在
设$t_{C}$交$t_{A}$和$t_{B}$分别于$P_{1}$和$P_{2}$，那么所求的问题也就分为了三个部分$S(P',a,c)$、$S(P',c,b)$和在$Delta TP_{1}P_{2}$内部以及严格在$P_{1}P_{2}$上的整点数，递归即可
但现在还有三个问题：
1.第三部分的计算
2.$P'$的选择以及递归结束条件
3.时间复杂度
问题1
首先，根据$P$的第2和3个条件，可以看作与$P$相交，具体来说再令$P_{1}$和$P_{2}$分别为$t_{C}$与$p_{a}$和$p_{b}$的交点，那么所求的也就是$Delta PP_{1}P_{2}$以及严格在$P_{1}P_{2}$上的整点数
对于上面的问题，考虑构建一个坐标系：
1.$P$为原点，$PP_{1}$和$PP_{2}$分别为$x$轴和$y$轴（同时方向也按照这个，即$P_{1}$和$P_{2}$都在$x$和$y$轴正半轴上）
2.令$vec{i}=(a_{0},-a_{1})$，$vec{j}=(-b_{0},b_{1})$，则单位长度即为两向量模长
不难发现，$x$轴方向上单位长度的向量即为$vec{i}$，类似地$y$轴上是$j$，那么对于$(x,y)$，其所对应的新坐标系中的点$(u,v)$即满足$uvec{i}+vvec{j}+P=(x,y)$
由于是恰好是两个方程两个未知数，且显然不矛盾或等价，因此两坐标系中的点一一对应
接下来，我们还要证明两坐标系中的整点一一对应，换言之若$(u,v)$是整点则$(x,y)$是整点，若$(x,y)$是整点则$(u,v)$也是整点
设$P(x_{0},y_{0})$，由于$P$是整点，$x_{0},y_{0}$也是整数，那么前者显然成立
后者考虑暴力解出方程组，即
$$
egin{cases}ua_{0}-vb_{0}+x_{0}=x\-ua_{1}+vb_{1}+y_{0}=yend{cases}iffegin{cases}u=frac{b_{1}(x-x_{0})+b_{0}(y-y_{0})}{a_{0}b_{1}-a_{1}b_{0}}\v=frac{a_{1}(x-x_{0})+a_{0}(y-y_{0})}{a_{0}b_{1}-a_{1}b_{0}}end{cases}
$$
由于$a$和$b$是Farey neighbor，根据定义即$a_{0}b_{1}-a_{1}b_{0}=1$，因此$u$和$v$也是整数
由于整点一一对应，我们将$Delta PP_{1}P_{2}$放到这个新的坐标系中去，考虑$P_{1}P_{2}$斜率：
$$
-ua_{1}+vb_{1}+y_{0}=-frac{a_{1}+b_{1}}{a_{0}+b_{0}}(ua_{0}-vb_{0}+x_{0})+b
$$
简单化简，即
$$
v=frac{(a_{0}+b_{0})a_{1}-(a_{1}+b_{1})a_{0}}{(a_{0}+b_{0})b_{1}-(a_{1}+b_{1})b_{0}}u+b'
$$
展开后可以发现，这个斜率恰好为-1，同时$P_{1}$与$x$轴有交点且在$P_{1}P_{2}$上，因此在新坐标系中有$P_{1}(b',0)$，类似地$P_{2}(0,b')$，那么不难求出整点数为$frac{lfloor b'
floor(lfloor b'
floor-1)}{2}$
具体式子：$C(sqrt{frac{N}{c}},sqrt{Nc})$，$b=ccdot x_{C}+y_{C}=2y_{C}$，$b'=(a_{0}+b_{0})(b-y_{0})-(a_{1}+b_{1})x_{0}$
问题2
选择$P'(lfloor b' floor,0)$和$P'(0,lfloor b' floor)$即可，正确性显然
边界条件需要考虑曲线$L$在上面所构造的新坐标系中的一些性质，首先就是对应的方程：
$$
(ua_{0}-vb_{0}+x_{0})(-ua_{1}+vb_{1}+y_{0})=N
$$
构造$w_{0}$和$w_{1}$，使得
$$
egin{cases}x_{0}=a_{0}w_{0}-b_{0}w_{1}\y_{0}=b_{1}w_{1}-a_{1}w_{0}end{cases}iffegin{cases}w_{0}=frac{b_{1}x_{0}+b_{0}y_{0}}{a_{0}b_{1}-a_{1}b_{0}}\w_{1}=frac{a_{1}x_{0}+a_{0}y_{0}}{a_{0}b_{1}-a_{1}b_{0}}end{cases}
$$
由于$a_{0}b_{1}-a_{1}b_{0}=1$，即
$$
egin{cases}w_{0}=b_{1}x_{0}+b_{0}y_{0}\w_{1}=a_{1}x_{0}+a_{0}y_{0}end{cases}
$$
接下来，代入即可得到
$$
(a_{0}(u+w_{0})-b_{0}(v+w_{1}))(b_{1}(v+w_{1})-a_{1}(u+w_{0}))=N
$$
化简后，即
$$
(a_{0}b_{1}+a_{1}b_{0})(u+w_{0})(v+w_{1})-a_{0}a_{1}(u+w_{0})^{2}-b_{0}b_{1}(v+w_{1})^{2}=N
$$
这个形式已经比较简洁了，直接暴力使用求根公式，即
$$
v=frac{(a_{0}b_{1}+a_{1}b_{0})(u+w_{0})pmsqrt{(a_{0}b_{1}+a_{1}b_{0})^{2}(u+w_{0})^{2}-4b_{0}b_{1}(a_{0}a_{1}(u+w_{0})^{2}-N)}}{2b_{0}b_{1}}-w_{1}
$$
对后面的根号中提出$(u+w_{0})^{2}$的系数，恰好是$(a_{0}b_{1}-a_{1}b_{0})^{2}=1$，即
$$
v=frac{(a_{0}b_{1}+a_{1}b_{0})(u+w_{0})pm sqrt{(u+w_{0})^{2}-4b_{0}b_{1}N}}{2b_{0}b_{1}}-w_{1}
$$
对于正负号，联系图像不难发现另一个解是在$L$上$B$的左侧，其的$v$较大，因此应该取负，即
$$
V(u)=frac{(a_{0}b_{1}+a_{1}b_{0})(u+w_{0})-sqrt{(u+w_{0})^{2}-4b_{0}b_{1}N}}{2b_{0}b_{1}}-w_{1}
$$
对其求导，即
$$
V'(u)=frac{a_{0}b_{1}+a_{1}b_{0}-frac{u+w_{0}}{sqrt{(u+w_{0})^{2}-4b_{0}b_{1}N}}}{2b_{0}b_{1}}
$$
考虑后式，即
$$
frac{u+w_{0}}{sqrt{(u+w_{0})^{2}-4b_{0}b_{1}N}}=sqrt{1+frac{4b_{0}b_{1}N}{(u+w_{0})^{2}-4b_{0}b_{1}N}}
$$
因此$V(u)$的导数单调递增，换言之其下凸
进一步的，考虑$A$所对应的点必然是$v$最小的点，根据下凸那么$A$所对应的$u$之前的点都是单调递减的，这也就是我们所关心的范围
这个范围中，由于$x$轴和$y$轴是不算的，因此若$u=1$时$V'(u)>0$或$V(u)<1$，结果即为0
在这个基础上，我们再设置一个边界条件：
若$a=frac{0}{1}$或$b=frac{1}{0}$（以$a=frac{0}{1}$为例），当$b_{0}>N^{frac{1}{3}}$，暴力枚举$yin [1,lfloor y_{B} floor]$进行统计
问题3
下面我们进行复杂度分析，同样按照上面的边界条件分为两类
1.当$a=frac{0}{1}$或$b=frac{1}{0}$（以$a=frac{0}{1}$为例），$b_{0}$是不断增加1的，且由于$b_{0}>N^{frac{1}{3}}$时即结束，那么状态只有$o(N^{frac{1}{3}})$个，同时枚举$y$的时间复杂度为$y_{B}=frac{N}{x_{B}}=sqrt{bN}=sqrt{frac{N}{b_{0}}}$，计算复杂度也是$o(N^{frac{1}{3}})$
2.当$a e frac{0}{1}$且$b e frac{1}{0}$，将所有状态，再分为两类：
(1)直接返回的状态，假设有$N_{0}$个
(2)递归下去的状态，假设有$N_{2}$个
考虑这棵搜索树，不难证明$N_{0}=N_{2}+1$，即我们的复杂度是$o(N_{0}+N_{2})=o(N_{2})$
对于$N_{2}$，由于每一个$(a,b)$最多被搜索一次，也就是有多少对$(a,b)$是第二类
考虑$(a,b)$是第二类状态的条件，必要条件为新坐标系中$(u,v)=(0,1)或(1,0)$都在新坐标系中$V(u)$与$x$轴、$y$轴内部
这两个点所对应的点的$x$坐标即为$x_{0}-b_{0}$和$x_{0}+a_{0}$，其恰好在$PP_{1}$或$PP_{2}$上，该图形的$x$坐标范围为$[x_{P_{2}},x_{P_{1}}]$，因此即有
$$
x_{P_{2}}le x_{0}-b_{0}<x_{0}+a_{0}le x_{P_{1}}
$$
根据$P$的第二个性质，有$x_{P_{1}}le x_{A}+1$以及$x_{B}-1le x_{P_{2}}$，否则必然存在整点，再根据$x_{A}=sqrt{frac{N}{a}}$以及$x_{B}=sqrt{frac{N}{b}}$，代入即
$$
sqrt{frac{N}{b}}-1le x_{0}-b_{0}<x_{0}+a_{0}le sqrt{frac{N}{a}}+1
$$
由于不确定$x_{0}$，不能直接判断，但可以推出$sqrt{frac{N}{a}}-sqrt{frac{N}{b}}+2ge a_{0}+b_{0}$，状态数即为
$$
sum_{a_{0}b_{1}-a_{1}b_{0}=1}[sqrt{frac{N}{a}}-sqrt{frac{N}{b}}+2ge a_{0}+b_{0}]
\=sum_{a_{0}b_{1}-a_{1}b_{0}=1}[frac{sqrt{a_{0}b_{1}}-sqrt{a_{1}b_{0}}}{sqrt{a_{1}b_{1}}}ge frac{a_{0}+b_{0}-2}{sqrt{N}}]
$$
注意到$a_{0}b_{1}=a_{1}b_{0}+1$，有$sqrt{w+1}-sqrt{w}=frac{1}{sqrt{w}+sqrt{w+1}}$，即
$$
sum_{a_{0}b_{1}-a_{1}b_{0}=1}[frac{sqrt{a_{0}b_{1}}-sqrt{a_{1}b_{0}}}{sqrt{a_{1}b_{1}}}ge frac{a_{0}+b_{0}-2}{sqrt{N}}]\le sum_{a_{0}b_{1}-a_{1}b_{0}=1}[frac{1}{sqrt{a_{1}^{2}b_{0}b_{1}}}ge frac{a_{0}+b_{0}-2}{sqrt{N}}]\=sum_{a_{0}b_{1}-a_{1}b_{0}=1}[a_{1}^{2}b_{0}b_{1}(a_{0}+b_{0}-2)^{2}le N]
$$
由于$a_{0},b_{0}ge 1$，因此$(a_{0}+b_{0}-2)^{2}ge 2a_{0}b_{0}$，即
$$
sum_{a_{0}b_{1}-a_{1}b_{0}=1}[a_{1}^{2}b_{0}b_{1}(a_{0}+b_{0}-2)^{2}le N]\le sum_{a_{0}b_{1}-a_{1}b_{0}=1}[2a_{0}a_{1}^{2}b_{0}^{2}frac{1+a_{1}b_{0}}{a_{0}}le N]\le sum_{a_{0}b_{1}-a_{1}b_{0}=1}[2a_{1}^{3}b_{0}^{3}le N]\le sum_{i=1}^{lfloor N^{frac{1}{3}}
floor}sigma(i)cdot sigma(i+1)
$$
又注意到有$xyle frac{x^{2}+y^{2}}{2}$，即
$$
sum_{i=1}^{lfloor N^{frac{1}{3}}
floor}sigma(i)cdot sigma(i+1)le sum_{i=1}^{lfloor N^{frac{1}{3}}
floor}sigma(i)^{2}
$$
后面的这个东西有结论，$sum_{i=1}^{n}sigma(i)^{2}$是$o(nlog^{3}n)$级别的，因此总复杂度即$o(N^{frac{1}{3}}log^{3}N)$
似乎还可以更精确地证明其是$o(N^{frac{1}{3}}log N)$的
代码实现
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define ll long long 4 #define lll __int128 5 #define eps 1e-7 6 int t,a[105]; 7 ll n; 8 lll ans; 9 bool check(lll k){ 10 return k>(lll)1<<63; 11 } 12 lll sqr(lll k){ 13 return k*k; 14 } 15 void write(lll ans){ 16 a[0]=0; 17 while (ans){ 18 a[++a[0]]=ans%10; 19 ans/=10; 20 } 21 if (!a[0])printf("0"); 22 for(int i=a[0];i;i--)printf("%d",a[i]); 23 printf(" "); 24 } 25 void dfs(lll x0,lll y0,int a1,int a0,int b1,int b0){ 26 if ((!a1)&&(a0==1)){ 27 if ((lll)b0*b0*b0>n){ 28 ll b=floor(2*sqrt(n)*sqrt(b0)+eps); 29 for(int i=1;(lll)i*i*b0<=n;i++)ans+=n/i-(b-1LL*i*b0); 30 return; 31 } 32 } 33 else{ 34 ll s=(ll)a0*b1+(ll)a1*b0; 35 lll w0=b1*x0+b0*y0+1,w1=a1*x0+a0*y0+1; 36 if (check(w0))return; 37 lll D=sqr(w0)-(lll)4*b0*b1*n; 38 if ((D<0)||(D>sqr(w0)/sqr(s))||(s*w0<(lll)2*b0*b1*w1))return; 39 if ((!check(s*w0-(lll)2*b0*b1*w1))&&(sqr(s*w0-(lll)2*b0*b1*w1)<D))return; 40 } 41 int c1=a1+b1,c0=a0+b0; 42 double c=1.0*c1/c0,b=2*sqrt(n)*sqrt(c); 43 lll bb=floor((a0+b0)*b+eps)-(a0+b0)*y0-(a1+b1)*x0; 44 ans+=bb*(bb-1)/2; 45 dfs(bb*a0+x0,-bb*a1+y0,a1,a0,c1,c0); 46 if ((b1==1)&&(b0==0))ans=ans*2-bb*(bb-1)/2; 47 else dfs(-bb*b0+x0,bb*b1+y0,c1,c0,b1,b0); 48 } 49 int main(){ 50 scanf("%d",&t); 51 while (t--){ 52 scanf("%lld",&n); 53 ans=0; 54 dfs(0,0,0,1,1,0); 55 write(ans); 56 } 57 }
这份代码具有较大的浮点误差，只能保证在$Nle 10^{14}$的正确性，且递归常数较大
优化
考虑上面这个过程中所有结束状态的$p_{a}$和$p_{b}$，若将其以此连在一起，其在曲线$L$的下方（可以有重合的点），且之间没有任何整点
同时，如果能够确定这条折线，由于所有点都是整点，根据皮克定理，不难确定这个折线下方的整点数，也就是所求的答案
事实上，这也就是上面这个算法的本质，我们不妨考虑直接构造这条直线
另外，由于$L$是下凸的，在$L$上方两点连线一定在$L$的外部，而内部没有这个性质，因此考虑构造在$L$外部的折线
具体实现可以参考代码，代码之后会再说明一下
（代码大概就是抄了一下论文中所给的代码）
1 #include<bits/stdc++.h> 2 using namespace std; 3 #define ll long long 4 #define lll __int128 5 struct frac{ 6 int x,y; 7 frac operator +(const frac &a)const{ 8 return frac{x+a.x,y+a.y}; 9 } 10 }; 11 int t,a[105]; 12 ll n; 13 lll x,y,ans; 14 stack<frac>st; 15 bool check(lll x,lll y){ 16 return x*y>n; 17 } 18 void write(lll ans){ 19 a[0]=0; 20 while (ans){ 21 a[++a[0]]=ans%10; 22 ans/=10; 23 } 24 if (!a[0])printf("0"); 25 for(int i=a[0];i;i--)printf("%d",a[i]); 26 printf(" "); 27 } 28 int main(){ 29 scanf("%d",&t); 30 while (t--){ 31 scanf("%lld",&n); 32 ll K=floor(sqrt(n)); 33 ans=0; 34 x=n/K,y=n/x+1; 35 st.push(frac{1,0}); 36 st.push(frac{1,1}); 37 while (!st.empty()){ 38 frac k=st.top(); 39 st.pop(); 40 while (check(x+k.x,y-k.y)){ 41 ans+=x*k.y+(k.y+1)*(k.x-1)/2; 42 x+=k.x,y-=k.y; 43 } 44 if (y*y*y<=n)break; 45 frac o; 46 while (1){ 47 o=k,k=st.top(); 48 if (check(x+k.x,y-k.y))break; 49 st.pop(); 50 } 51 while (1){ 52 frac s=k+o; 53 if (check(x+s.x,y-s.y))st.push(k=s); 54 else{ 55 if ((x+s.x)*(x+s.x)*k.y>=(lll)n*k.x)break; 56 o=s; 57 } 58 } 59 } 60 for(int i=1;i<y;i++)ans+=n/i; 61 ans=ans*2-K*K; 62 write(ans); 63 } 64 }
这个过程其实与前面的过程是一致的，实际上可以看作是一个将dfs的栈手动实现的过程
从dfs的角度考虑，栈中相邻两个元素即为之前的一次dfs中的两个分数
之后相当于不断递归，其中有一个比较特殊的地方就是关于
if ((x+s.x)*(x+s.x)*k.y>=(lll)n*k.x)break;
这里判定的其实是$-f'(x+s.x)ge frac{k.y}{k.x}$，这个判定是类似于前面判定$U(1)ge 1$的，也就是说没有整点就结束
关于没有整点，其实就等价于没有斜率（绝对值）比我更大的直线，否则那一个点就在其下方且是整点，因此要求$-f'(x+s.x)ge frac{k.y}{k.x}$
另外还有一些细节问题：
1.从$K=lfloorsqrt{N} floor$，因为根据对称性有$sum_{i=1}^{N}lfloorfrac{N}{i} floor=2sum_{i=1}^{K}lfloorfrac{N}{i} floor-K^{2}$
2.做第一个优化可以使得我们只需要判定$yle N^{frac{1}{3}}$的时候（不需要判定$xle N^{frac{1}{3}}$），根据前面这里要暴力处理
这样的时间复杂度与之前的复杂度相同，即为$o(N^{frac{1}{3}}log^{3}N)$，且没有浮点误差、常数更小