2.1.7 内存结构

下面让我们看看数组对象是如何在内存中存储的。如图2-3所示，数组的描述信息保存在一个数据结构中，这个结构引用两个对象：用于保存数据的存储区域和用于描述元素类型的dtype对象。

图2-3 ndarray数组对象在内存中的存储方式

数据存储区域保存着数组中所有元素的二进制数据，dtype对象则知道如何将元素的二进制数据转换为可用的值。数组的维数和形状等信息都保存在ndarray数组对象的数据结构中。图2-3中显示的是下面的数组a的内存结构：

    a = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32)

数组对象使用strides属性保存每个轴上相邻两个元素的地址差，即当某个轴的下标增加1时，数据存储区中的指针所增加的字节数。例如图2-3中的strides为(12,4)，即第0轴的下标增加1时，数据的地址增加12个字节。也就是a[1,0]的地址比a[0,0]的地址大12，正好是3个单精度浮点数的总字节数。第1轴的下标增加1时，数据的地址增加4个字节，正好是一个单精度浮点数的字节数。

如果strides属性中的数值正好和对应轴所占据的字节数相同，那么数据在内存中是连续存储的。通过切片下标得到的新数组是原始数组的视图，即它和原始数组共享数据存储区域，但是新数组的strides属性会发生变化：

    b = a[::2, ::2]
         b      b.strides
    ------------  ---------
    [[ 0.,  2.],  (24, 8)  
     [ 6.,  8.]]           

由于数组b和数组a共享数据存储区，而数组b中的第0轴和第1轴都是从a中隔一个元素取一个，因此数组b的strides变成了(24, 8)，正好都是数组a的两倍。对照前面的图2-3很容易看出数据0和2的地址相差8个字节，而数据0和6的地址相差24个字节。

元素在数据存储区中的排列格式有两种：C语言格式和Fortran语言格式。在C语言中，多维数组的第0轴是最上位的，即第0轴的下标增加1时，元素的地址增加的字节数最多；而Fortran语言中的多维数组的第0轴是最下位的，即第0轴的下标增加1时，地址只增加一个元素的字节数。在NumPy中默认以C语言格式存储数据，如果希望改为Fortran格式，只需要在创建数组时，设置order参数为"F"：

    c = np.array([[0,1,2],[3,4,5],[6,7,8]], dtype=np.float32, order="F")
    c.strides
    (4, 12)

了解了数组的内存结构，就可以解释使用下标取得数据时的复制和引用问题：

●当下标使用整数和切片时，所取得的数据在数据存储区域中是等间隔分布的。因为只需要修改图2-3所示的数据结构中的dim count、dimensions、stride等属性以及指向数据存储区域的指针data，就能实现整数和切片下标，所以新数组和原始数组能够共享数据存储区域。

●当使用整数序列、整数数组和布尔数组时，不能保证所取得的数据在数据存储区域中是等间隔的，因此无法和原始数组共享数据，只能对数据进行复制。

数组的flags属性描述了数据存储区域的一些属性，直接查看flags属性将输出各个标志的值，也可以单独获得其中的某个标志值：

    print a.flags
    print "c_contiguous:", a.flags.c_contiguous
      C_CONTIGUOUS : True
      F_CONTIGUOUS : False
      OWNDATA : True
      WRITEABLE : True
      ALIGNED : True
      UPDATEIFCOPY : False
    c_contiguous: True

下面是几个比较重要的标志：

●C_CONTIGUOUS：数据存储区域是否是C语言格式的连续区域。

●F_CONTIGUOUS：数据存储区域是否是Fortran语言格式的连续区域。

●OWNDATA：数组是否拥有此数据存储区域，当一个数组是其他数组的视图时，它不拥有数据存储区域。

由于数组a是通过array()直接创建的，因此它的数据存储区域是C语言格式的连续区域，并且它拥有数据存储区域。下面我们看看数组a的转置标志，数组的转置可以通过其T属性获得，转置数组将其数据存储区域看作Fortran语言格式的连续区域，并且它不拥有数据存储区域。

    a.T.flags
    C_CONTIGUOUS : False
      F_CONTIGUOUS : True
      OWNDATA : False
      WRITEABLE : True
      ALIGNED : True
      UPDATEIFCOPY : False

下面查看数组b的标志，它不拥有数据存储区域，其数据也不是连续存储的。通过视图数组的base属性可以获得保存数据的原始数组：

    b.flags
    C_CONTIGUOUS : False
      F_CONTIGUOUS : False
      OWNDATA : False
      WRITEABLE : True
      ALIGNED : True
      UPDATEIFCOPY : False
    id(b.base)    id(a)  
    ----------  ---------
    119627272   119627272

