数学研发论坛

 找回密码
 欢迎注册
查看: 541|回复: 34

[擂台] 求 n 使得 n! 以 2021 开头

[复制链接]
发表于 2021-1-27 17:37:27 | 显示全部楼层 |阅读模式

马上注册,结交更多好友,享用更多功能,让你轻松玩转社区。

您需要 登录 才可以下载或查看,没有帐号?欢迎注册

x
某书证明了存在 n,使得 n! 可以以任意一串数字作为开头。现求正整数 n,使得 n! 以2021 开头。
2021-01-27_173701.png
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-27 18:03:42 | 显示全部楼层
16979
17791
24764
26619
27727
34990
42004
57216
71869
74494
77981
81253
87457
98508

点评

应该都是取对数之后再比较大小的吧?  发表于 2021-1-27 19:44
服!老大就是老大!  发表于 2021-1-27 18:28
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2021-1-27 20:22:16 | 显示全部楼层
使用 Mathematica 硬算,竟然也好使:
  1. n = 11; prod = 10!;
  2. While[n < 100000, prod = prod*n;
  3. If[Take[RealDigits[prod][[1]], 4] == {2, 0, 2, 1}, Print[n]];
  4. n = n + 1
  5. ]
复制代码
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-27 20:23:21 | 显示全部楼层
  1. {16979,17791,24764,26619,27727,34990,42004,57216,71869,74494,77981,81253,87457,98508,101646,106752,115082,117536,123759,124963,125097,128930,129496,143176,145430,145608,147954,150084,159119,160231,161163,168079,170851,175620,179179,186731,194395,197814,200531,202408,205634,220996,229949,232031,232165,232359,233538,233812,236569,239429,241831,244379,256844,258008,260086,266522,274849,279710,282887,283345,283429,287093,297038,313125,320693,324548,327849,339662,343693,344403,347492,347553,348061,355546,357943,358607,359862,366076,369742,375529,384520,385526,411111,414986,427070,428915,442399,446228,446642,446722,452782,453888,456909,459849,463837,465062,465878,470988,474609,478127,480111,490862,496341,499823,502028,506166,507946,510457,517036,528199,532074,543024,551055,554680,555384,556545,562137,562545,573616,574616,591316,597692,603304,606644,608358,615043,616439,621847,624892,630419,631890,641452,654517,661103,661844,665506,667020,675154,676942,678469,679987,681658,682597,682784,687780,687940,689781,699900,708254,709378,710487,722016,727589,728845,732825,736823,743496,749697,754686,765342,766274,770955,777918,783395,785148,786380,793419,797545,797657,799510,807757,819293,843222,845303,847777,849418,858218,867231,872432,878030,888901,889824,894485,906501,908758,908999,913891,916410,916542,920553,922326,922440,930970,937076,946271,949542,953722,972282,972778,981548,982177,988996}
复制代码

点评

