每一个好的抽象都能解决一个问题，而这篇文章将涵盖我迄今为止所知道的一个杰出的数学抽象——张量。从简单的 2 层多层感知机到 GPT-5，神经网络归根结底都是浮点数在操作图中流动。这篇文章将从头开始用 C 语言构建一个完整的加速张量库。它深受 Bellard 的 libnc 启发，但不幸的是该库尚未开源。张量只不过是一个扁平的数字数组，加上一些元数据，告诉你如何将这些数字解释为多维对象。我们都学过 2D 数组可以更好地表示为一个 1D 数组加上行/列的数量——这本质上就是张量。但超越两个维度，我们可能需要一些其他的元数据，例如一个广义形状：float data [ 32 * 3 * 28 * 28 ]; // 32 张图像，3 个通道，28x28 像素 int shape [ 4 ] = { 32 , 3 , 28 , 28 }; // 张量的形状 int ndim = 4 ; // 维度的数量 拥有形状后，我们可以得知在 4D 张量中位置 [n,c,h,w] 的元素在数据中的偏移量为 data[n*(3*28*28)+c*(28*28)+h*28+w] ，如果我们将张量保存在按行优先的 C 顺序格式中（通常是目前大多数张量库的默认格式）。现在，每次计算元素索引都会低效，因此一旦知道形状，我们可以预先计算步幅。步幅指示如果我们想在给定维度上推进一个元素时要跳过多少个元素。我们还可以将所有与形状相关的字段组合在一起： struct ut_shape { int ndim ; // 维度的数量 int nelem ; // 元素的数量 int shape [ 4 ]; // 张量的形状 int strides [ 4 ]; // 每个维度的步幅 }; struct ut_tensor { struct ut_shape shape ; // 张量的形状 float * data ; // 指向数据的指针 ... // 之后会有更多字段 }; 我们可以添加一些辅助函数来创建形状并从多维索引获取扁平索引： ut_shape ut_shape_new ( int ndims , int * dims ) { ut_shape s = {. ndims = ndims , . nelems = 1 }; for ( int i = 0 ; i < ndim ; i ++ ) { s . shape [ i ] = dims [ i ]; s . nelem *= dims [ i ]; } return s ; } int ut_index ( ut_shape s , const int * idx ) { int flat = 0 , stride = 1 ; for ( int i = s . ndim - 1 ; i >= 0 ; i -- ) { flat += idx [ i ] * stride ; stride *= s . shape [ i ]; } return flat ; } ut_shape s = ut_shape_new ( 3 , ( int []){ 2 , 3 , 4 }); assert ( s . nelem == 24 ); assert ( s . ndim == 3 ); // 元素 [1][2][3] 应在偏移量 1*12 + 2*4 + 3 = 23 assert ( ut_index ( s , ( int []){ 1 , 2 , 3 }) == 23 ); 张量通常是动态分配的，因此我们应该提供方式来创建和销毁它们，仅仅是在 malloc / free 之上的一个封装。有时我们想要创建一个张量，该张量与另一个张量共享相同的数据，例如获得张量某一部分的“视图”，或在不复制数据的情况下转置张量，或者修改其形状（扁平化）。在这种情况下，我们需要跟踪数据的所有权。保留张量也是有用的，以便它们可以超出其原始范围，例如在迭代训练中以避免相同分配的重复发生。因此，我们在张量结构中添加了一个引用计数字段： struct ut_tensor { struct ut_shape shape ; // 张量的形状 float * data ; // 指向数据的指针 int refcount ; // 共享所有权的引用计数，达到 0 时释放 struct ut_tensor * owner ; // 如果该张量是视图，则指向拥有张量的指针 }; 我们可以进一步优化内存管理，添加区域分配器或内存池以避免频繁的 malloc / free 调用，但我们已经拥有的东西是一个良好的开端。然而，如果没有对它们的操作，我们的张量几乎没有用处。逐元素 操作张量的最基本操作是逐元素操作——对所有元素进行单次循环，对每个元素（单目）或一对元素（双目）应用函数。我们可以像这样实现它们： static void ew_neg ( float * out , const float * a , int n ) { for ( int i = 0 ; i < n ; i ++ ) out [ i ] = - a [ i ]; } // ...更多单目操作... static void ew_relu ( float * out , const float * a , int n ) { for ( int i = 0 ; i < n ; i ++ ) out [ i ] = fmaxf ( 0.f , a [ i ]); } static void ew_add ( float * out , const float * a , const float * b , int n ) { for ( int i = 0 ; i < n ; i ++ ) out [ i ] = a [ i ] + b [ i ]; } // ...更多双目操作... static ut_tensor * ew_unary ( ut_tensor * a , void ( * fn )( float * , const float * , int )) { ut_tensor * out = ut_alloc ( a -> shape . ndim , a -> shape . shape ); fn ( out -> data , a -> data , a -> shape . nelem ); return out ; } static ut_tensor * ew_binary ( ut_tensor * a , ut_tensor * b , void ( * fn )( float * , const float * , const float * , int )) { ut_tensor * out = ut_alloc ( a -> shape . ndim , a -> shape . shape ); fn ( out -> data , a -> data , b -> data , a -> shape . nelem ); return out ; } ut_tensor * ut_neg ( ut_tensor * a ) { return ew_unary ( a , ew_neg ); } ut_tensor * ut_exp ( ut_tensor