惊奇数学事实

Tupper自我指涉公式：图象里竟然包含式子本身

/ 32 条评论

    你认为，一个函数图象里是否有可能包含这个函数本身的“图象”？难以置信的是，还真有人构造了这样一个东西。2001年，Jeff Tupper发表的一篇论文里提到了这样一个有趣的不等式：
  
    在0 <= x <= 105，n <= y <= n + 16的范围内，这个不等式对应的图象是这个样子：
  

其中，n = 96093937991895888497167296212785275471500433966012930665150551927170280239526642
46896428421743507181212671537827706233559932372808741443078913259639413377234878
57735749823926629715517173716995165232890538221612403238855866184013235585136048
82869333790249145422928866708109618449609170518345406782773155170540538162738096
76025656250169814820834187831638491155902256100036523513703438744618483787372381
98224849863465033159410054974700593138339226497249461751545728366702369745461014
655997933798537483143786841806593422227898388722980000748404719

    你会觉得这个很神奇吗？你也许会想，天哪，这个是怎么构造出来的啊！但仔细思考之后，你会发现这个一点都不神奇。事实上明白了道理之后你可以构造出无数个这样的式子来。现在给你一些时间让你思考一下，你能否看出其中的奥秘？

    就像魔术揭秘一样，说穿了真相后上面的这些东西就一点意思都没有了。在这个式子里，涉及到x和y的变量时都加上了取整符号，因此整个图象都是一格一格的。这样，不等式右边的式子就简化为y div 17 * 2^(-17x – y mod 17) mod 2，其中x和y都为整数。接着观察，一个数乘以2的负k次方相当于对应的二进制数右移k位，那么x * 2^(-k) mod 2实质上就是二进制数x右起第k位上的数字。对于某个自然数t，当17t <= y < 17(t+1)时，指数-17x – y mod 17恰好对应所有的负整数，于是位于y=17t和y=17t+16之间的图象的每个像素和t的二进制中的每一位数字一一对应。随着t值的增加，图形的像素会一点一点地变化。当纵坐标足够大时，必然会出现一段高度为17的图象，图象的样子和不等式本身的样子相同。当然，你也可以在里面“找到”任何你想要的图象，只需要把图象还原为二进制数并转换为十进制即可。你甚至可以告诉你的MM，说你发现了一个函数，函数在某个位置的图象正好是某某某我爱你的字样。

Matrix67原创
转贴请注明出处
最近发现了一些很不厚道的人，希望大家注意哦！

神奇的分形艺术（四）：Julia集和Mandelbrot集

/ 44 条评论

    考虑函数f(z)=z^2-0.75。固定z0的值后，我们可以通过不断地迭代算出一系列的z值：z1=f(z0), z2=f(z1), z3=f(z2), …。比如，当z0 = 1时，我们可以依次迭代出：

z1 = f(1.0) = 1.0^2 – 0.75 = 0.25
z2 = f(0.25) = 0.25^2 – 0.75 = -0.6875
z3 = f(-0.6875) = (-0.6875)^2 – 0.75 = -0.2773
z4 = f(-0.2773) = (-0.2773)^2 – 0.75 = -0.6731
z5 = f(-0.6731) = (-0.6731)^2 – 0.75 = -0.2970
…

    可以看出，z值始终在某一范围内，并将最终收敛到某一个值上。
    但当z0=2时，情况就不一样了。几次迭代后我们将立即发现z值最终会趋于无穷大：

z1 = f(2.0) = (2.0)^2 – 0.75 = 3.25
z2 = f(3.25) = (3.25)^2 – 0.75 = 9.8125
z3 = f(9.8125) = (9.8125)^2 – 0.75 = 95.535
z4 = f(95.535) = (95.535)^2 – 0.75 = 9126.2
z5 = f(9126.2) = (9126.2)^2 – 0.75 = 83287819.2
…

    经过计算，我们可以得到如下结论：当z0属于[-1.5, 1.5]时，z值始终不会超出某个范围；而当z0小于-1.5或大于1.5后，z值最终将趋于无穷。
    现在，我们把这个函数扩展到整个复数范围。对于复数z0=x+iy，取不同的x值和y值，函数迭代的结果不一样：对于有些z0，函数值约束在某一范围内；而对于另一些z0，函数值则发散到无穷。由于复数对应平面上的点，因此我们可以用一个平面图形来表示，对于哪些z0函数值最终趋于无穷，对于哪些z0函数值最终不会趋于无穷。我们用深灰色表示不会使函数值趋于无穷的z0；对于其它的z0，我们用不同的颜色来区别不同的发散速度。由于当某个时候|z|>2时，函数值一定发散，因此这里定义发散速度为：使|z|大于2的迭代次数越少，则发散速度越快。这个图形可以编程画出。和上次一样，我用Pascal语言，因为我不会C的图形操作。某个MM要过生日了，我把这个自己编程画的图片送给她^_^

