关于单位分数的一些难题

mathe

采用25#附件中给出的算法,让N=8,k=6
对于序列$n_1,n_2,...,n_8$,定义$R_0=1,R_1=n_1-1,R_t=n_tR_{t-1}-\prod_{i=1}^{t-1}n_i$
对于候选序列$n_1,n_2,...,n_k$,我们需要判断
$(\prod_{i=1}^kn_i)^2+R_k$是否有关于$R_k$和$-\prod_{i=1}^kn_i$同余的因子
由此在${\prod_{i=1}^k n_i}/{R_k}$不是太大的情况(比如小于1000),我们可以穷举所有这样可能的因子.
但是如果合格比例比较大的时候,可能采用因子分解$(\prod_{i=1}^kn_i)^2+R_k$的方法比较.
不过我们可以用C语言写一个程序,将比例比较小的情况全部处理,然后余下的部分输出到一个文件,在用其它工具对这些数据进行因子分解,然后继续处理.
当然文章里面还有一个重要的筛选条件是对于任意$1<=t<=k$,
${\prod_{i=1}^tn_i}/{R_t}<n_{t+1}<(N-t){\prod_{i=1}^tn_i}/{R_t}$
谁有空可以先将这部分功能实现一下,然后根据处理结果我们就可以看出到底还有多少工具量遗留.

medie2005

注意到，搜索过程中能搜到最后一层的几乎都会产生解。
是否可以由这个来提前剪枝？

数学星空

以下附件对于mathe研究计算s(8) 可能会有帮助...
S0025-5718-99-01088-1.pdf (269.38 KB, 下载次数: 2)

2010-1-14 09:28 上传

点击文件名下载附件

数学星空

下面附件是有关如何计算单位分数的一些算法研究及计算程序...
或许对编程有所帮助
EgyptianFractions.pdf (130.49 KB, 下载次数: 1)

2010-1-14 09:46 上传

点击文件名下载附件

mathe

注意到，搜索过程中能搜到最后一层的几乎都会产生解。
是否可以由这个来提前剪枝？
medie2005 发表于 2010-1-14 09:23

什么叫搜索到最后一层,直接搜索到8层吗?

我现在写了个程序,统计了对s(8)搜索到6层产生的数据分布情况
发现计算出来的$\prod_{i=1}^6n_i$没有超过14位数的,这是一个很好的信息.
也就是说$(\prod_{i=1}^6n_i)^2+R_6$不会超过28位数.
也就是我们需要一个能够对不超过28位数进行确定性素性判断和因子分解的程序.
而搜索到6层产生的所有候选解只有3382746个,不算太大.
但是其中${\prod_{i=1}^6n_i}/{R_6}$通常都比较大,主要分布在7位数,所以如果如25#论文采用试初法的确效率不高,而采用因子分解的方法应该更加有效

mathe

统计前6层的代码,大家看看有没有错误,完全按照那个论文上的算法写的

1. 
   
2. #include "stdafx.h"
   
3. #include <math.h>
   
4. #define L 40
   
5. #define N 8
   
6. #define K 6
   
7. unsigned n[K];
   
8. unsigned long long M[K];
   
9. unsigned long long R[K];
   
10. unsigned long long maxM;
    
11. int A[4];
    
12. long long count[K];
    
13. long long gc[L];
    
14. void output(int level, int overflow)
    
15. {
    
16.     if(level==K){
    
17.         int a3=A[3]%4;
    
18.         int a1=A[1]%4;
    
19.         if(A[2]==1&&a1==1)
    
20.             return;
    
21.         if(A[0]==1&&a3==0)
    
22.             return;
    
23.         if(A[0]==0&&A[2]==0){
    
24.             if((a1==1||a1==2)&&a3<2)
    
25.                 return;
    
26.         }
    
27.         int c=(int)ceil(log10((double)(M[level-1])));
    
28.         if(c>=L)c=L-1;
    
29.         if(c<0)c=0;
    
30.         gc[c]++;
    
31.     }
    
32.     if(M[level-1]>maxM)maxM=M[level-1];
    
33.     count[level-1]++;
    
34. }
    
35. 
    
36. unsigned long long gcd(unsigned long long x,unsigned long long y)
    
37. {
    
38.     while(y!=0){
    
39.         unsigned long long r=x%y;
    
40.         x=y;
    
41.         y=r;
    
42.     }
    
43.     return x;
    
44. }
    
45. 
    
46. void search(int level)
    
47. {
    
48.     if(level>=K){
    
49.         output(level,0);
    
50.         return;
    
51.     }
    
52.     double ll=(double)M[level-1]/R[level-1];
    
53.     double uu=ll*(N-level);
    
54.     if(uu>=4294967295.5||uu>=18446744073709551615.0/M[level-1]){///overflow
    
55.         output(level,1);
    
56.         return;
    
57.     }
    
58.     unsigned li=(unsigned)ceil(ll);
    
59.     unsigned ui=(unsigned)uu;
    
60.     if(li<=n[level-1])li=n[level-1]+1;
    
61.     unsigned i;
    
62.     for(i=li;i<=ui;i++){
    
63.         if(gcd(M[level-1],i)!=1)
    
64.             continue;
    
65.         n[level]=i;
    
66.         M[level]=M[level-1]*i;
    
67.         R[level]=R[level-1]*i-M[level-1];
    
68.         A[i%4]++;
    
69.         search(level+1);
    
70.         A[i%4]--;
    
71.     }
    
72. }
    
73. 
    
74. void search0()
    
75. {
    
76.     int i;
    
77.     for(i=0;i<4;i++)A[i]=0;
    
78.     for(i=2;i<=7;i++){
    
79.         n[0]=i;
    
80.         R[0]=i-1;
    
81.         M[0]=i;
    
82.         A[i%4]++;
    
83.         search(1);
    
84.         A[i%4]--;
    
85.     }
    
86. }
    
87. 
    
88. int _tmain(int argc, _TCHAR* argv[])
    
89. {
    
90.     search0();
    
91.     int i;
    
92.     for(i=0;i<K;i++){
    
93.         printf("count[%d]=%lld\n",i+1,count[i]);
    
94.     }
    
95.     for(i=0;i<L;i++)
    
96.     printf("gc[%d]=%lld\n",i,gc[i]);
    
97.     printf("maxM=%llud\n",maxM);
    
98.         return 0;
    
99. }

复制代码

medie2005

什么叫搜索到最后一层,直接搜索到8层吗?
==============================
我是指搜索到第6层。后两个解不定方程就可以（看s(7)的结果，要解的方程是很少的）。

mathe

那可以s(8)的情况不同,我上面估计出来还有
3382746
个候选的情况

mathe

你是如何过滤的,我试了一下,s(7)情况余下1378个方程,好像解远远没有这么多

medie2005

呵呵，我是猜测的，没实际算过。最后的不定方程有解的可能性比较大。