搜寻互质点中的空白

KeyTo9_Fans

在趣题妙解板块《由互质点织成的一块布》一贴中的第二问，mathe给出了漂亮的结论以及巧妙的构造：

http://bbs.emath.ac.cn/thread-1946-2-1.html

为此特意摆设一个编程擂台，对于n=2、n=3、n=4的情况，找对应的(x,y)，使得整数对

(x,y)，(x+1,y)，...，(x+n-1,y)
(x,y+1)，(x+1,y+1)，...，(x+n-1,y+1)
............
(x,y+n-1)，(x+1,y+n-1)，...，(x+n-1,y+n-1)

均不互质。

在满足上面的条件的前提下，(x+y)的值越小越好。

看看对于n=2、n=3、n=4的情况，找出来的整数对到底有多大。

northwolves

N=2: (14,20)

northwolves

N=3: (1274 ,1308)

数学星空

这个问题的计算很费时啊,有谁能给出N=4,5,6,7,8,9,10
的值呢??
可能需要借助数论知识来编程,缩短运算时间,有谁能给出答案及程序???

mathe

N=4我们可以试着搜索公因子模式如下的结构，其中?用互素而且不小于5的素数。
6 ?  2 3
? ?  ? ?
2 ?  2 ?
3 ?  ? 3

mathe

Found x=2038235252988,y=2891651852154
Found x=721235275761,y=383376261639
Found x=329601462366,y=322399219389
Found x=8197745490,y=507009915048
Found x=8197745490,y=41984229963
Found x=3333736029,y=394313397
Found x=1668280158,y=642103992
Found x=939158010,y=770462691

mathe

1. 
   
2. // gcds.cpp : Defines the entry point for the console application.
   
3. //
   
4. 
   
5. #include "stdafx.h"
   
6. #include <algorithm>
   
7. #include <assert.h>
   
8. using namespace std;
   
9. #define ASSERT(x) assert(x)
   
10. int primes[]={5,7,11,13,19,23,29,31,37,41,43};
    
11. int cand[9];
    
12. int buf[4][4];
    
13. long long minr=-1;
    
14. long long mx, my;
    
15. long long gcd(long long x, long long y)
    
16. {
    
17.     if(y==0)return x;
    
18.     do{
    
19.         long long r=x%y;
    
20.         x=y;y=r;
    
21.     }while(y!=0);
    
22.     return x;
    
23. }
    
24. 
    
25. long long extGCD(long long a, long long b, long long & ar, long long & br)
    
26. {
    
27.         long long x = 0, y = 1;
    
28.     long long old_a=a,old_b=b;
    
29.         ar = 1, br = 0;
    
30.         unsigned long long t, q;
    
31.         while (b != 0)
    
32.         {
    
33.                 t = b; q = a / b; b = a % b; a = t;
    
34.                 t = x; x = ar - q * x; ar = t;
    
35.                 t = y; y = br - q * y; br = t;
    
36.         }
    
37.     if(ar<0)ar+=old_b;
    
38.     if(br<0)br+=old_a;
    
39.         return (a);
    
40. }
    
41. long long lcm(long long x, long long y)
    
42. {
    
43.     return x/gcd(x,y)*y;
    
44. }
    
45. 
    
46. long long chinese(long long c[4])
    
47. {
    
48.     long long s,t;
    
49.     long long d1=extGCD(c[0],c[1],s,t);
    
50.     ASSERT(d1==1);
    
51.     if(s!=0)s=c[1]-s;
    
52.     long long m1=c[0]*c[1];
    
53.     long long r1=c[0]*s;
    
54.     d1=extGCD(c[2],c[3],s,t);
    
55.     ASSERT(d1==1);
    
56.     long long a,b;
    
57.     a=((c[3]-3)*s)%c[3];
    
58.     b=((c[2]-2)*t)%c[2];
    
59.     long long m2=c[2]*c[3];
    
60.     long long r2=(a*c[2]+b*c[3])%m2;
    
61.     d1=extGCD(m1,m2,s,t);
    
62.     ASSERT((r2-r1)%d1==0);
    
63.     b=(r1*t)%m1;
    
64.     a=(r2*s)%m2;
    
65.     long long d2=gcd(d1,r1);
    
66.     m1/=d1;
    
67.     m2/=d2;
    
68.     d1/=d2;
    
69.     r1%=m1;r2%=m2;
    
70.     d1=extGCD(m1,m2,s,t);
    
71.     ASSERT(d1==1);
    
72.     a=(r1*t)%m1;
    
73.     b=(r2*s)%m2;
    
74.     long long r3=(a*m2+b*m1)%(m1*m2);
    
75.     return lcm(r3,d2);
    
76. }
    