{$ASSERTIONS+}

uses graph;

type
   complex=record
      re:real;
      im:real;
   end;

operator * (a:complex; b:complex) c:complex;
begin
   c.re := a.re*b.re - a.im*b.im;
   c.im := a.im*b.re + a.re*b.im;
end;

operator + (a:complex; b:complex) c:complex;
begin
   c.re := a.re + b.re;
   c.im := a.im + b.im;
end;

var
   z,c:complex;
   gd,gm,i,j,k:integer;
begin
   gd:=D8bit;
   gm:=m640x480;
   InitGraph(gd,gm,'');
   Assert(graphResult=grOk);

   c.re:=-0.75;
   c.im:=0;
   for i:=-300 to 300 do
   for j:=-200 to 200 do
   begin
      z.re:=i/200;
      z.im:=j/200;
      for k:=0 to 200 do
      begin
         if sqrt(z.re*z.re + z.im*z.im) >2 then break
         else z:=(z*z)+c;
      end;
      PutPixel(i+300,j+200,k)
   end;

   readln;
   CloseGraph;
end.

    代码在Windows XP SP2，FPC 2.0下通过编译，麻烦大家帮忙报告一下程序运行是否正常（上次有人告诉我说我写的绘图程序不能编译）。在我这里，程序运行的结果如下：

    这个美丽的分形图形表现的就是f(z)=z^2-0.75时的Julia集。考虑复数函数f(z)=z^2+c，不同的复数c对应着不同的Julia集。也就是说，每取一个不同的c你都能得到一个不同的Julia集分形图形，并且令人吃惊的是每一个分形图形都是那么美丽。下面的六幅图片是取不同的c值得到的分形图形。你可能不相信这样一个简单的构造法则可以生成这么美丽的图形，这没什么，你可以改变上面程序代码中c变量的值来亲自验证。

c = 0.45, -0.1428
  

c = 0.285, 0.01
  

c = 0.285, 0
  

c = -0.8, 0.156
  

c = -0.835, -0.2321
  

c = -0.70176, -0.3842
  

    类似地，我们固定z0=0，那么对于不同的复数c，函数的迭代结果也不同。由于复数c对应平面上的点，因此我们可以用一个平面图形来表示，对于某个复数c，函数f(z)=z^2+c从z0=0开始迭代是否会发散到无穷。我们同样用不同颜色来表示不同的发散速度，最后得出的就是Mandelbrot集分形图形：
    

    前面说过，分形图形是可以无限递归下去的，它的复杂度不随尺度减小而消失。Mandelbrot集的神奇之处就在于，你可以对这个分形图形不断放大，不同的尺度下你所看到的景象可能完全不同。放大到一定时候，你可以看到更小规模的Mandelbrot集，这证明Mandelbrot集是自相似的。下面的15幅图演示了Mandelbrot集的一个放大过程，你可以在这个过程中看到不同样式的分形图形。

网上可以找到很多小程序实现Mandelbrot集的放大过程。把上面给出的代码改一改，你也可以写出一个这样的程序来。

Update：2011 年 8 月 31 日，我对这个话题做了更进一步的讨论 http://www.matrix67.com/blog/archives/4570

令人敬畏的十维空间

/ 50 条评论


    我们把一个边长为2的正方形划分成4个小正方形，每个小正方形里作一个内切圆，然后在原来的大正方形中间作一个同时外切于这4个圆的小圆（红色标注）。我们把这个小圆叫做“中心圆”。你怎么来求这个中心圆的半径？
    仔细观察其中一个小正方形，思路就出来了：红色的中心圆变成了一个90度扇形，它的中心位于单位正方形的一角，并且外切于直径为1的圆。可以看到扇形半径加上圆的半径等于单位正方形对角线的一半，这样我们就得出，中心圆的半径等于(sqrt(2)-1)/2。
    对于一个立方体同样如此。我们把立方体切成8个小立方体，得到的8个球体中间夹住的那个中心球半径就应该为(sqrt(3)-1)/2。你会发现一个惊人的事实，在超立方体中，位于16个四维球体间的中心球半径为(sqrt(4)-1)/2 = 1/2，它竟然与那16个小球一样大。真正可怕的事情发生在九维立方体中，此时的九维中心球半径为(sqrt(9)-1)/2 = 1，竟然内切于最初的九维立方体！而到了十维空间后，中心球的直径将超过十维立方体的边长，这个中心球将突破立方体的边界！被围在里面的中心球居然比原来的N维立方体还大，这显然违反了大多数人的直觉；如果你能想象出这个画面来，你就牛B了。科幻小说中把对十维空间的感知能力作为文明发达程度的标准，除了一些相关的宇宙模型外，这可能也是其中一个原因吧。

