关于维度

Posted 2020-07-03

1.什么是维度。

其实这个话题是欧氏几何的一个延伸。我们称零维的点，一维的线，二维的面，三维的体，四维的时空。你要注意到，这里0，1，2，3，4都是整数。你有没有想过，到底什么是维度？有没有分数维？比如3.1415926维。讨论这个的数学分支被称为分形数学。事实上分形数学已经广泛应用于物理，化学，地质，金融，社会科学等的方方面面，甚至到艺术及时尚。那么什么叫分形，什么是维度？先从一组图看起。

 

                                 

 图0： 大自然中分形结构的例子  

这组图特点在于，每个图中，某一部分挑出来都跟整体有类似性：比如，每根树枝结构都类似于整棵树。早在1860年，Rubkin就写到：”…大自然的鬼斧神工，可以把大尺度的山脉缩小成小尺度，在一块石头上你可以找到自然界的各种构造的变化。石块上的苔藓像森林，晶粒像悬崖，岩石表面像山脉…”。具有自相似性或标度不变性的几何对象我们成为它是分形的。严格的分形数学基础是测度论及公度拓扑学，艰深难懂。分形理论的奠基人Mandelbrot给出了更直观的定义：部分以某种形式与整体相似的集合，称作分形。

欧式几何的法则之一是整数维，这里把维度的概念做一拓展。考虑：把一个正方形的每个边变为原来的三倍，得到一个9倍面积的正方形：3^2=9; 类似的，把一个正方体每个边变为原来的三倍，得到一个27倍体积的新正方体：3^3=27。推而广之，把一个d维的物体，将每个独立方向变为原来的n倍，结果得到N个原来的对象。这三个数之间的关系是n^d=N。对此关系式取对数，还可以改写为：

d=ln(N)/ln(n)。

这个认识已经是一个质的飞跃。因为当我们不知道维度，用右边进行计算的时候，很可能得到的不是整数（但是在欧式几何的公理空间里，必为整数），也因此扩展了欧式几何。下面以一个简单而著名的Koch曲线的例子，来描述如何构造分形，以及如何计算其维数。

如图1，从一段单位长度的线段开始，第一步将线段三等分，将中间部分用两个1/3长度线段取代（图1.b），第二步，将现有每条线段三等分，中间有两条1/3^2的线段取代。如此反复直至无穷。Koch曲线是分形的，因为它有自相似性。同时要注意，Koch曲线无限长，处处连续，处处不可微商。那么它的维度是多少呢。将图2中部图形放大三倍后，得到图二下部的曲线，它是由原来的四个图形组成的，因此，Koch曲线的分形维数为d=ln4/ln3=1.26。

 

除此以外，常见的分形还有康托集，Sierpinski垫片。此类分形称作规则分形或者决定论的分形，还有些无规则的分形，大多物理现象中的分形为此类，其特点是不具备严格的自相似性，而是统计意义上的自相似。后者更加复杂，以后会单独分析布朗运动的分形特性。

刚才介绍的分形维数的计算只适合规则分形，还有其他不同的方式计算复杂分形的维数，其中之一为Hausdoff维度：

很吓人对不对？其实通俗来讲，计算hausdoff维度的方法是（以嵌入二维平面的分形为例）：

设想有一个由2维空间内具有有限大小的点组成的集合，N是用来覆盖这个集合内所有点所需的半径为R的圆形的最少个数，则这个最小数N是R的一个函数，记作N(R)。显然R越小则N越大，假设N(R)和R之间存在一个反比的关系，我们把这个关系记作

当R趋向于0时，我们得到Hausdoff维度为

 

Hausdoff维数的计算并不容易，尤其是想通过数值计算。通常只会给出一个范围。以下程序中有具体实现。 

 2.著名的分形集

前边列举了有具体几何形状的分形几何体，这里列举2个著名分形集合（当然，严格说几何体也是集合）。

2.1 Mandelbrot/Julia集合

定义复二次多项式f(Z)=Z^2+c，考虑无穷次迭代数列f(z),f(f(z))…直至无穷f(f(f…f(z))),假如此序列的模收敛，z为集合中一元素。当c为特定常数的情况下为Julia集合；假如c非定常，而是当前点的坐标的情况下，我们称为Mandelbrot集合。

