可表示为连续正整数平方和的四次方数

数论爱好者

northwolves 发表于 2025-5-3 22:46
我打开了10个进程计算到大约10^7了，目前就发现10个解：
13, 295, 330, 364, 1085, 5005, 6305, 15516, 415 ...

经过对比,wayne的代码不是最快.我用Wolfram Mathematica 8,验证区间199万至200万,两分钟后得出:{1990368,{{47761921,1146120288,549179184}}},先给出一个答案,继续运行,平均每秒跑3个整数.5分钟后跑到1990836,距离200万还很远.但是,我的maple代码在450秒跑完199万至200万,然后一次性给出答案.

1. restart;
   
2. with(numtheory):
   
3. with(Threads):
   
4. 
   
5. findN_highspeed := proc(m_start, m_end)
   
6.     local m, mm, div, d, S, p, solutions;
   
7.     solutions := [];
   
8.     for m from m_start to m_end do
   
9.         mm := 12*m^4;
   
10.         div := select(d -> d^3 - d < mm, divisors(mm)[2..-2]);
    
11.         
    
12.         for d in div do
    
13.             S := (mm/d + 1 - d^2)/3;
    
14.             if S <= 0 then next; end if;
    
15.             
    
16.             if issqr(S) and type(sqrt(S), integer) then
    
17.                 p := (sqrt(S) + 1 - d)/2;
    
18.                 if type(p, posint) and (d - sqrt(S)) mod 2 = 1 then
    
19.                     printf("Found m=%d, d=%d, p=%d\n", m, d, p);
    
20.                     solutions := [op(solutions), m];
    
21.                     break;
    
22.                 end if;
    
23.             end if;
    
24.         end do;
    
25.     end do;
    
26.     return solutions;
    
27. end proc;
    
28. 
    
29. # 并行执行设置（分割搜索范围）
    
30. search_range := 1990000..2000000:
    
31. num_threads := 4;
    
32. total_range := rhs(search_range) - lhs(search_range) + 1;
    
33. step := iquo(total_range, num_threads, 'remainder');
    
34. 
    
35. ranges := [seq(
    
36.     lhs(search_range) + (i-1)*step ..
    
37.     `if`(i <= remainder, lhs(search_range) + i*step, lhs(search_range) + i*step - 1),
    
38.     i = 1..num_threads
    
39. )];
    
40. 
    
41. # 启动并行计算
    
42. results := Threads:-Map(r -> findN_highspeed(op(1,r), op(2,r)), ranges);
    
43. final_result := ListTools:-Flatten(results);

复制代码

wayne

1) 你说的不快，是啥意思，是你改进了Mathematica代码吗，改进的地方在哪，
还是你用2010年发布的 Mathematica 8 跟 Maple比吗，你的 Maple也是2010年的吗。先不说拿15年前的A产品跟现在的B产品比较是否合理，就算都是同一年发布的，不同的软件比较起来也没啥好讨论的余地。
更何况我粗略看了你的Maple代码，算法逻辑好像是一样的。

2) 要打破记录不难。 换成C/Rust肯定可以很快。
再不济，用大家的代码稍微改改，随机指定几个区间，凭借运气也能碰到几个， 关键是是要没有遗漏的找到所有的解。

数论爱好者

wayne 发表于 2025-5-3 23:21
1) 你说的不快，是啥意思，是你改进了Mathematica代码吗，改进的地方在哪，
还是你用2010年发布的 Mathemat ...

你说的也对,参考你的代码,并无实质性改进,maple版本2013年的....

数论爱好者

数论爱好者 发表于 2025-5-4 00:16
你说的也对,参考你的代码,并无实质性改进,maple版本2013年的....

二者确实不能够比较,Mathematica要好得多.
Mathematica 8:CPU占用50%,物理内存占用27%,Mathematica 8软件内存100MB左右
maple 17:CPU占用70%-100%,物理内存占用60%,maple 17软件内存1.2GB左右

northwolves

对于已有的前10个解：
m = {13, 295, 330, 364, 1085, 5005, 6305, 15516, 415151,  1990368};
n = {2, 177, 352, 1536, 2401, 40898, 60625, 185761,  19512097, 47761921}

m因数分解：{{{13,1}},{{5,1},{59,1}},{{2,1},{3,1},{5,1},{11,1}},{{2,2},{7,1},{13,1}},{{5,1},{7,1},{31,1}},{{5,1},{7,1},{11,1},{13,1}},{{5,1},{13,1},{97,1}},{{2,2},{3,2},{431,1}},{{11,2},{47,1},{73,1}},{{2,5},{3,2},{6911,1}}}

n因数分解：{{{2,1}},{{3,1},{59,1}},{{2,5},{11,1}},{{2,9},{3,1}},{{7,4}},{{2,1},{11,2},{13,2}},{{5,4},{97,1}},{{431,2}},{{11,2},{47,2},{73,1}},{{6911,2}}}

wayne

northwolves 发表于 2025-5-4 07:57
对于已有的前10个解：
m = {13, 295, 330, 364, 1085, 5005, 6305, 15516, 415151,  1990368};
n = {2, 17 ...

计算了一下，发现只有$\frac{12m^4}{n}$是整数，而$\frac{m}{n}，\frac{m^2}{n}，\frac{m^3}{n}，\frac{m^4}{n}$都不能保证。

northwolves

m = 9491625, n = 284731876 是一个很接近的解，可惜p不是整数。

northwolves

A001032 is helpful filting the divisors of 12*m^4.

wayne

也就是说，咱们应该先讨论 【可表示为连续正整数平方和的平方数】

northwolves

https://oeis.org/draft/A383359

David A. Corneth 用很短时间就刷新到2*10^7了。我的机子现在还在计算1.4*10^7