找回密码
 欢迎注册
楼主: wayne

[擂台] projectEuler 454 -- Diophantine reciprocals

[复制链接]
发表于 2014-4-18 09:42:14 | 显示全部楼层
几何意义: 1/A+1/B=1/C
要是有人能把这个图片用latex实现,那就把这个图片删除了吧,
不过删之前先发短信息告诉我一下


QQ截图20140418093547.png
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-4-18 09:46:01 | 显示全部楼层
几何意义就是三个纵坐标都是整数,然后还满足上面的图形要求
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-5-23 13:17:48 | 显示全部楼层
因为
\(\frac1x+\frac1y=\frac1n\),\(0\lt x\lt y\leq L\)
所以 \(n\lt x\lt y\leq L\)
设 \(x=n+s\),\(y=n+t\),则\(1\leq s\lt t\leq L-n\)
\(\frac1{n+s}+\frac1{n+t}=\frac1{n} \iff n^2=st\)
所以
\(n^2=st\geq1\times 2=2 \Rightarrow n\geq 2\)
\(s=n^2/t\geq \frac{n^2}{L-n}\)
\(\frac{n^2}{L-n}\leq s\lt n\lt t\leq L-n\)
\(\frac1{n}=\frac1x+\frac1y\gt \frac1L+\frac1L=\frac2L \Rightarrow n\lt \frac{L}{2}\)
所以
\(f(L)\)
\(=\sum_{n=2}^{\lfloor \frac{L-1}{2} \rfloor} \sum_{\begin{subarray}{cc} s\geq\frac{n^2}{L-n} \\ s\,|\, n^2 \end{subarray}}^{n-1} 1\)
\(=\sum_{n=2}^{\lfloor \frac{L-1}{2} \rfloor} \sum_{\begin{subarray}{cc} t=n+1 \\ t\,|\, n^2 \end{subarray}}^{L-n} 1\)
\(=\sum_{n=2}^{\lfloor \frac{L-1}{2} \rfloor} ((\sum_{\begin{subarray}{cc} u\geq\frac{n^2}{L-n} \\ u\,|\, n^2 \end{subarray}}^{L-n} 1)-1)/2\)

关键的是
给定 \(n\),\(L\),求满足
\(t\,|\, n^2\),\(n\lt t\leq L-n\)
的\(t\)的数目
这一步的复杂度与\(n^2\)的质因数分解相同

评分

参与人数 1威望 +3 金币 +3 贡献 +3 经验 +3 鲜花 +3 收起 理由
wayne + 3 + 3 + 3 + 3 + 3 赞一个!

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-23 22:41:48 | 显示全部楼层
可否使用多线程?
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
 楼主| 发表于 2014-6-23 23:24:53 | 显示全部楼层

这个是求和性质的,所以可以使用多线程提升速度。但也只是常数优化。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-24 08:38:38 | 显示全部楼层
算这个的精确值困难的话,算算估计值呢。就像n!   当然有公式算精确值最好了
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-24 13:00:28 | 显示全部楼层
sunwukong 发表于 2014-5-23 13:17
因为
\(\frac1x+\frac1y=\frac1n\),\(0\lt x\lt y\leq L\)
所以 \(n\lt x\lt y\leq L\)


其实就是涉及大数分解,虽然大数分解本身没有什么高效率的算法,但是这个题目只是统计解的个数,而非找到每一组解。
联想到给出不大于`x`的素数个数`\pi(x)`的算法目前最高效的算法中,并不需要找出所有小于x的素数然后统计个数(见M. DELEGLISE AND J. RIVAT.COMPUTING  (x): THE MEISSEL, LEHMER, LAGARIAS,MILLER, ODLYZKO METHOD)
因而,`n^2`满足一定限制要求的分解方案数,应该有相应不需要找出每种具体方案的高效算法。
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-25 08:15:54 | 显示全部楼层
还有种几何关系满足这个等式:平面过一点的三射线中间射线与两边射线夹角都是60度。一条直线截这3射线为三段,这三段长度满足这个等式

点评

大哥,画张图上来吧  发表于 2014-6-25 11:52
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-25 19:45:24 | 显示全部楼层
本帖最后由 282842712474 于 2014-6-25 20:02 编辑

有没有已知的数据进行对比?
我用python算得$f(10^7)=30093331,\quad f(10^8)=349446716$,前者用时1.7秒,后者用时17秒。
照$O(N)$复杂度算,$f(10^12)$需要200000秒,也就是六七个小时左右(实际可能更多)。使用多个线程,可以把时间减半吧~

我是直接用
http://bbs.emath.ac.cn/forum.php ... 2667&fromuid=70
的第一个公式算的。

补充:
算了$f(10^9)=3979600400$,用时196.748994145秒...

补充内容 (2014-6-26 08:38):
我看错了...200000秒是差不多3天了~
毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
发表于 2014-6-25 19:49:14 | 显示全部楼层
本帖最后由 282842712474 于 2014-6-25 20:06 编辑
282842712474 发表于 2014-6-25 19:45
有没有已知的数据进行对比?
我用python算得$f(10^7)=30093331,\quad f(10^8)=349446716$,前者用时1.7秒 ...


下面的代码在python2.7和python3均可以运行,17秒的结果是在pypy下运行得到
由于python是号称运行速度最慢的编程语言,如果转为C++,应该有优化空间。但是我不知道C++上面如何处理大整数($10^12$对于C的整数型来说已经溢出了)。


  1. import time
  2. start=time.clock()
  3. import math
  4. l=10000000 #上限
  5. s=0 #统计个数
  6. rl=int(math.sqrt(l))
  7. rrl=int(math.sqrt(int(math.sqrt(l))))

  8. #弄出个素数表来
  9. prime=[i for i in range(1,rl+1)]#定义整数表
  10. j=2
  11. while j<=rrl:
  12.   if prime[j-1] != 0:
  13.     m=int(rl/j)
  14.     for i in range(j,m+1):
  15.       prime[i*j-1]=0
  16.     j=j+1
  17.   else:
  18.     j=j+1
  19. prime.sort() #从小到大重新排列(让0位于前面)
  20. z=prime.count(0) #统计0的个数
  21. prime=prime[z+1:rl]
  22. #素数表生成完毕

  23. for m in range(1,int(rl/1.4)):
  24.   #找出与m互素的数,存在p里
  25.   p=[i for i in range(1,rl+1)]
  26.   m2=m
  27.   for i in prime:
  28.     if i>m:
  29.       break
  30.     if m%i==0:
  31.       m=m/i
  32.       for j in range(m2//i+1,int(rl/i)+1): #!!!
  33.         p[i*j-1]=0
  34.   p.sort() #从小到大重新排列(让0位于前面)
  35.   z=p.count(0) #统计0的个数
  36.   p=p[z+m2:rl]
  37.   #找出与m互数的完毕。
  38.   for n in p:
  39.     s=s+l//((n+m2)*n)
  40. print(s)
  41. end=time.clock()
  42. print(end-start)
复制代码

评分

参与人数 1威望 +5 金币 +5 贡献 +5 收起 理由
wayne + 5 + 5 + 5 赞一个!

查看全部评分

毋因群疑而阻独见  毋任己意而废人言
毋私小惠而伤大体  毋借公论以快私情
您需要登录后才可以回帖 登录 | 欢迎注册

本版积分规则

小黑屋|手机版|数学研发网 ( 苏ICP备07505100号 )

GMT+8, 2024-3-29 06:52 , Processed in 0.053760 second(s), 22 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表