三圆相切的解法难题

northwolves

northwolves 发表于 2025-3-8 18:43
mathe版主一语道破天机

1/2提到括号之外还能快一倍：

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

复制代码

{0.5625, 48058365099774112022844947}

四来

根据以上公式，嵌套循环求解, 取 N=2000

简单暴力粗犷的 Java 代码:

1. private void printSn(int N) {
   
2.         int[] S = new int[N + 2]; //结果保存数组
   
3.         int fz, fm;  //对应公式右端的分子，分母
   
4.         int r1, r2, r3;  //三个半径
   
5.         int root; //r1r2 的平方根
   
6.         int n;  //S(n) 的 n
   
7.         for (r2 = 1; r2 <= N; r2++) { //r2 从 1 到 N
   
8.             for (r1 = 1; r1 <= N; r1++) { //r1 从 1 到 N
   
9.                 fz = r1 * r2;  //计算分子 = r1r2
   
10.                 root = sqrtNumber(fz);//尝试开方，如果 fz 不是完全平方数则返回 -1
    
11.                 if (-1 != root) {  //如果是完全平方数
    
12.                     fm = r1 + r2 + 2 * root; //则计算分母
    
13.                     if (0 == fz % fm) { //如果分母可以整除分子
    
14.                         r3 = fz / fm;  //则计算 r3
    
15.                         n = FastMath.max(r1, r2) + 1;//取r1,r2大者+1作为 n
    
16.                         S[n] += r1 + r2 + r3; //累加
    
17.                     }
    
18.                 }
    
19.             }
    
20.         }
    
21.         //打印结果(略)
    
22. }

复制代码

嵌套循环结果如下，S(100)根本无解，还是我理解有误？