能否给个代码欣赏一下?  发表于 2021-1-27 20:24
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-27 20:34:14 | 显示全部楼层
无聊的跑了10分钟,刷出了10^8的所有结果,突破了我就给代码,
  1. {16979,17791,24764,26619,27727,34990,42004,57216,71869,74494,77981,81253,87457,98508,101646,106752,115082,117536,123759,124963,125097,128930,129496,143176,145430,145608,147954,150084,159119,160231,161163,168079,170851,175620,179179,186731,194395,197814,200531,202408,205634,220996,229949,232031,232165,232359,233538,233812,236569,239429,241831,244379,256844,258008,260086,266522,274849,279710,282887,283345,283429,287093,297038,313125,320693,324548,327849,339662,343693,344403,347492,347553,348061,355546,357943,358607,359862,366076,369742,375529,384520,385526,411111,414986,427070,428915,442399,446228,446642,446722,452782,453888,456909,459849,463837,465062,465878,470988,474609,478127,480111,490862,496341,499823,502028,506166,507946,510457,517036,528199,532074,543024,551055,554680,555384,556545,562137,562545,573616,574616,591316,597692,603304,606644,608358,615043,616439,621847,624892,630419,631890,641452,654517,661103,661844,665506,667020,675154,676942,678469,679987,681658,682597,682784,687780,687940,689781,699900,708254,709378,710487,722016,727589,728845,732825,736823,743496,749697,754686,765342,766274,770955,777918,783395,785148,786380,793419,797545,797657,799510,807757,819293,843222,845303,847777,849418,858218,867231,872432,878030,888901,889824,894485,906501,908758,908999,913891,916410,916542,920553,922326,922440,930970,937076,946271,949542,953722,972282,972778,981548,982177,988996,1010616,1013325,1031953,1033245,1035162,1043571,1044638,1045835,1046346,1050392,1050579,1051136,1055948,1058351,1083839,1087892,1093823,1101817,1115838,1115922,1121485,1127931,1128313,1128808,1128865,1130864,1130995,1134176,1135142,1136426,1136498,1140001,1144434,1144792,1145302,1148876,1162926,1168305,1170702,1179688,1183182,1186279,1188371,1210638,1220499,1221272,1225837,1227651,1230031,1231528,1232411,1240899,1249064,1260279,1278026,1282555,1286628,1287048,1292226,1294330,1298675,1301238,1303339,1305928,1306066,1308378,1310036,1310760,1312754,1321407,1322257,1324025,1328826,1337609,1340299,1340739,1341547,1341876,1354060,1363894,1364539,1370147,1373860,1376015,1376923,1388870,1395426,1395944,1398237,1400512,1406591,1407508,1410128,1427201,1429410,1435955,1437645,1458822,1466927,1469535,1473386,1479761,1483737,1484227,1501671,1504186,1508340,1512522,1526343,1526490,1526958,1529786,1531090,1542711,1550132,1553888,1555125,1557504,1559678,1561127,1568230,1572536,1578871,1579118,1594172,1606070,1620370,1624220,1632212,1632865,1632973,1633372,1636627,1648284,1650580,1652558,1664419,1665196,1666279,1667230,1677012,1680541,1692051,1697981,1702549,1705120,1707107,1717127,1734340,1736366,1743935,1746321,1746544,1752082,1753436,1763439,1768286,1770465,1770731,1785561,1788888,1800624,1809649,1816780,1823025,1823669,1824033,1827842,1828487,1830275,1831284,1839497,1839565,1844469,1849187,1854251,1854512,1857855,1860380,1879112,1881730,1890020,1909653,1912142,1914704,1917052,1921280,1926200,1934752,1936969,1938257,1944328,1946615,1949719,1955414,1962540,1967418,1973681,1977543,1989783,1997106,2007688,2013520,2016028,2022836,2029244,2039315,2042981,2045200,2047602,2052106,2054436,2055312,2056443,2060675,2061385,2070363,2077371,2077582,2078432,2079530,2089561,2094605,2095539,2102805,2105603,2106858,2114324,2116681,2117338,2121657,2133299,2143013,2150166,2152198,2155736,2156671,2162049,2163980,2179251,2185073,2186089,2201412,2202780,2229301,2233627,2234452,2242992,2243505,2261649,2264581,2266034,2267978,2275084,2277925,2279172,2279390,2280407,2284450,2289805,2290333,2294524,2295447,2298739,2300235,2303436,2308617,2310724,2311642,2318442,2318672,2322239,2323015,2334488,2338181,2348019,2352029,2353528,2355808,2356117,2359824,2363294,2363551,2378075,2386796,2394091,2402438,2404452,2404906,2410315,2417816,2423827,2426421,2437562,2439662,2439915,2440888,2447054,2450029,2451981,2455016,2455521,2457842,2458825,2468282,2470857,2470969,2472338,2486814,2492161,2493119,2494064,2495978,2496197,2517703,2519281,2519737,2525555,2532127,2532915,2536856,2537363,2540257,2542200,2544849,2548689,2553196,2555450,2561372,2566449,2571569,2575941,2579501,2581534,2582221,2583017,2584841,2586259,2587991,2589337,2592445,2594320,2595129,2595636,2602037,2605925,2609214,2615348,2615985,2618706,2622296,2626768,2631056,2633554,2636565,2636736,2639748,2641031,2641645,2647233,2656613,2660115,2660468,2667234,2669709,2672266,2672875,2673917,2680137,2690638,2691780,2694713,2697485,2699918,2703333,2704583,2705534,2709762,2727072,2730210,2732978,2739847,2741021,2744077,2745557,2749175,2755059,2758583,2760205,2762302,2766437,2769341,2772732,2778443,2790100,2790735,2794672,2797375,2804194,2804424,2813057,2814991,2825757,2832364,2835388,2836265,2839627,2862173,2870980,2873385,2881496,2884597,2889707,2895794,2896803,2902513,2910235,2911653,2912603,2912631,2912659,2916669,2929883,2936621,2937478,2939533,2940560,2944460,2952452,2958616,2971266,2973684,2986650,2989515,2994444,2999810,3001566,3007002,3024162,3033342,3033425,3035518,3039087,3046809,3055398,3056171,3058215,3058767,3068618,3069503,3074226,3078275,3085472,3088628,3092544,3095337,3099907,3101785,3106436,3107086,3118387,3126771,3143066,3147316,3148063,3149125,3152632,3154749,3159880,3172660,3176254,3179098,3179737,3184022,3185485,3192103,3201348,3209581,3211853,3220441,3220756,3222061,3229090,3229255,3229911,3230074,3232736,3235604,3239614,3241915,3244825,3256753,3263989,3271431,3274905,3280493,3283804,3284048,3284989,3286013,3294796,3297572,3299374,3302380,3306258,3309986,3318593,3321993,3322739,3323389,3326055,3327327,3333243,3333396,3335826,3339389,3353935,3358385,3367663,3385637,3388687,3389053,3396289,3400874,3416494,3416747,3420660,3422588,3422908,3427163,3438536,3444751,3447630,3449609,3450323,3469782,3470476,3472463,3479366,3481755,3485443,3487425,3495068,3495436,3497967,3499561,3505621,3513652,3519517。。。。。。。。。。。,99352598,99354373,99360445,99369117,99370574,99371669,99372400,99379388,99379758,99383470,99384588,99385708,99393985,99395882,99397023,99399694,99405072,99405458,99409332,99415584,99415977,99421110,99423893,99427491,99435569,99453762,99455868,99458406,99460956,99464375,99466523,99472584,99476079,99480481,99482696,99486260,99489399,99489849,99491653,99492105,99493917,99500771,99501692,99504002,99513352,99525299,99526753,99531142,99534091,99539048,99542550,99546079,99548108,99549126,99552707,99553221,99556317,99565738,99576452,99581363,99591925,99595886,99596455,99598167,99600461,99601613,99603927,99611543,99612728,99615109,99620521,99629717,99637876,99648832,99651460,99657447,99659466,99662858,99676072,99683975,99686165,99689109,99698876,99700407,99704269,99705828,99711351,99716982,99723558,99724391,99732005,99741616,99744300,99745201,99746105,99747923,99748837,99750675,99764918,99774914,99781131,99783243,99787530,99790803,99794128,99814091,99817841,99821670,99824269,99832308,99842204,99849671,99852764,99854335,99860797,99862460,99864143,99867574,99869323,99871096,99884265,99886271,99903812,99908723,99916630,99919438,99922348,99925371,99931818,99935282,99938942,99942836,99947016,99951555,99956566,99956567,99962237,99968926,99977524,99977525,99977526,99993316,99993317,99993318,99993319,99993320,99993321,99993322}