77. 
    
78. void test_buf()
    
79. {
    
80.     long long c[4];
    
81.     int i;
    
82.     for(i=0;i<4;i++){
    
83.         long long l1=lcm(buf[i][0],buf[i][1]);
    
84.         long long l2=lcm(buf[i][2],buf[i][3]);
    
85.         c[i]=lcm(l1,l2);
    
86.     }
    
87.     long long x=chinese(c);
    
88.     for(i=0;i<4;i++){
    
89.         long long l1=lcm(buf[0][i],buf[1][i]);
    
90.         long long l2=lcm(buf[2][i],buf[3][i]);
    
91.         c[i]=lcm(l1,l2);
    
92.     }
    
93.     long long y=chinese(c);
    
94.     if(minr<0||minr>x+y){
    
95.         minr=x+y;
    
96.         mx=x;
    
97.         my=y;
    
98.         printf("Found x=%lld,y=%lld\n",mx,my);
    
99.     }
    
100. }
     
101. 
     
102. int _tmain(int argc, _TCHAR* argv[])
     
103. {
     
104.     int i,j,k,t;
     
105.     for(i=0;i<11;i++)for(j=i+1;j<11;j++){
     
106.         for(k=0,t=0;k<11;k++){
     
107.             if(k!=i&&k!=j){
     
108.                 cand[t++]=primes[k];
     
109.             }
     
110.         }
     
111.         do{
     
112.             buf[0][0]=6;buf[0][1]=cand[0];buf[0][2]=2;buf[0][3]=3;
     
113.             buf[1][0]=cand[1];buf[1][1]=cand[2];buf[1][2]=cand[3];buf[1][3]=cand[4];
     
114.             buf[2][0]=2;buf[2][1]=cand[5];buf[2][2]=2;buf[2][3]=cand[6];
     
115.             buf[3][0]=3;buf[3][1]=cand[7];buf[3][2]=cand[8];buf[3][3]=3;
     
116.             test_buf();
     
117.             buf[0][0]=3;buf[0][1]=2;buf[0][2]=cand[0];buf[0][3]=6;
     
118.             buf[1][0]=cand[1];buf[1][1]=cand[2];buf[1][2]=cand[3];buf[1][3]=cand[4];
     
119.             buf[2][0]=cand[5];buf[2][1]=2;buf[2][2]=cand[6];buf[2][3]=2;
     
120.             buf[3][0]=3;buf[3][1]=cand[7];buf[3][2]=cand[8];buf[3][3]=3;
     
121.             test_buf();
     
122.             buf[0][0]=3;buf[0][1]=cand[0];buf[0][2]=cand[1];buf[0][3]=3;
     
123.             buf[1][0]=cand[2];buf[1][1]=2;buf[1][2]=cand[3];buf[1][3]=2;
     
124.             buf[2][0]=cand[4];buf[2][1]=cand[5];buf[2][2]=cand[6];buf[2][3]=cand[7];
     
125.             buf[3][0]=3;buf[3][1]=2;buf[3][2]=cand[8];buf[3][3]=6;
     
126.             test_buf();
     
127.             buf[0][0]=3;buf[0][1]=cand[0];buf[0][2]=cand[1];buf[0][3]=3;
     
128.             buf[1][0]=2;buf[1][1]=cand[2];buf[1][2]=2;buf[1][3]=cand[3];
     
129.             buf[2][0]=cand[5];buf[2][1]=cand[4];buf[2][2]=cand[6];buf[2][3]=cand[7];
     
130.             buf[3][0]=6;buf[3][1]=cand[8];buf[3][2]=2;buf[3][3]=3;
     
131.             test_buf();
     
132.         }while(next_permutation(cand,cand+9));
     
133.     }
     
134.         return 0;
     
135. }
     

复制代码

mathe

上面的程序我还没有运行完毕，谁来运行一下看看？

数学星空

请问：mathe
需要下载什么程序可以运行你写的程序呢？
因为我用的是matlabR2009

mathe

一个更优的解:x=726156780,y=68838054