S(5)=9	S(9)=18	S(13)=27	S(17)=36	S(21)=45
S(25)=54	S(29)=63	S(33)=72	S(37)=179	S(41)=90
S(45)=99	S(49)=108	S(53)=117	S(57)=126	S(61)=135
S(65)=144	S(69)=153	S(73)=358	S(77)=171	S(81)=180
S(85)=189	S(89)=198	S(93)=207	S(97)=216	S(101)=225
S(105)=234	S(109)=537	S(113)=252	S(117)=261	S(121)=270
S(125)=279	S(129)=288	S(133)=297	S(137)=306	S(141)=315
S(145)=1054	S(149)=333	S(153)=342	S(157)=351	S(161)=360
S(165)=369	S(169)=378	S(173)=387	S(177)=396	S(181)=895
S(185)=414	S(189)=423	S(193)=432	S(197)=441	S(201)=450
S(205)=459	S(209)=468	S(213)=477	S(217)=1074	S(221)=495
S(225)=504	S(226)=722	S(229)=513	S(233)=522	S(237)=531
S(241)=540	S(245)=549	S(249)=558	S(253)=1253	S(257)=576
S(261)=585	S(265)=594	S(269)=603	S(273)=612	S(277)=621
S(281)=630	S(285)=639	S(289)=2108	S(293)=657	S(297)=666
S(301)=675	S(305)=684	S(309)=693	S(313)=702	S(317)=711
S(321)=720	S(325)=1611	S(329)=738	S(333)=747	S(337)=756
S(341)=765	S(345)=774	S(349)=783	S(353)=792	S(357)=801
S(361)=1790	S(365)=819	S(369)=828	S(373)=837	S(377)=846
S(381)=855	S(385)=864	S(389)=873	S(393)=882	S(397)=1969
S(401)=1782	S(405)=909	S(409)=918	S(413)=927	S(417)=936
S(421)=945	S(425)=954	S(429)=963	S(433)=3162	S(437)=981
S(441)=990	S(445)=999	S(449)=1008	S(451)=1444	S(453)=1017
S(457)=1026	S(461)=1035	S(465)=1044	S(469)=2327	S(473)=1062
S(477)=1071	S(481)=1080	S(485)=1089	S(489)=1098	S(493)=1107
S(497)=1116	S(501)=1125	S(505)=2506	S(509)=1143	S(513)=1152
S(517)=1161	S(521)=1170	S(525)=1179	S(529)=1188	S(533)=1197
S(537)=1206	S(541)=2685	S(545)=1224	S(549)=1233	S(553)=1242
S(557)=1251	S(561)=1260	S(565)=1269	S(569)=1278	S(573)=1287
S(577)=4216	S(581)=1305	S(585)=1314	S(589)=1323	S(593)=1332
S(597)=1341	S(601)=1350	S(605)=1359	S(609)=1368	S(613)=3043
S(617)=1386	S(621)=1395	S(625)=1404	S(629)=1413	S(633)=1422
S(637)=1431	S(641)=1440	S(645)=1449	S(649)=3222	S(653)=1467
S(657)=1476	S(661)=1485	S(665)=1494	S(669)=1503	S(673)=1512
S(676)=2166	S(677)=1521	S(681)=1530	S(685)=3401	S(689)=1548
S(693)=1557	S(697)=1566	S(701)=1575	S(705)=1584	S(709)=1593
S(713)=1602	S(717)=1611	S(721)=5270	S(725)=1629	S(729)=1638
S(733)=1647	S(737)=1656	S(741)=1665	S(745)=1674	S(749)=1683
S(753)=1692	S(757)=3759	S(761)=1710	S(765)=1719	S(769)=1728
S(773)=1737	S(777)=1746	S(781)=1755	S(785)=4502	S(789)=1773
S(793)=3938	S(797)=1791	S(801)=3564	S(805)=1809	S(809)=1818
S(813)=1827	S(817)=1836	S(821)=1845	S(825)=1854	S(829)=4117
S(833)=1872	S(837)=1881	S(841)=1890	S(845)=1899	S(849)=1908
S(853)=1917	S(857)=1926	S(861)=1935	S(865)=6324	S(869)=1953
S(873)=1962	S(877)=1971	S(881)=1980	S(885)=1989	S(889)=1998
S(893)=2007	S(897)=2016	S(901)=9285	S(905)=2034	S(909)=2043
S(913)=2052	S(917)=2061	S(921)=2070	S(925)=2079	S(929)=2088
S(933)=2097	S(937)=4654	S(941)=2115	S(945)=2124	S(949)=2133
S(953)=2142	S(957)=2151	S(961)=2160	S(965)=2169	S(969)=2178
S(973)=4833	S(977)=2196	S(981)=2205	S(985)=2214	S(989)=2223
S(993)=2232	S(997)=2241	S(1001)=2250	S(1005)=2259	S(1009)=7378
S(1013)=2277	S(1017)=2286	S(1021)=2295	S(1025)=2304	S(1029)=2313
S(1033)=2322	S(1037)=2331	S(1041)=2340	S(1045)=5191	S(1049)=2358
S(1053)=2367	S(1057)=2376	S(1061)=2385	S(1065)=2394	S(1069)=2403
S(1073)=2412	S(1077)=2421	S(1081)=5370	S(1085)=2439	S(1089)=2448
S(1093)=2457	S(1097)=2466	S(1101)=2475	S(1105)=2484	S(1109)=2493
S(1113)=2502	S(1117)=5549	S(1121)=2520	S(1125)=2529	S(1126)=3610
S(1129)=2538	S(1133)=2547	S(1137)=2556	S(1141)=2565	S(1145)=2574
S(1149)=2583	S(1153)=8432	S(1157)=2601	S(1161)=2610	S(1165)=2619
S(1169)=2628	S(1173)=2637	S(1177)=2646	S(1181)=2655	S(1185)=2664
S(1189)=5907	S(1193)=2682	S(1197)=2691	S(1201)=5346	S(1205)=2709
S(1209)=2718	S(1213)=2727	S(1217)=2736	S(1221)=2745	S(1225)=6086
S(1226)=3042	S(1229)=2763	S(1233)=2772	S(1237)=2781	S(1241)=2790
S(1245)=2799	S(1249)=2808	S(1253)=2817	S(1257)=2826	S(1261)=6265
S(1265)=2844	S(1269)=2853	S(1273)=2862	S(1277)=2871	S(1281)=2880
S(1285)=2889	S(1289)=2898	S(1293)=2907	S(1297)=9486	S(1301)=2925
S(1305)=2934	S(1309)=2943	S(1313)=2952	S(1317)=2961	S(1321)=2970
S(1325)=2979	S(1329)=2988	S(1333)=6623	S(1337)=3006	S(1341)=3015
S(1345)=3024	S(1349)=3033	S(1351)=4332	S(1353)=3042	S(1357)=3051
S(1361)=3060	S(1365)=3069	S(1369)=6802	S(1373)=3087	S(1377)=3096
S(1381)=3105	S(1385)=3114	S(1389)=3123	S(1393)=3132	S(1397)=3141
S(1401)=3150	S(1405)=6981	S(1409)=3168	S(1413)=3177	S(1417)=3186
S(1421)=3195	S(1425)=3204	S(1429)=3213	S(1433)=3222	S(1437)=3231
S(1441)=10540	S(1445)=3249	S(1449)=3258	S(1453)=3267	S(1457)=3276
S(1461)=3285	S(1465)=3294	S(1469)=3303	S(1473)=3312	S(1477)=7339
S(1481)=3330	S(1485)=3339	S(1489)=3348	S(1493)=3357	S(1497)=3366
S(1501)=3375	S(1505)=3384	S(1509)=3393	S(1513)=7518	S(1517)=3411
S(1521)=3420	S(1525)=3429	S(1529)=3438	S(1533)=3447	S(1537)=3456
S(1541)=3465	S(1545)=3474	S(1549)=7697	S(1553)=3492	S(1557)=3501
S(1561)=3510	S(1565)=3519	S(1569)=9004	S(1573)=3537	S(1576)=5054
S(1577)=3546	S(1581)=3555	S(1585)=11594	S(1589)=3573	S(1593)=3582
S(1597)=3591	S(1601)=11930	S(1605)=3609	S(1609)=3618	S(1613)=3627
S(1617)=3636	S(1621)=8055	S(1625)=3654	S(1629)=3663	S(1633)=3672
S(1637)=3681	S(1641)=3690	S(1645)=3699	S(1649)=3708	S(1653)=3717
S(1657)=8234	S(1661)=3735	S(1665)=3744	S(1669)=3753	S(1673)=3762
S(1677)=3771	S(1681)=3780	S(1685)=3789	S(1689)=3798	S(1693)=8413
S(1697)=3816	S(1701)=3825	S(1705)=3834	S(1709)=3843	S(1713)=3852
S(1717)=3861	S(1721)=3870	S(1725)=3879	S(1729)=12648	S(1733)=3897
S(1737)=3906	S(1741)=3915	S(1745)=3924	S(1749)=3933	S(1753)=3942
S(1757)=3951	S(1761)=3960	S(1765)=12469	S(1769)=3978	S(1773)=3987
S(1777)=3996	S(1781)=4005	S(1785)=4014	S(1789)=4023	S(1793)=4032
S(1797)=4041	S(1801)=18570	S(1805)=4059	S(1809)=4068	S(1813)=4077
S(1817)=4086	S(1821)=4095	S(1825)=4104	S(1829)=4113	S(1833)=4122
S(1837)=9129	S(1841)=4140	S(1845)=4149	S(1849)=4158	S(1853)=4167
S(1857)=4176	S(1861)=4185	S(1865)=4194	S(1869)=4203	S(1873)=13702
S(1877)=4221	S(1881)=4230	S(1885)=4239	S(1889)=4248	S(1893)=4257
S(1897)=4266	S(1901)=4275	S(1905)=4284	S(1909)=9487	S(1913)=4302
S(1917)=4311	S(1921)=4320	S(1925)=4329	S(1929)=4338	S(1933)=4347
S(1937)=4356	S(1941)=4365	S(1945)=9666	S(1949)=4383	S(1953)=4392
S(1957)=4401	S(1961)=4410	S(1965)=4419	S(1969)=4428	S(1973)=4437
S(1977)=4446	S(1981)=9845	S(1985)=4464	S(1989)=4473	S(1993)=4482
S(1997)=4491	S(2001)=8910

