三圆相切的解法难题

northwolves

只求和的情况下，可以提速一些：

1. f[n_]:=(c=0;Do[r=k/(1+k-2Sqrt[k]);If[IntegerQ[r3*r],c+=r3(1+k+r)],{r3,1,n/4},{k,Select[Range@n,IntegerQ[Sqrt[r3*#]]&&#>=4r3&]/r3}];c);f@1234

复制代码

王守恩

基本数量关系这样就可以。\(\D(r_1+r_2)^2-(r_1-r_2)^2=\big(\sqrt{(r_1+r_3)^2-(r_1-r_3)^2\ \ \ \ }+\sqrt{(r_2+r_3)^2-(r_2-r_3)^2\ \ \ \ }\ \big)^2\)

三个圆的半径都是整数时,三个半径之和——由上可得——Union@Flatten@Table[a (1 + b + b^2)^2 , {a, 9}, {b, 9}]——这个算式不好。还是请 northwolves 指正。

三个圆的半径都是整数时,三个半径之和是这样一串数——{9, 18, 27, 36, 45, 49, 54, 63, 72, 81, 98, 147, 169, 196, 245, 294, 338, 343, 392, 441, 507, 676, 845, 882, 961,
1014, 1183, 1323, 1352, 1521, 1764, 1849, 1922, 2205, 2646, 2883, 3087, 3249, 3528, 3698, 3844, 3969, 4805, 5329, 5547, 5766, 6498,——OEIS没有。

补充内容 (2025-3-10 10:33):
这样可能好一些()。Union@Flatten@Table[a (1 + b + b^2)^2, {b, 5}, {a, 625/(1 + b + b^2)^2}]

northwolves

观察比对了一会儿半径比值的规律，得到一个相对高效的代码

1. a[n_,r_]:=(b=SelectFirst[Reverse@Divisors@r,IntegerQ@Sqrt@#&,1];Range[2Sqrt@b,Sqrt@(n*b/r)]^2/b);
   
2. f[n_]:=(c=0;Do[r=k/(1+k-2 Sqrt[k]);If[IntegerQ[r3*r],c+=r3 (1+k+r)],{r3,1,n/4},{k,a[n,r3]}];c);
   
3. f[123456]

复制代码

4800044591

这样也可计算大一点的数字了：f[2^20]=346532068711

四来

mathe

根据关系式\(r_1r_2\)是完全平方数，可以设\(r_1=da^2,r_2=db^2\), 其中\((a,b)=1,a\le b\)
带入关系得到\(r_3=\frac{da^2b^2}{(a+b)^2}\),由此得到\((a+b)^2|d\),于是我们可以设\(d=(a+b)^2h\)
所以代入得到\(r_1=(a+b)^2a^2h,r_2=(a+b)^2b^2h,r_3=a^2b^2h\)
这说明\(r_2\)因子分解必然含平方因子。
题目说明的不清楚，假设题目是给定n时，求最大的S(n)
于是我们需要\((a+b)^2b^2h\le n\)时求\(((a+b)^2(a^2+b^2)+a^2b^2)h\)的最大值。

mathe

现在如果假设有一个解满足\(b\ge 1\), 我们选择\(a'=(a+b)b-1,b'=1\),h不变，那么同样满足约束条件，对应的S(n)值变为
于是\(a'+b'=(a+b)b\ge a+b, a'b'=ab+b^2-1\gt ab, a'^2+b'^2=(ab+b^2-1)^2+1\gt (2a-1+b^2)^2\ge (a+b^2)^2\ge(a+b)^2\gt a^2+b^2\).
所以S(n)必然变大，也就是我们总可以选择b=1
所以我们需要寻找满足
\((a+1)^2 h\le n\)时，求\( ((a+1)^2(a^2+1)+a^2)h\)的最大值。
记\(c=a+1\ge 2\), 我们要求\(c^2h\le n\)时求\((c^2-c+1)^2 h\)的最大值。
同样容易证明h中含有平方因子时，可以将平方因子转移到c里面结果会更大。所以最大情况h必然不含平方因子。

mathe

如果时求所有满足条件的情况之和，那么我们可以记
\(c=(a+b)b\ge 2\)
于是我们可以按c来分类
\(2\le c\le \sqrt{n}\), 在此条件下可以得到\(1\le h\le h_c=\lfloor \frac{n}{c^2}\rfloor\)
  然后我们可以对每个范围内的c将它表示为两个不同整数a+b,b的乘积，而且要满足\(b\ge a\),也就是\(b\lt a+b \le 2b\),也就是这两个因子相差不能太远
  对每个可能的a,b计算\((a+b)^2(a^2+b^2)+a^2b^2=(a^2+ab+b^2)^2\)并且求和。然后对和再乘上\(\frac{h_c(h_c+1)}2\)。
最后对于不同的c累加上面的结果即可。

northwolves

13楼代码改进版：

1. 
   
2. a[n_,r_]:=(b=Sqrt@SelectFirst[Reverse@Divisors@r,IntegerQ@Sqrt@#&,1];Range[2b,Sqrt@(n/r)*b]/b);
   
3. f[n_]:=(c=0;Do[r=k/(k-1);If[IntegerQ[r3*r^2],c+=r3 (1+k^2+r^2)],{r3,1,n/4},{k,a[n,r3]}];c);
   
4. f[123456]//Timing

复制代码

{0.328125, 4800044591}

northwolves

参照mathe版主15楼的思路，提速很可观：

1. f[n_]:=Total@Flatten@Table[{(a+b)^2*(a^2)*h,(a+b)^2(b^2)*h,a^2b^2*h},{a,1,n^(1/4)},{b,Select[Range@a,GCD[a,#]==1&]},{h,1,n/((a+b)a)^2}];f[1234567]//Timing

复制代码

{0.25, 480386841332}

northwolves

mathe版主一语道破天机

1. f[n_]:=Sum[h=Floor[n/((a+b)a)^2];1/2 (a^2+a b+b^2)^2 h (1+h),{a,1,n^(1/4)},{b,Select[Range@a,GCD[a,#]==1&]}];f[12345678901234]//Timing

复制代码

{1.04688, 480583650997741120228449478}