我们还可以通过view()方法从同一块数据区创建不同的dtype的数组对象，也就是使用不同的数值类型查看同一段内存中的二进制数据：

由于数组a的元素类型是单精度浮点数，占用4个字节，通过a.view(np.uint32)，我们创建了一个新的数组，它和数组a使用同一段数据内存，但是它将每4个字节的数据当作无符号32位整数处理。而a.view(np.uint8)将每个字节都当作一个单字节的无符号整数，因此得到一个形状为(3, 8)的数组。通过view()方法获得的新数组与原数组共享内存，当a[0, 0]被修改时，b[0, 0]和c[0, :4]都会改变：

    a[0, 0] = 3.14
    b[0, 0]          c[0, :4]      
    ----------  --------------------
    1078523331  [195, 245,  72,  64]

下面我们看一个使用view()方法的有趣的例子。在《雷神之锤III：竞技场》的C语言源代码中有这样一个神奇的计算平方根倒数的函数Q_rsqrt()。代码中使用牛顿迭代法计算平方根倒数，这并没有任何神奇之处，但是其中包含了一个神奇的数字0x5f3759df，并将单精度浮点数当作32位的整数进行了一次令人毫无头绪的运算：

    float Q_rsqrt( float number )
    {
        long i;
        float x2, y;
        const float threehalfs = 1.5F;
    
        x2 = number * 0.5F;
        y  = number;
        i  = * ( long * ) &y;                          // 对浮点数的邪恶位级hack
        i  = 0x5f3759df - ( i >> 1 );            // 这到底是怎么回事？
        y  = * ( float * ) &i;
        y  = y * ( threehalfs - ( x2 * y * y ));       // 第一次牛顿迭代
        return y;
    }

下面我们用NumPy实现同样的计算：

    number = np.linspace(0.1, 10, 100)
    y = number.astype(np.float32) ❶
    x2 = y * 0.5
    i = y.view(np.int32) ❷
    i[:] = 0x5f3759df - (i >> 1) ❸
    y = y * (1.5 - x2 * y * y) ❹
    np.max(np.abs(1 / np.sqrt(number) - y)) ❺
    0.0050456140410597428

❶由于linspace()创建的数组的类型为双精度浮点数，因此这里首先通过astype()方法将其转换成单精度浮点数数组y。❷通过view()方法创建一个与y共享内存的32位整数数组i。❸对整数数组i进行那段完全摸不着头脑的运算，并且将结果重新写入数组i中。由于i和y共享内存，此时y中的值也发生了变化。注意这里的赋值不能使用i=0x5f3759df -(i >> 1)，如果这样写，那么数组i就是一个全新的数组了。❹进行一次牛顿迭代运算，这里由于使用y=...的写法，因此y将变成一个全新的数组，和原来的i不再共享内存。在这段代码中有很多数组运算，关于这方面的内容将在下一节进行详细说明。❺最后输出真实值和近似值之间的最大误差。图2-4显示了绝对误差与自变量的关系，当number很小时绝对误差较大，但此时的函数值也较大，因此相对误差的变化并不大。

图2-4 《雷神之锤》中计算平方根倒数算法的绝对误差

除了使用切片从同一块数据区创建不同的shape和strides的数组对象之外，还可以直接设置这些属性，从而得到用切片实现不了的效果，例如：

    from numpy.lib.stride_tricks import as_strided
    a = np.arange(6)
    b = as_strided(a, shape=(4, 3), strides=(4, 4))
            a                b     
    ------------------  -----------
    [0, 1, 2, 3, 4, 5]  [[0, 1, 2],
    [1, 2, 3],
    [2, 3, 4],
    [3, 4, 5]]

这个例子中，我们从NumPy的辅助模块中载入了as_strided()函数，并使用它从一个长度为6的一维数组a创建了一个shape为(4, 3)的二维数组b。由于通过strides参数直接指定了数组b的strides属性，因此不仅数组b和数组a共享数据区，而且b中的前后两行有两个元素是重合的。例如下面修改a[2]的值，b中的前三行中对应的元素也发生改变：

    a[2] = 20
    b
    array([[ 0,  1, 20],
           [ 1, 20,  3],
           [20,  3,  4],
           [ 3,  4,  5]])

在对数据进行处理时，可能经常需要对数据进行分块处理，而且为了保持平滑，每块数据之间需要有一定的重叠部分。这时可以使用上面介绍的方法对数据进行带重叠的分块。需要注意的是，使用as_strided()时NumPy不会进行内存越界检查，因此shape和strides设置不当可能会发生意想不到的错误。