又是黑色星期五！为什么总是黑色星期五呢？

/ 12 条评论

    你是否还记得，今年的4月13日是黑色星期五。短短三个月后，黑色星期五再次现身！为什么这一天老是出现呢？
    恐怕心理原因是最好的解释。人们对黑色星期五的出现记忆更深刻，给人一种黑色星期五常常出现的错觉。有趣的是，仔细算一算，你会发现13日是星期五的次数真的要多一些。
    很多人以为现在实行的历法是4年一循环，这是不对的。现在实行的历法以400年为一个循环。大家很容易忽略整百年的问题。一个很有意思的智力题就是问一个人是否可能连续5年不过生日。有个笑话说文科MM感叹她等了她男友4年，整整1460天，学理科的好友脱口而出“难道你男友是1900年的”。我们可以算一下在这400年中共有多少天：365*303+366*97=146097。这个数正好能被7整除。换句话说，现在与400年后的星期数不变，日历完全相同。在这400年里一共有4800个月，利用Zeller公式（见这里的最后一小节）可以编程统计出13日是星期几的次数最多，这对于OIer们再熟悉不过了，因为USACO有一道题就是干这种无聊的事情。下面就是程序运行后的结果：

  • Sunday     687
  • Monday     685
  • Tuesday    685
  • Wednesday  687
  • Thursday   684
  • Friday     688
  • Saturday   684


    可以看到，事实上13日是星期五的概率确实是最高的。
    另外，注意到了么，利用“400年一周期”这个结论我们可以对USACO的那个题进行扩展，出一个Friday the Thirteenth数据加强版。
    做人要厚道，转贴不注明出处者将受到黑色星期五的诅咒。

神奇的分形艺术（二）：一条连续的曲线可以填满整个平面

