题目108：求解丢番图方程1/x + 1/y = 1/n

欧拉计划

Diophantine reciprocals I

In the following equation x, y, and n are positive integers.

For n = 4 there are exactly three distinct solutions:

What is the least value of n for which the number of distinct solutions exceeds one-thousand?

NOTE: This problem is an easier version of Problem 110; it is strongly advised that you solve this one first.

题目：

在如下方程中，x，y 和 n 皆为正整数。

对于 n=4，一共有三个不同的解：

要使方程拥有超过一千个不同的解，最小的 n 是多少？

注意：该题目是题目 110 的一个简单版，强烈建议你先做这一个题目。

jerryxjr1220

这题还是要找规律，通过观察前几个n（n>60），可以发现，基本上每个不同解的个数+1，n都是60的倍数。

1. n = 60
   
2. while True:
   
3.         count = 0
   
4.         for x in range(n+1,2*n+1):
   
5.                 if x*n % (x-n) == 0:
   
6.                         count += 1
   
7.         if count > 1000:
   
8.                 break
   
9.         n += 60
   
10. print (n, count)

复制代码

答案180180

guosl

这个方程可以通过代换x=n+x1，y=n+y1。得到x1*y1=n^2。所以方程的解的个数等于n^2的因数分解的个数。

永恒的蓝色梦想

guosl 发表于 2019-4-26 10:50
这个方程可以通过代换x=n+x1，y=n+y1。得到x1*y1=n^2。所以方程的解的个数等于n^2的因数分解的个数。

好像并不对哦

guosl

永恒的蓝色梦想 发表于 2019-7-24 17:35
好像并不对哦

这个方程可以通过代换x=n+x1，y=n+y1。得到x1*y1=n^2。所以方程的解的个数等于n^2的因数分解中小于等于n的因数个数。

代号-K

110 的 answer is ： 1260
由于 n 的 可拆解 为  n^2 中 两个因数的乘积
所以可以编写代码如下：

1. void getDiophantineReciprocals(int flag){
   
2.     int num = 1;
   
3.     int count;
   
4.     int i;
   
5.     ElementType numPow;
   
6.     while(true){
   
7.         count = 1;
   
8.         i = 2;
   
9.         num++;
   
10.         numPow = (ElementType) pow(num, 2);
    
11.         while(i <= num){
    
12.             if(numPow % i == 0){
    
13.                 count++;
    
14.             }
    
15.             i++;
    
16.         }
    
17.         if(count >= flag){
    
18.             printf("get it\n");
    
19.             printf("%d has %d\n", num,count);
    
20.             return;
    
21.         }
    
22.         printf("%d has %d\n",num,count);
    
23.     }
    
24. }
    
25. 
    
26. int main(int argc, char *argv[]){
    
27.     getDiophantineReciprocals(110);
    
28.     return 0;
    
29. }

复制代码

guosl

1. /*
   
2. 答案：180180
   
3. 耗时：0.0027585秒
   
4. */
   
5. #include <iostream>
   
6. #include <algorithm>
   
7. #include <omp.h>
   
8. using namespace std;
   
9. 
   
10. const int N = 1000000;
    
11. const int M = 1000;
    
12. const int nSteps = 100;
    
13. int p[N];//p[i]为i的一个素因数，用于i的素因数分解
    
14. 
    
15. int main(void)
    
16. {
    
17.   int nK = 5;
    
18.   volatile bool bContinue = true;
    
19.   int nMinK = 0x7fffffff;
    
20.   double t = omp_get_wtime();
    
21. #pragma omp parallel shared(nK,bContinue) reduction(min:nMinK)
    
22.   {
    
23.     //应用筛法求每个数的一个素因数
    
24. #pragma omp for
    
25.     for (int i = 1; i < N; ++i)
    
26.       p[i] = i;
    
27. #pragma omp single
    
28.     for (int i = 2; i <= M; ++i)
    
29.     {
    
30.       if (p[i] < i)
    
31.         continue;
    
32.       for (int j = i * i; j < N; j += i)
    
33.         p[j] = i;
    
34.     }
    
35.     while (bContinue)
    
36.     {
    
37.       int k = 0;
    
38. #pragma omp critical(c1)
    
39.       {
    
40.         k = nK;
    
41.         nK += nSteps;
    
42.       }
    
43.       for (int i = k; i < k + nSteps; ++i)
    
44.       {
    
45.         //对i^2求其因数的分解
    
46.         int nCount = 1;//记录i的因数个数
    
47.         int x = i;
    
48.         while (x > 1)
    
49.         {
    
50.           int y = 0, z = p[x];
    
51.           while ((x % z) == 0)
    
52.           {
    
53.             ++y;
    
54.             x /= z;
    
55.           }
    
56.           nCount *= (2 * y + 1);
    
57.         }
    
58.         if (nCount > 2000)//此时至少有1000个因数小于i
    
59.         {
    
60.           nMinK = min(nMinK, i);
    
61. #pragma omp critical(c2)
    
62.           bContinue = false;
    
63.           break;
    
64.         }
    
65.       }
    
66.     }
    
67.   }
    
68.   t = omp_get_wtime() - t;
    
69.   cout << nMinK << endl << t << endl;
    
70.   return 0;
    
71. }
    

复制代码