18 Dec

迟到一年的建模:再探碎纸复原

前言:一年前国赛的时候,很初级地做了一下B题,做完之后还写了个《碎纸复原:一个人的数学建模》。当时就是对题目很有兴趣,然后通过一天的学习,基本完成了附件一二的代码,对附件三也只是有个概念。而今年我们上的数学建模课,老师把这道题作为大作业让我们做,于是我便再拾起了一年前的那份激情,继续那未完成的一个人的数学建模...

与去年不同的是,这次将所有代码用Python实现了,更简洁,更清晰,甚至可能更高效~~以下是论文全文。

研究背景

2011年10月29日,美国国防部高级研究计划局(DARPA)宣布了一场碎纸复原挑战赛(Shredder Challenge),旨在寻找到高效有效的算法,对碎纸机处理后的碎纸屑进行复原。[1]该竞赛吸引了全美9000支参赛队伍参与角逐,经过一个多月的时间,有一支队伍成功完成了官方的题目。

近年来,碎纸复原技术日益受到重视,它显示了在碎片中“还原真相”的可能性,表明我们可以从一些破碎的片段中“解密”出原始信息来。另一方面,该技术也和照片处理领域中的“全景图拼接技术”有一定联系,该技术是指通过若干张不同侧面的照片,合成一张完整的全景图。因此,分析研究碎纸复原技术,有着重要的意义。

阅读剩余部分...

28 Oct

在Python中使用GMP(gmpy2)

之前笔者曾写过《初试在Python中使用PARI/GP》,简单介绍了一下在Python中调用PARI/GP的方法。PARI/GP是一个比较强大的数论库,“针对数论中的快速计算(大数分解,代数数论,椭圆曲线...)而设计”,它既可以被C/C++或Python之类的编程语言调用,而且它本身又是一种自成一体的脚本语言。而如果仅仅需要高精度的大数运算功能,那么GMP似乎更满足我们的需求。

了解C/C++的读者都会知道GMP(全称是GNU Multiple Precision Arithmetic Library,即GNU高精度算术运算库),它是一个开源的高精度运算库,其中不但有普通的整数、实数、浮点数的高精度运算,还有随机数生成,尤其是提供了非常完备的数论中的运算接口,比如Miller-Rabin素数测试算法、大素数生成、欧几里德算法、求域中元素的逆、Jacobi符号、legendre符号等[来源]。虽然在C/C++中调用GMP并不算复杂,但是如果能在以高开发效率著称的Python中使用GMP,那么无疑是一件快事。这正是本文要说的gmpy2

阅读剩余部分...

17 Oct

两百万素数之和与“电脑病”

原则上来讲,同样的算法,如果分别在Python和C++上实现,那么Python的速度肯定比不上C++的。但是Python还被称为“胶水语言”,它允许我们把主要计算的部分用C或C++等高效的语言编写好,然后它作为“粘合剂”把两者粘合在一起。正因为如此,Python才有了各种各样的扩展库,这些库中有不少是用C语言编写的。因此,我们在编写Python程序的时候,如果可以用这些现成的库,速度会快很多。本文就是用Numpy来改进之前的《两百万前素数之和与前两百万素数之和》的计算。

算法本身是没有变的,只是用了Numpy来处理数组计算,代码如下:

阅读剩余部分...

7 Oct

班门弄斧:Python的代码能有多简洁?

此文很水,高手略过...

Python以它的开发效率而闻名,优秀的开发效率自然意味着它能够用更少的代码实现更多的功能。那么,对于同一个问题,Python的代码能有多简洁?而我们怎么平衡开发效率和运行效率?笔者学了几个月Python,略懂一点,在此班门弄斧一翻。

在此,我们来编程计算
$$\sum_{n=0}^{10} n^2$$
这当然是一道非常简单的习题。按照一般思路,写出来的最自然的代码就是:

s=0
for i in range(11):
  s=s+i**2
print(s)

阅读剩余部分...

22 Jul

初试在Python中使用PARI/GP

BoJone很喜欢Python,也很喜欢数论,所以就喜欢利用Python玩数论了。平时也喜欢自己动手写一些数论函数,毕竟Python支持大整数高精度运算,这点是非常好的;但是,在很多实际应用中,还是希望能有一个现成的数论函数库来调用。之前尝试过数学研发网的HugeCalc库,但是由于各种不熟悉不了了之。后来论坛上的无心老兄推荐了PARI/GP,小试一下,居然在Python上成功调用了。以后再也不用担心Python上的数论计算问题了,呵呵~

阅读剩余部分...

27 Jun

Project Euler 454 :五天攻下“擂台”

进入期末了,很多同学都开始复习了,这学期我选的几门课到现在还不是很熟悉,本想也在趁着这段时间好好看看。偏生五天前我在浏览数学研发论坛的编程擂台时看到了这样的一道题目

设对于给定的$L$,方程
$$\frac{1}{x}+\frac{1}{y}=\frac{1}{n}$$
满足$0 < x < y \leq L$的正整数解共有$f(L)$种情况。比如$f(6)=1,f(12)=3,f(1000)=1069$,求$f(10^{12})$。

这道题目的来源是Project Euler的第454题:Diophantine reciprocals III(丢潘图倒数方程),题目简短易懂,但又不失深度,正符合我对理想题目的定义。而且最近在学习Python学习得不亦乐乎,看到这道题目就跃跃欲试。于是乎,我的五天时间就没有了,而且过程中几乎耗尽了我现在懂的所有编程技巧。由于不断地测试运行,我的电脑发热量比平时大了几倍,真是辛苦了我的电脑。最后的代码,自我感觉已经是我目前写的最精彩的代码了。在此与大家共享和共勉~

上述表达式是分式,不利于编程,由于$n=\frac{xy}{x+y}$,于是上述题目也等价于求$(x+y)|xy$(意思是$x+y$整除$xy$)的整数解。

阅读剩余部分...

11 Jun

用PyPy提高Python脚本执行效率

《两百万前素数之和与前两百万素数之和》中,我们用Python求了前两百万的素数和以及两百万前的素数和,并且得到了在Python 3.3中的执行时间如下:

两百万前的素数之和:
142913828922
time: 2.4048174478605646

前两百万的素数之和:
31381137530481
time: 46.75734807838953

于是想办法提高python脚本的执行效率,我觉得在算法方面,优化空间已经比较小了,于是考虑执行器上的优化。在搜索的无意间我看到了一个名词——Psyco!这是python的一个外部模块,导入后可以加快.py脚本的执行。网上也有《用 Psyco 让 Python 运行得像 C一样快》、《利用 psyco 让 Python 程序执行更快》之类的文章,说明Psyco确实是一个可行的选择,于是就跃跃欲试了,后来了解到Psyco在2012年已经停止开发,只支持到Python 2.4版本,目前它由 PyPy所接替。于是我就下载了PyPy

阅读剩余部分...

10 Jun

两百万前素数之和与前两百万素数之和

标题说了两道比较好玩的编程题,如果读者觉得标题绕的让人眩晕的话,那么让我再说得清晰一点:

两百万前素数之和指的是所有不超过两百万的素数的和;
前两百万素数之和指的是前两百万个素数的和。

我是从子谋的blog中看到这道题目的,前一道题目是Project Euler的第10题,后一道则是我跟子谋探索着玩的。关于子谋的研究和代码,大家可以去他的blog上学习。本文分享一下我自己的想法。

阅读剩余部分...