/ 16 条评论

    虽然有些东西似乎是显然的，但一个完整的定义仍然很有必要。比如，大多数人并不知道函数的连续性是怎么定义的，虽然大家一直在用。有人可能会说，函数是不是连续的一看就知道了嘛，需要定义么。事实上，如果没有严格的定义，你很难把下面两个问题说清楚。
    你知道吗，除了常函数之外还存在其它没有最小正周期的周期函数。考虑一个这样的函数：它的定义域为全体实数，当x为有理数时f(x)=1，当x为无理数时f(x)=0。显然，任何有理数都是这个函数的一个最小正周期，因为一个有理数加有理数还是有理数，而一个无理数加有理数仍然是无理数。因此，该函数的最小正周期可以任意小。如果非要画出它的图象，大致看上去就是两根直线。请问这个函数是连续函数吗？如果把这个函数改一下，当x为无理数时f(x)=0，当x为有理数时f(x)=x，那新的函数是连续函数吗？
    Cauchy定义专门用来解决这一类问题，它严格地定义了函数的连续性。Cauchy定义是说，函数f在x=c处连续当且仅当对于一个任意小的正数ε，你总能找到一个正数δ使得对于定义域上的所有满足c-δ< x <c+δ的x都有f(c)-ε<f(x)<f(c)+ε。直观地说，如果函数上有一点P，对于任意小的ε，P点左右一定范围内的点与P的纵坐标之差均小于ε，那么函数在P点处连续。这样就保证了P点两旁的点与P无限接近，也就是我们常说的“连续”。这又被称作为Epsilon-Delta定义，可以写成“ε-δ定义”。
    有了Cauchy定义，回过头来看前面的问题，我们可以推出：第一个函数在任何一点都不连续，因为当ε< 1时，δ范围内总存在至少一个点跳出了ε的范围；第二个函数只在x=0处是连续的，因为此时不管ε是多少，只需要δ比ε小一点就可以满足ε-δ定义了。
    在拓扑学中，也有类似于ε-δ的连续性定义。假如一个函数f(t)对应空间中的点，对于任意小的正数ε，总能找到一个δ使得定义域(t-δ,t+δ)对应的所有点与f(t)的距离都不超过ε，那么我们就说f(t)所对应的曲线在点f(t)处连续。

    回到我们的话题，如何构造一条曲线使得它可以填满整个平面。在这里我们仅仅说明存在一条填满单位正方形的曲线就够了，因为将此单位正方形平铺在平面上就可以得到填满整个平面的曲线。大多数人可能会想到下面这种构造方法：先画一条单位长的曲线，然后把它变成一个几字形，接着把每一条水平的小横线段变成一个几字形，然后不断迭代下去，最后得到的图形一定可以填满整个单位正方形。我们甚至可以递归地定义出一个描述此图形的函数：将定义域平均分成五份，第二和第四份对应两条竖直线段上的点，并继续对剩下的三个区间重复进行这种操作。这个函数虽然分布得有些“不均匀”，但它确实是一个合法的函数。最后的图形显然可以填充一个正方形，但它是不是一条曲线我们还不知道呢。稍作分析你会发现这条“曲线”根本不符合前面所说的ε-δ定义，考虑任何一个可以无限细分的地方（比如x=1/2处），只要ε<1/2，δ再小其范围内也有一条竖线捅破ε的界线。这就好像当n趋于无穷时sin(nx)根本不是一条确定的曲线一样，因为某个特定的函数值根本不能汇聚到一点。考虑到这一点，我们能想到的很多可以填满平面的“曲线”都不是真正意义上的连续曲线。为了避免这样的情况出现，这条曲线必须“先把自己周围填满再延伸出去”，而填满自己周围前又必须先填满“更小规模的周围”。这让我们联想到分形图形。

    德国数学家David Hilbert发现了这样一种可以填满整个单位正方形的分形曲线，他称它为Hilbert曲线。我们来看一看这条曲线是怎么构造出来的。首先，我们把一个正方形分割为4个小正方形，然后从左下角的那个小正方形开始，画一条线经过所有小正方形，最后到达右下角。现在，我们把这个正方形分成16个小正方形，目标同样是从左下角出发遍历所有的格子最后到达右下角。而在这之前我们已经得到了一个2×2方格的遍历方法，我们正好可以用它。把两个2×2的格子原封不动地放在上面两排，右旋90度放在左下，左旋90度放在右下，然后再补三条线段把它们连起来。现在我们得到了一种从左下到右下遍历4×4方格的方法，而这又可以用于更大规模的图形中。用刚才的方法把四个4×4的方格放到8×8的方格中，我们就得到了一条经过所有64个小方格的曲线。不断地这样做下去，无限多次地迭代后，每个方格都变得无穷小，最后的图形显然经过了方格上所有的点，它就是我们所说的Hilbert曲线。下图是一个迭代了n多次后的图形，大致上反映出Hilbert曲线的样子。
        

    根据上面这种方法，我们可以构造出函数f(t)使它能映射到单位正方形中的所有点。Hilbert曲线首先经过单位正方形左下1/4的所有点，然后顺势北上，东征到右上角，最后到达东南方的1/4正方形，其中的每一个阶段都是一个边长缩小了一半的“小Hilbert曲线”。函数f(t)也如此定义：[0,1/4]对应左下角的小正方形中所有的点，[1/4,1/2]就对应左上角，依此类推。每个区间继续划分为四份，依次对应面积为1/16的正方形，并无限制地这么细分下去。注意这里的定义域划分都是闭区间的形式，这并不会发生冲突，因为所有m/4^n处的点都是两个小Hilbert曲线的“交接处”。比如那个f(1/4)点就是左上左下两块1/4正方形共有的，即单位正方形正左边的那一点。这个函数是一条根正苗红的连续曲线，完全符合ε-δ定义，因为f(t-δ)和f(t+δ)显然都在f(t)的周围。
    Hilbert曲线是一条经典的分形曲线。它违背了很多常理。比如，把Hilbert曲线平铺在整个平面上，它就成了一条填满整个平面的曲线。两条Hilbert曲线对接可以形成一个封闭曲线，而这个封闭曲线竟然没有内部空间。回想我们上次介绍的Hausdorff维度，你会发现这条曲线是二维的，因为它包含自身4份，每一份的一维大小都是原来的一半，因此维度等于log(4)/log(2)。这再一次说明了它可以填满整个平面。

    Hilbert曲线的价值在于建立一维空间与二维空间一一对应的关系。Hilbert曲线可以看作是一个一维空间到二维空间的映射，也就是说我们证明了直线上的点和平面上的点一样多（集合的势相同）。Hilbert曲线也是一种遍历二维格点的方法，它同样可以用来证明自然数和有理数一样多。如果你已经知道此结论的Cantor证明，你会立刻明白Hilbert遍历法的证明，这里就不再多说了。当然，Hilbert曲线也可以扩展到三维空间，甚至更高维的空间，从而建立一维到任意多维的映射关系。下图就是一个三维Hilbert曲线（在迭代