题目44：找出最小的和与差都是五角数的五角数对

欧拉计划

Pentagon numbers

Pentagonal numbers are generated by the formula, . The first ten pentagonal numbers are:

1, 5, 12, 22, 35, 51, 70, 92, 117, 145, ...

It can be seen that . However, their difference, 70 − 22 = 48, is not pentagonal.

Find the pair of pentagonal numbers, , for which their sum and difference are pentagonal and is minimised; what is the value of D?

题目：

五角数通过如下公式定义： 。前十个五角数是：

1, 5, 12, 22, 35, 51, 70, 92, 117, 145, ...

可以看出 . 但是它们的差 70 − 22 = 48 却不是五角数。

找出最小的五角数对 ， 使得它们的和与差都是五角数，并且 取到最小。这时 D 的值是多少?

迷雾少年

看不懂。。

jerryxjr1220

Find: P1:7042750,P2:1560090,D:5482660,S:8602840

1. 
   
2. def isRoot(x):        #根的求解方程，判别是否是整数
   
3.     if (1+(1+24*x)**0.5) % 6 == 0:
   
4.         return 1
   
5.     else:
   
6.         return 0
   
7. 
   
8. for i in range(1,10000):
   
9.     P1 = int(i*(3*i-1)/2)        
   
10.     for j in range(i,1,-1):
    
11.         P2 = int(j*(3*j-1)/2)        
    
12.         if isRoot(P1+P2) and isRoot(P1-P2):
    
13.             print ('Find: P1:%d,P2:%d,D:%d,S:%d' % (P1,P2,P1-P2,P1+P2))
    
14.             exit()
    

复制代码

jerryxjr1220

利用set（）元组去搜索，更快速 0.3 sec出结果

1. from time import time
   
2. start = time()
   
3. 
   
4. def generatePentagonals(n):
   
5.     pentagonals = [0]
   
6.     x=1
   
7.     while pentagonals[len(pentagonals)-1] < n:
   
8.         pentagonals.append(x*(3*x-1)/2)
   
9.         x+=1
   
10.     pentagonals.pop(0)
    
11.     return pentagonals
    
12. 
    
13. 
    
14. def findPair(p):
    
15.     pSet = set(p)
    
16. 
    
17.     for i in range (1,len(p)):
    
18.         for j in range(0,i):
    
19.             if (p[i]-p[j] in pSet) and (p[i]+p[j] in pSet):
    
20.                 return p[i] - p[j]
    
21.                
    
22. 
    
23.     return "Pair not in range"
    
24. 
    
25.    
    
26. print (findPair(generatePentagonals(10000000)))
    
27. 
    
28. print("Time: {0} secs".format(time()-start))

复制代码

飘飞的白杨

1. def progrem44():
   
