寻找雷劈数

2019-10-02  本文已影响0人  丁嵩冰

戴公子预初,看到个故事:
有位数学家卡普利加在一次旅行中,遇到猛烈的暴风雨,电闪雷鸣过后,他看到路边一块里程碑,被雷电劈成两半,一半上刻著30,另一半刻著25。这时,卡普利加的脑际中忽然发现了一个绝妙的数学关系——把劈成两半的数加起来,再平方,正好是原来的数位。除此之外,还有没有别的数,也具有这样的性质呢?
熟悉速算的人很快就找到了另一个数:2025 按照第一个发现者的名字,这种怪数被命名为“卡普利加数”,又称“雷劈数”。
据说这样的书有无穷多个,怎么找到他们呢?
----------------------------------------------------------------------------------------------->
首先想到的是excel,用row()^2,用Log找到结果是偶数位的行,对半拆开,相加平方后比较,轻松找到十来个五位数以内的对劈的雷劈数。
不过,戴公子发现在excel中公式拖到底运算量太大电脑差点死机,总共能算出几十个,再多用excel就不行了。
到此为止已经可以证明有无穷多个雷劈数了:N个9,接一个8,接N个0,再接一个1,肯定是是N+1个9的平方。
如果包含非对劈的,1后面接偶数个0的都是。
其实也不是excel不行,算法得改进,不过要更好解决,还是交给python吧。
---------------------------------------------------------------------------------------------->

import math
for n in range(1,9):
    L=pow(10,n-1)
    H=pow(10,n)
    for a in range(L,H):
        if a%10 in [0,4,5,8]: #加014569后还是完全平方的数肯定以0458结尾
            for b in range(1,H):
                if b%10 in [0,1,4,5,6,9]: # 完全平方数肯定以0124569结尾
                    if pow(a+b,2)==a*pow(10,n)+b:
                        print(a,b,a*pow(10,n)+b)

转啊转算出十来个,好像比较慢诶......
作为解释型语言,python确实比较慢,不过...我们有numba!
---------------------------------------------------------------------------------------------->

from numba import jit
import math

@jit(nopython=True,fastmath = True  )
def dai():
    for n in range(1,18):
        L=pow(10,n-1)
        H=pow(10,n)
        for aa in range(L,H,10):
            for k in [0,4,8]: #简单分析可知a肯定是048结尾,运算量一下子少了70%
                a=aa+k
                rootp=math.sqrt(H)
                Lb=int(math.sqrt(a)*rootp)-a #关键优化2,只需遍历几个可能的b
                Hb=int(math.sqrt(a+1)*rootp)-a+1
                for b in range(Lb,Hb):
                    t1=b%100 #b的结尾只有几种可能性,运算量再减少近八成
                    t2=pow(a%100 + t1,2)%100 #同余
                    if t1==t2 and (t1 in [1,4,9,16,25,36,49,64,81,0,21,44,69,96,56,89,24,61,41,84,29,76]):
                        if pow(a+b,2)==a*H+b:
                            m=int(n-math.log10(b+1)) #补0,优化输出格式
                            print('a=',a,'\t\tb=',b,'\ta+b=',a+b,'\t\tab=',a,'0'*m,b)

if __name__ == '__main__':
    dai()

结果是竹筒倒豆子,噼里啪啦输出一大堆,速度提高100倍都不止
网上查了一圈,好多种不同编程语言都有找雷劈数的代码,但基本上都是通过第一种算法,找到万亿级后速度就很慢了。
利用python的无限制长整型找雷劈数,再用numba优化编译成机器码加速,轻松找到亿亿亿级雷劈数。

补充:本来想用取模平方再取模的同余算法优化判断长整型平方效率,debug过程中发现,jit后对大整数会溢出,好在平方和ab同样溢出,不影响等值判断,但没法用同余算法优化。

睡一觉起来,笔记本i7-6700HQ已经得到了200来个对劈数,最大的是这个28位数,约3809亿亿亿:
ab= 38098433988110 23625494911556

上一篇下一篇

猜你喜欢

热点阅读