王守恩

S(04)=A(1)=1+4+4=9,
S(08)=A(2)=2+8+8+9=27,
S(12)=A(3)=3+12+12+27=54,
S(16)=A(4)=4+16+16+54=90,
S(20)=A(5)=5+20+20+90=135,
S(24)=A(6)=6+24+24+135=189,
S(28)=A(7)=7+28+28+189=252,
S(32)=A(8)=8+32+32+252=324,
S(36)=A(9)=9+36+36+4+9+36=454,
S(72)=A(18)=18+72+72+8+18+72+1426=1686,
S(96)=A(24)=24+96+96+2631=2847,
S(100)=A(25)=25+100+100=3072,

{9, 27, 54, 90, 135, 189, 252, 324, 454, 544, 643, 751, 868, 994, 1129, 1273, 1426, 1686, 1857, 2037, 2226, 2424, 2631, 2847, 3072, 3306, 3696, 3948, 4209, 4479}

Table[Sum[Floor[(4 k)/(a (a + b))^2 + 1] Floor[(4 k)/(a (a + b))^2] (a^2 + b^2 + a b)^2/2, {a, Power[4 k, (4)^-1]}, {b, Select[Range@a, CoprimeQ[a, #] &]}], {k, 30}]

S(1234)=A(308)==475425,

王守恩

四来 发表于 2025-3-8 21:11
根据以上公式，嵌套循环求解, 取 N=2000

简单暴力粗犷的 Java 代码:

F(04)=B(1)=1+4+4=9,
F(08)=B(2)=2+8+8=18,
F(12)=B(3)=3+12+12=27,
F(16)=B(4)=4+16+16=36,
F(20)=B(5)=5+20+20=45,
F(24)=B(6)=6+24+24=54,
F(28)=B(7)=7+28+28=63,
F(32)=B(8)=8+32+32=72,
F(36)=B(9)=9+36+36+4+9+36=130,

{9, 18, 27, 36, 45, 54, 63, 72, 130, 90, 99, 108, 117, 126, 135, 144, 153, 260, 171, 180, 189, 198, 207, 216, 225, 234, 390, 252——这通项公式我也没有。

王守恩

三个圆的半径都是整数,三个半径之和=1521={169,676,676}+{81,144,1296},我们说这个是2组解。

三个圆的半径都是整数,三个半径之和=441={49,196,196}+{36,81,324}+{16,25,400},我们说这个是3组解。

会有更多组解的吗？！

northwolves

王守恩 发表于 2025-3-10 15:14
三个圆的半径都是整数,三个半径之和=1521={169,676,676}+{81,144,1296},我们说这个是2组解。

三个圆的半径 ...

1亿内的解(原来的计算有问题，重发，？ 暂未确认)：

前15项：
a = {9, 2401, 441, 8281, 21609, ? , 405769, ?, 74529, 19882681, ?, 68574961, 2989441, ?, 3651921}

a(22)=146482609

a(27)=26904969

a(40)=2872852801

northwolves

测试发现：

723877356890961有125组解

45594931810762921041有225组解

942822584306490332808429866283371017321532378401 有7654组解

王守恩

northwolves 发表于 2025-3-10 16:06
1亿内的解：

前15项：

小心翼翼的问:  我只有8个解？a(23)=102030201,

1. Solve[{z == a^2*b, y == (1 + a)^2*b, x == b (1 + a)^2*a^2, x + y + z == b (1 + a + a^2)^2 == 102030201, a > 0, b > 0}, {a, b, z, y, x}, Integers]

复制代码

{a -> 1, b -> 11336689, z -> 11336689, y -> 45346756, x -> 45346756},
{a -> 2, b -> 2082249, z -> 8328996, y -> 18740241, x -> 74960964},
{a -> 3, b -> 603729, z -> 5433561, y -> 9659664, x -> 86936976},
{a -> 4, b -> 231361, z -> 3701776, y -> 5784025, x -> 92544400},
{a -> 9, b -> 12321, z -> 998001, y -> 1232100, x -> 99800100},
{a -> 10, b -> 8281, z -> 828100, y -> 1002001, x -> 100200100},
{a -> 16, b -> 1369, z -> 350464, y -> 395641, x -> 101284096},
{a -> 100, b -> 1, z -> 10000, y -> 10201, x -> 102010000}}

northwolves

{27,102030201,{{10000,10201,102010000},{116964,125316,101787921},{350464,395641,101284096},{828100,1002001,100200100},{998001,1232100,99800100},{1368900,1758276,98903025},{2047761,2802276,97180164},{2775556,4036081,95218564},{3168400,4752400,94109401},{3701776,5784025,92544400},{3980025,6350400,91699776},{5433561,9659664,86936976},{5721664,10394176,85914361},{6426225,12320100,83283876},{6708100,13147876,82174225},{8076964,17757796,76195441},{8328996,18740241,74960964},{8803089,20738916,72488196},{9144576,22325625,70560000},{9759376,25603600,66667225},{10310521,29246464,62473216},{10601536,31629376,59799289},{10732176,32867289,58430736},{11088900,37271025,53670276},{11256025,40704400,50069776},{11329956,43996689,46703556},{11336689,45346756,45346756}}}

mathe

\(102030201=3^2\times 7^2\times 13^2\times 37^2\)
这个数字选择的比较好。
根据15#，\(r_1+r_2+r_3=(a^2+ab+b^2)^2 h\)
所以我们需要寻找这些数字有更多的上面表达形式即可。
而\(a^2+ab+b^2\)有个特点，其所有素因子要么是3，要么模3余1.
而且记\(\omega\)为三次单位根，那么\(a^2+ab+b^2=(a-\omega b)(a-\bar{\omega}b)\), 也就是是它总是爱森斯坦整环中整数的模长平方形式。
所以由于上面3,7,13,37都是在爱森斯坦整环中可以分解，比如\(3=(1-\omega)(1-\bar{\omega}), 7=(2-\omega)(2-\bar{\omega})=(1-2\omega)(1-2\bar{\omega}), ...\)
由此我们可以选择把部分上面表达形式相乘，得到一些复合形式，如
\(21=(1-\omega)(2-\omega)(1-\bar{\omega})(2-\bar{\omega})=(1-4\omega)(1-4\bar{\omega})\)
也可以表示为
\(21=(1-\omega)(1-2\omega)(1-\bar{\omega})(1-2\bar{\omega})=(-1-5\omega)(-1-5\bar{\omega})=(1+5\omega)(1+5\bar{\omega})\)
其中第一种形式对应a=1,b=4符合都是正数要求，但是第二种a=1,b=-5不符合要求。
另外对于每个\(a-b\omega\)，我们分别将它继续乘上\(\omega,\omega^2\)还可以得到另外两种表达形式，如
\(b+(a+b)\omega, (-a - b)+(-a)\omega\), 再加上它们的之间a,b轮换可以有6种不同形式。
所以最终需要穷举的数目如下
i) 3可以选择\(1-\omega\)或者不参与(h*=9)，两种选择
ii) 7可以选择\(1-2\omega\)或者\(2-\omega\)或者不参与，三种选择
iii)13可以选择\(1-3\omega\)或者\(3-\omega\)或者不参与，三种选择
iv)37可以选择\(3-4\omega\)或者\(4-3\omega\)或者不参与，三种选择。
所以上面共2*3*3*3=54种选择，对每个结果，再可以继续乘上\(\omega,\omega^2\)，所以共需枚举54*3=162种情况，选择其中a,b同号的情况。