2.2 Newton 分形

 

3. 分形的程序实现（python）

以下程序包含刚刚提到的分形集，程序截图如下：

程序除了plot Julia(稍作修改就可plot Mandelbrot集)/Newton集外，还提供了Housdoff维度的计算函数。

程序上色采用了最简单的方法：计算逃逸时间。比如4次迭代后发散，此点数值为4，收敛点（在最大迭代次数到达前未发散）此点数值为itermax，这里的contour仅仅根据此数值来做图。

"""
   _author_ = 杜鹏飞, Iowa state university
   _Python_version_= 2.7 
"""
import numpy as np
import matplotlib.pyplot as plt
import itertools as iter

class Fractal(object):
    def __init__(self, name, color, iteration, pixelNum, domain):
        self.name = name               #name of the fractal plot
        self.colorScheme = color       #which color scheme to use, current only 1
        self.maxIter = iteration       #maximum iteration
        self.pixelNum = pixelNum       #number of pixels along each dimension
        self.domain = domain           #(xmin,xmax,ymin,ymax) square domain
        #create the mesh grid
        self.x = np.linspace(self.domain[0],self.domain[1],pixelNum)
        self.y = np.linspace(self.domain[2],self.domain[3],pixelNum)
    self.xx, self.yy = np.meshgrid(self.x, self.y, sparse=True)


    def plot(self,matrix):
        #ravel() returns a contiguous flattened array;
        self.pic = np.reshape(matrix,[self.pixelNum,self.pixelNum])
        plt.imshow(self.pic)
        plt.title(self.name, fontsize=18, y=1.03)
        plt.tight_layout()
        plt.show()

    @staticmethod
    def abs(x, y):
        return np.sqrt(x*x + y*y)

    @staticmethod
    def square(x,y):
        return (x*x-y*y,2*x*y)

#subclass for Newton fractal drawing
class Newton(Fractal):
    def __init__(self, name, color, iteration, pixelNum, domain):
        Fractal.__init__(self, name, color, iteration, pixelNum, domain)
        #------construct the one variable function
        self.f = np.poly1d([1,0,0,-1]) # x^3 - 1
        #------define the derivative of f
        self.fp = np.polyder(self.f)

    def newton(self, i, guess):
        a = np.empty(guess.shape,dtype=int)
        a[:] = i
        j = np.abs(self.f(guess))>.00001
        if np.any(j):
            a[j] = self.newton(i+1, guess[j] - np.divide(self.f(guess[j]),self.fp(guess[j])))
        return a

    def plot(self):
        self.matrix = self.newton(0,np.ravel(self.xx+self.yy*1j)) 
        Fractal.plot(self,self.matrix)
 

#subclass for Julia fractal drawing
class Julia(Fractal):
    def __init__(self, name, color, iteration, pixelNum, domain, c0, c1):
        Fractal.__init__(self, name, color, iteration, pixelNum, domain)
        self.c_re = c0
    self.c_im = c1
        self.iter = 0
    self.matrix = np.zeros((pixelNum,pixelNum), dtype=np.int)

    def selfIterate(self,a,b):
        if self.abs(a,b) > 2:
            return 0
        else:
            a,b = self.square(a,b) 
            a += self.c_re
        b += self.c_im
        if self.iter < self.maxIter:
        self.iter += 1
                self.selfIterate(a,b)
          if self.iter >= self.maxIter:
             return self.iter
            return self.iter

    def julia(self):
        #wrong" self.matrix = [self.selfIterate(self.x[i],self.x[j]) for i in range(self.pixelNum) for j in range(self.pixelNum)]
        for i,j in iter.product(range(self.pixelNum), range(self.pixelNum)):
        self.matrix[i,j] = self.selfIterate(self.x[i],self.x[j])
        self.iter = 0
    return 0

    def plot(self):
        self.julia()
        Fractal.plot(self,self.matrix)

newton = Newton("Newton set",1,1000,1000,(-10,10,-10,10))
newton.plot()
#j = Julia("Julia set",1,50,1000,(-2,2,-2,2),-0.4,0.6)
#j.plot()