2.     s = set(i*(3*i-1)//2 for i in range(1,2501))
   
3.     for i in s:
   
4.         for j in s:
   
5.             if i < j: continue
   
6.             if i+j in s and i+2*j in s:
   
7.                 return j,i+j
   
8. print(progrem44())

复制代码

芒果加黄桃

没有好办法，学习4楼老大代码

愤怒的大头菇

1. def Pentagon(n):
   
2.     numbers = n*(3*n-1)/2
   
3.     return numbers
   
4. 
   
5. 
   
6. def Judg(x, i):
   
7.     for n in range(1, 2*i-1):
   
8.         if x == n*(3*n-1)/2:
   
9.             return True
   
10.     return False
    
11. 
    
12. 
    
13. i = 2
    
14. temps = True
    
15. 
    
16. while temps:
    
17.     for n in range(1, i):
    
18.         temp = Pentagon(i) - Pentagon(n)
    
19.         tmp = Pentagon(i) + Pentagon(n)
    
20.         if Judg(temp, i) and Judg(tmp, i):
    
21.             print(temp)
    
22.             temps = False
    
23. 
    
24.     i += 1

复制代码

能算出来，但是速度很慢

渡风

%% Problem44.m
%最后编辑时间：2017-06-13 9:21 版本1.0
%找出 两个五角星数，使得它们的和为五角星数，差为五角星数且为最小
% Problem44所用时间为28.294秒
% Problem44的答案为5482660

1. %% Problem44.m
   
2. %最后编辑时间：2017-06-13 9:21 版本1.0
   
3. %找出 两个五角星数，使得它们的和为五角星数，差为五角星数且为最小
   
4. % Problem44所用时间为28.294秒
   
5. % Problem44的答案为5482660
   
6. function Output = Problem44()
   
7. tic
   
8. 
   
9. Start = 3;
   
10. Meet = 0;
    
11. while (Meet == 0)
    
12.     Bignum = (Start * (3*Start - 1))/2;
    
13.     for ii = Start - 1: -1 : 1
    
14.         Smallnum = (ii * (3*ii - 1))/2;
    
15.         if IsPentagon(Bignum - Smallnum) == 1 && IsPentagon(Bignum + Smallnum) == 1
    
16.             Meet = 1;
    
17.             break
    
18.         end
    
19.     end
    
20.     Start = Start + 1;
    
21.     disp(Start)
    
22. end
    
23. 
    
24. Output = Bignum - Smallnum;
    
25. toc
    
26. disp('此代码使用matlab编程')
    
27. disp(['Problem44所用时间为',num2str(toc),'秒'])
    
28. disp(['Problem44的答案为',num2str(Output)])
    
29. end
    
30. %% IsPentagon.m
    
31. % 测试一个数是否为五角星
    
32. % 最后编辑时间：17-06-13 9:03 版本：1.0
    
33. % 是五角型数则返回 1，否则返回 0
    
34. function Judge = IsPentagon(n)
    
35. if nargin == 0
    
36. n = 22;
    
37. end
    
38. if n == 1
    
39.     Judge = 1;
    
40. else
    
41.     front = floor(sqrt(2*n)/2);
    
42.     after = ceil(sqrt(2*n)/1.7);
    
43.     Judge = 0;
    
44.     for ii = after:-1:front
    
45.         if n == ((3*ii - 1) * ii)/2
    
46.             Judge = 1;
    
47.             break
    
48.         end
    
49.     end
    
50. end
    
51. end

复制代码

cq_xuke

1. import time
   
2. s=time.clock()
   
3. from itertools import takewhile
   
4. def pentagon(init=None):
   
5.     if init is None:
   
6.         init=1
   
7.     while not ispentagon(init):
   
8.         init+=1
   
9.     n=((1+24*init)**0.5+1)//6
   
10.     while True:
    
11.         yield n*(3*n-1)//2
    
12.         n+=1
    
13. 
    
14. def ispentagon(n):
    
15.     flag=False
    
16.     if int(((1+24*n)**0.5+1)/6)==((1+24*n)**0.5+1)/6:
    
17.         flag=True
    
18.     return flag
    
19. 
    
20. flag=False
    
21. for p1 in pentagon(3):
    
22.     if flag==True:
    
23.         break
    
24.     for p2 in takewhile(lambda x:x<p1,pentagon(1)):
    
25.         if ispentagon(p1+p2) and ispentagon(p1-p2):
    
26.             print(int(p1-p2))
    
27.             flag=True
    
28.             break
    
29. print('Runing time: {0}s'.format(round(time.clock()-s,3)))

复制代码

5482660
Runing time: 3.316s

najin

用的matlab
结果是：
>> Untitled
     7042750

     1560090

     5482660

时间已过 0.167725 秒。
>>

王小召

最小值是：442.0
用时：0.7488047999999999秒

import math
import time


def pentagon_num(num):
    if int(math.sqrt(24*num + 1)) == math.sqrt(24*num + 1):
        return True


def cal_result():
    result = []
    #  2000的五角数范围已经达到 1000 * (2999)/2 = 1499500
    for i in range(2, 1000):
        for j in range(1, i):
            if pentagon_num(i*(3*i-1)/2 + j*(3*j-1)/2) and pentagon_num(i*(3*i-1)/2 - j*(3*j-1)/2):
                result.append(i*(3*i-1)/2 - j*(3*j-1)/2)
    return result

print("最小值是：{}\n用时：{}秒".format(min(cal_result()), time.process_time()))

永恒的蓝色梦想

抄了一下3L的求根公式：

1. isRoot=lambda x:not(1+(1+24*abs(x))**0.5)%6
   
2. pentagon=lambda x:(3*x-1)*x//2
   
3. for k in range(1,10000):
   
4.   i=pentagon(k)
   
5.   for j in range(k,10000):
   
6.     j=pentagon(j)
   
7.     if isRoot(i-j)and isRoot(i+j):
   
8.       print(f'运行结果：{abs(i-j)}\n运行时间：{t()-t1}s')
   
9.       exit()

复制代码

debuggerzh

5482660

Process returned 0 (0x0)   execution time : 1.176 s
Press any key to continue.
纯暴力解法。不妨设Pj>Pk
先枚举D，再枚举Pk，二分查找
D小于P[k+1] - P[k]时枚举下一个D

1. #include<iostream>
   
2. #include<algorithm>
   
3. using namespace std;
   
4. 
   
5. const int M = 3e3;
   
6. int p[M] = {0,1};
   
7. 
   
8. void ini(){
   
9.   for (int i = 2;i < M;i++){
   
10.     p[i] = p[i-1] + 3*i - 2;
    
11.   }
    
12. }
    
13. 
    
14. void print(int n){
    
15.   for (int i = 1;i < n;i++){
    
16.     cout << p[i] << endl;
    
17.   }
    
18. }
    
19. 
    
20. int main(){
    
21.   ini();
    
22. 
    
23.   for (int a = 2;;a++){
    
24.     for (int k = 1;k+1 < M && p[k+1] - p[k] <= p[a];k++){
    
25.       if (!binary_search(p, p+M, p[a] + p[k]))  continue;
    
26.       if (!binary_search(p, p+M, p[a] + 2*p[k]))  continue;
    
27.       cout << p[a] << endl;
    
28.       return 0;
    
29.     }
    
30.   }
    
31. 
    
32.   return 0;
    
33. }
    

复制代码

a1351468657

1. #include <stdio.h>
   
2. #include <math.h>
   
3. 
   
4. int check_Pentagonal(int);
   
5. int check_Pentagonal(int num)
   
6. {
   
7. 
   
8.         if (num == 0)
   
9.         {
   
10.                 return 0;
    
11.         }
    
12.         float j;
    
13.         int i;
    
14.         j = (sqrt(24 * num + 1) + 1) / 6;
    
15.         i = (sqrt(24 * num + 1) + 1) / 6;
    
16. 
    
17.         if ((float)i == j)
    
18.         {
    
19.                 return 1;
    
20.         }
    
21.         else
    
22.         {
    
23.                 return 0;
    
24.         }
    
25. }
    
26. 
    
27. main()
    
28. {
    
29.         int i, j, m, n;
    
30.        
    
31.         for (i = 2; i < 10000; i++)
    
32.         {
    
33.                 for (j = i; j >= 1; j--)
    
34.                 {
    
35.                         m = i * (i * 3 - 1) / 2;
    
36.                         n = j * (j * 3 - 1) / 2;
    
37.                         if (check_Pentagonal(m + n))
    
38.                         {
    
39.                                 if (check_Pentagonal(m - n))
    
40.                                 {
    
41.                                         goto Label;
    
42.                                 }
    
43.                         }
    
44.                 }
    
45.         }
    
46. Label:printf("Pj = %d, Pk = %d D = %d", m, n, m - n);
    
47. }

复制代码

Pj = 1560090, Pk = 7042750 D = 5482660

guosl

1. /*
   
2. 答案：5482660
   
3. 耗时：2.34444秒（单核）
   
4.       0.198174秒（八核）
   
5. */
   
6. #include <iostream>
   
7. #include <vector>
   
8. #include <algorithm>
   
9. #include <omp.h>
   
10. 
    
11. using namespace std;
    
12. 
    
13. int main(void)
    
14. {
    
15.   //打表
    
16.   vector<int> v;
    
17.   int k = 1;
    
18.   while (true)
    
19.   {
    
20.     int p = k * (3 * k - 1) / 2;
    
21.     if (p < 1000000000)
    
22.     {
    
23.       v.push_back(p);
    
24.       ++k;
    
25.     }
    
26.     else
    
27.       break;
    
28.   }
    
29.   int nMinD = 0x7fffffff, nK = int(v.size());
    
30.   volatile int i = 0;
    
31.   volatile bool bContinue = true;
    
32. #pragma omp parallel firstprivate(nK) shared(i,v,bContinue) reduction(min:nMinD)
    
33.   while (bContinue && i < nK)
    
34.   {
    
35.     int nD;
    
36. #pragma omp critical(c1)
    
37.     nD = v[i++];//枚举差
    
38.     bool bFind = false;
    
39.     for (int j = 0; j < nK; ++j)
    
40.     {
    
41.       auto itr = lower_bound(v.begin(), v.end(), nD + v[j]);
    
42.       if (itr != v.end() && *itr == nD + v[j])//找到两个五边形数使其差位nD
    
43.       {
    
44.         int k = int(itr - v.begin());
    
45.         int nD1 = v[k] + v[j];
    
46.         itr = lower_bound(v.begin(), v.end(), nD1);
    
47.         if (itr != v.end() && *itr == nD1)//找到两个五边形数使其和与差都是五边形数
    
48.         {
    
49.           bFind = true;
    
50.           break;
    
51.         }
    
52.       }
    
53.     }
    
54.     if (bFind)
    
55.     {
    
56.       nMinD = min(nD, nMinD);//记录最后结果
    
57. #pragma omp critical(c2)
    
58.       bContinue = false;
    
59.     }
    
60.   }
    
61.   cout << nMinD << endl;
    
62.   return 0;
    
63. }

复制代码

B1tetheDust

1. import time as t
   
2. import numpy as np
   
3. 
   
4. start = t.perf_counter()
   
5. 
   
6. 
   
7. def is_integer(a_num):
   
8.     delta = 0.00001
   
9.     if abs(a_num - int(a_num)) < delta:
   
10.         return True
    
11. 
    
12. 
    
13. def find_pentagon():
    
14.     given_range = 3000
    
15.     for j in range(1000, given_range):
    
16.         for k in range(j + 2, given_range + 2):
    
17.             k_add_j = np.sqrt(np.square(k) + np.square(j) - 1 / 3 * (k + j) + 1 / 36) + 1 / 6
    
18.             k_minus_j = np.sqrt(np.square(k) - np.square(j) - 1 / 3 * (k - j) + 1 / 36) + 1 / 6
    
19.             if is_integer(k_add_j) and is_integer(k_minus_j):
    
20.                 return [j, k, k_minus_j]
    
21. 
    
22. 
    
23. nums = find_pentagon()
    
24. print(nums[2] * (3 * nums[2] - 1) / 2)
    
25. print("It costs %f s" % (t.perf_counter() - start))
    

复制代码

5482660.0
It costs 0.448183 s