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$)的整数解。

点击阅读全文...

27 Jan

三个相切圆的公切圆

在学车的时候,我堂大哥曾问我一道作圆的问题:

三圆的外切圆和内切圆 (1)

三圆的外切圆和内切圆 (1)

平面上给出三个两两相切的圆以及它们的圆心,求作一个圆与这三个圆都相切(尺规作图)。

如果从纯几何的途径入手,我们甚至很难判断这样的圆是否存在。但是我之前似乎已经看过类似的题目,于是很快想到一个名词:反演。反演可以将圆反演成直线(圆过反演点),也可以将圆反演成圆(圆不过反演点),而其他的相切、相交等关系保持不变。对反演后的图形进行相同的反演,就变回原来的图形。本题的难点在于圆太多,利用反演,我们可以将它变为两条直线和一个圆的问题。

假设读者已经有了反演的基本知识,如果没有,请到
http://zh.wikipedia.org/wiki/反演

阅读相关内容。

点击阅读全文...

18 Jun

线性微分方程组:已知特解求通解

含有$n$个一阶常微分方程的一阶常微分方程组
$$\dot{\boldsymbol{x}}=\boldsymbol{A}\boldsymbol{x}$$
其中$\boldsymbol{x}=(x_1(t),\dots,x_n(t))^{T}$为待求函数,而$\boldsymbol{A}=(a_{ij}(t))_{n\times n}$为已知的函数矩阵。现在已知该方程组的$n-1$个线性无关的特解$\boldsymbol{x}_1,\boldsymbol{x}_2,\dots,\boldsymbol{x}_{n-1}$(解的列向量),求方程的通解。

这是我的一位同学在6月5号问我的一道题目,我当时看了一下,感觉可以通过李对称的方法很容易把解构造出来,当晚就简单分析了一下,发现根据李对称的思想,由上面已知的信息确实足以把通解构造出来。但是我尝试了好几天,尝试了几何、代数等思想,都没有很好地构造出相应的正则变量出来,从而也没有写出它的显式解,于是就搁置下来了。今天再分析这道题目时,竟在无意之间构造出了让我比较满意的解来~

点击阅读全文...

25 Feb

翻到新的维度,把积分解决!

一般来说,如果原函数容易找到的话,牛顿-莱布尼兹公式是定积分的通用方法。但是牛顿-莱布尼兹公式只适合连续函数的积分,如果积分区间含有奇点,那就不成立了。比如,我们考虑积分
$$\int_{-1}^1 \frac{1}{x^2}dx$$
当然,从严格的数学上来说,这种写法是不成立的,因为被积函数在原点没有意义。当然,从物理的角度来考虑,由于对称性,我们确信
$$\int_{-1}^1 \frac{1}{x^2}dx=2\int_{0}^1 \frac{1}{x^2}dx=\lim_{\varepsilon\to 0}2\int_{\varepsilon}^1 \frac{1}{x^2}dx$$
从而得出积分发散的结论。这种处理某种程度上是可以接受的,但是却不是让人满意的,因为它导致了分段。有什么办法可以直接处理这种情况呢?确实有的,同样引入参数,并且最终让参数为0,考虑带参数的积分
$$\int_{-1}^1 \frac{1}{x^2+\varepsilon^2}dx$$
只要参数为正,这个被积函数就在$\mathbb{R}$上处处连续了,也就是奇点消失了,这样子就可以用牛顿-莱布尼兹公式了
$$\int_{-1}^1 \frac{1}{x^2+\varepsilon^2}dx=\left.\frac{1}{\varepsilon}\arctan\left(\frac{x}{\varepsilon}\right)\right|_{-1}^{1}$$
考虑$\varepsilon\to 0$的情况,就自动得到了积分发散的结论。

点击阅读全文...

18 Mar

倒立单摆之分离频率

Mathieu方程

在文章《有质动力:倒立单摆的稳定性》中,我们分析了通过高频低幅振荡来使得倒立单摆稳定的可能性,并且得出了运动方程
$$l\ddot{\theta}+[h_0 \omega^2 \cos(\omega t)-g]\sin\theta=0$$

由此对单摆频率的下限提出了要求$\omega \gg \sqrt{\frac{g}{h_0}}$。然而,那个下限只不过是必要的,却不是充分的。如果要完整地分析该单摆的运动方程,最理想的方法当然是写出上述常微分方程的解析解。不过很遗憾,我们并没有办法做到这一点。我们只能够采取各种近似方法来求解。近似方法一般指数值计算方法,然后笔者偏爱的是解析方法,也就是说,即使是近似解,也希望能够求出近似的解析解。

点击阅读全文...

19 Mar

一本对称闯物理:相对论力学(一)

简单说说

《可畏的对称》

《可畏的对称》

笔者最近陶醉于从李对称的角度来理解力学和场论,并且计算得到一些比较有趣的结果,遂想在此与大家分享。我发现,仅仅需要一个描述对称的无穷小生成元和一些最基本的假设,几乎就可以完成地推导出整个相对论力学来,甚至推导出整个(经典)场论理论来。这确实是不可思议的,我现在能基本体会到当年徐一鸿大师写的《可畏的对称》的含义了。对称的威力如此之大,以至于我们真的不得不敬畏它。而在构思本文题目的时候,我也曾想到过用“可畏的对称”为题,但不免有抄袭和老套之嫌。后来想到曾有一部漫画叫《一本漫画闯天涯》,遂将“漫画”改成“对称”,“天涯”改成“物理”,似乎也能表达我对“对称”的感觉。

对称就是在某种变换下保持不变的性质,比如狭义相对论要求所有物理定律在所有惯性系中保持不变,这相对于要求描述物理定律的方程在匀速运动的坐标变换下保持不变,结合光速不变的要求,我们就可以推导出洛伦兹变换,从而完成地描述了狭义相对论里边的对称。然而,并不是任何时候都可以想推导洛伦兹变换那样,能够把一个完整的变换推导出来的。幸好,李对称的不需要完整的对称描述,它只需要“无穷小变换”(意味着我们可以忽略掉高阶项),对应地产生一个“无穷小生成元”,用这个无穷小生成元,就足以完整构建出我们所需要的物理来。这种“无穷小”决定“广泛”、“局部”决定“全局”的奇妙至今仍让我觉得不可思议。(关于李对称、无穷小生成元的基本概念,不妨先阅读:《求解微分方程的李对称方法》

点击阅读全文...

9 Apr

趣题:与橡皮绳赛跑的蚂蚁

这是一道流传很广的趣题,也许不少读者已经听说过它,然而广为人知却不一定“广为人‘解’”,在此把题目给出来,写下我自己的答案,并且谈谈我对答案的看法。题目是这样子的:

与橡皮绳赛跑的蚂蚁
一只蚂蚁沿着一条长$l=100$米的橡皮绳以每秒1厘米的匀速由一端向另一端爬行。每过1秒钟,橡皮绳就拉长100米,比如10秒后,橡皮绳就伸长了1000米。假设橡皮绳可以任意拉长,并且拉伸是均匀的。蚂蚁也会不知疲倦的一直往前爬,在绳子均匀拉长时,蚂蚁的位置理所当然的相对匀速向前挪动,问,如此下去,蚂蚁能否最终爬到橡皮绳的另一端?

点击阅读全文...

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

点击阅读全文...