复制代码

点评

开多核跑40秒钟  发表于 2021-1-27 21:44
我算的是 10^8以内是 21436个解  发表于 2021-1-27 21:44
find 21440 sols, use time = 90.092s  发表于 2021-1-27 20:48
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-27 21:46:16 | 显示全部楼层
单核版本的代码:
  1. n=1;Monitor[Reap[While[n<10^5,t=LogGamma[n+1.];If[Log[2021]<t-(Floor[t/Log[10]]-3)Log[10]<Log[2022],Print[Sow[n]]];n++]],n]
复制代码

多核版本的代码:
  1. CloseKernels[];LaunchKernels[12];
  2. Timing[ans=ParallelCombine[Select[#,Log[2021]<Last[#]-(Floor[Last[#]/Log[10]]-3)Log[10]<Log[2022]&]&,Table[{n,LogGamma[n+1.]},{n,10^8}]][[All,1]];]
复制代码
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2021-1-27 23:04:50 | 显示全部楼层
本帖最后由 uk702 于 2021-1-27 23:08 编辑
wayne 发表于 2021-1-27 21:46
单核版本的代码:

多核版本的代码:


我将你的多核版本的结果打印出来,似乎各有各的不准。
31419078!=2022000269906979489591697...,你的多了这个。
39600422!=2021000126106549851567511...,你的少了这个。
48400100!=2020999793431651954566836...,我的多了这个。
48844860!=2021000416575056021475393...,你的少了这个。
80495570!=2020999980052118200212907...,我的多了这个。

我实际给出的解是 21439 个(显示时多算了一个),多了 2 个,你给出了 21436 个,少了 2 个并多了 1 个。故检查后共同的解是 21437 个,暂不排除其中有我们都算错的。

点评

用高精度的跟你的对比,你的还少了三个,61543437,88289946,90711389  发表于 2021-1-28 00:04
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-27 23:52:16 | 显示全部楼层
这是有可能的,因为在做比较的时候用的是double精度,我稍微做了修改,用30位高精度计算,重新跑了下程序,速度慢了不少,算出正确的结果应该有21439 个。
  1. CloseKernels[];LaunchKernels[12];
  2. Timing[ans=ParallelCombine[Select[#,Log[2021]<Last[#]-(Floor[Last[#]/Log[10]]-3)Log[10]<Log[2022]&]&,Table[{n,LogGamma[n+1`30]},{n,10^8}]][[All,1]];]
复制代码
发现前面的double精度的版本产生了错误解31419078,60459888,还少了39600422,48844860,61543437,88289946,90711389。

点评

好!我也看到了,我的并行版本竟然为 52270002 生成了两条记录,不知那出乱子了。  发表于 2021-1-28 10:36
52270002在里面,因为你没有提及  发表于 2021-1-28 10:32
52270002!=20218440718601624628173414618456...,好像你的高精度版本还缺了这个,还是说漏写了。  发表于 2021-1-28 10:25
原来你也用vim  发表于 2021-1-28 08:38
我又跑了50位精度的版本,vimdiff之后,计算结果跟30位高精度的结果完全一致。  发表于 2021-1-28 00:09
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2021-1-28 07:11:27 | 显示全部楼层
我的单核版本,julia 代码,double 精度,使用的某位同学论文中的斯特林加强公式(https://wenku.baidu.com/view/77419c44336c1eb91a375d3d.html),直接计算 log(n) 的累加和在 double 精度下累积误差似乎相当大。



  1. t1=time();
  2. h(n) = (log(sqrt(2.0*pi*n)) + n*log(n)) - n * (1-1/(0.4+12*n^2));
  3. n=10; ct=1; u = log(2021.); v = log(2022.)
  4. while n < 10^8
  5.   global n, ct
  6.   a = (h(n) - u)/log(10.); t = floor(Int, a);
  7.   if h(n) < v + t*log(10) ct = ct+1;println(n) end
  8.   n=n+1
  9. end

  10. # 打印出来的 ct 比实际多 1 个
  11. println("find $ct sols, use time = $(time() - t1)")
复制代码
2021-01-28_070731.png
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2021-1-28 08:47:32 | 显示全部楼层
斯泰林公式可以采用
\(\log(n!) \approx \frac{\log(2\pi n)}2+n\log(n)-n+\frac1{12n}+O(\frac{1}{n^3})\)

点评

n大了就不行了。  发表于 2021-1-28 09:13
有更精确的:log(n!) = log(sqrt(2*pi*n)) + log(n)*n - n + log(1+1/(12*n)+1/(288*n^2)-139/(51840*n^3)-571/(2488320*n^4)),但初步验证似论文给的公式比这个还精确,不过我只验证 n 取小值时,充分大没验证  发表于 2021-1-28 09:11
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2021-3-8 23:22 , Processed in 0.070263 second(s), 20 queries .

Powered by Discuz! X3.4

© 2001-2017 Comsenz Inc.

快速回复 返回顶部 返回列表