kgr*_*yte 5 arrays algorithm math numpy multidimensional-array
是否存在将线性索引转换为支持负步幅的下标列表的算法？
背景
MATLAB、Julia 等环境和 NumPy 等库提供对跨步数组（也称为 ndarray）的支持。跨步数组由线性内存（例如，单个底层缓冲区）支持，这与嵌套数组形成鲜明对比，其中每个嵌套数组对应一个维度。例如，考虑以下 2x2 矩阵
[ 1 2
3 4 ]
作为数组的数组来实现
[ 1 2
3 4 ]
其中（使用从零开始的索引）
A = [ [ 1, 2 ], [ 3, 4 ] ]
我们可以将相同的 2x2 矩阵表示为跨步数组，如下所示（假设行优先）
a01 = A[0][1] = 2
在哪里
A = [ 1, 2,
3, 4 ]
一般来说，对于跨步 NxM 矩阵，(i,j)可以通过以下方式访问元素
a01 = A[ 2*0 + 1*1 ] = 2
其中buffer是底层数据缓冲区，si和分别对应于沿和维度sj的步幅。假设一个行主跨步数组，对于上面的 2x2 矩阵，并且（省略元素字节长度）。ijsi = 2sj = 1
一般来说，步长可以根据数组形状计算如下：
function get( i, j ) {
return buffer[ si*i + sj*j ];
}
为了方便使用跨步数组，环境/库通常提供方便的函数，允许在线性索引和下标之间轻松转换。例如，在MATLAB中，从下标转换为线性索引
function shape2strides( shape, order ) {
var out = new Array( shape.length );
var s = 1;
var i;
if ( order === 'column-major' ) {
for ( i = 0; i < shape.length; i++ ) {
out[ i ] = shape[ i ];
s *= shape[ i ];
}
return out;
} else { // row-major
for ( i = shape.length-1; i >= 0; i-- ) {
out[ i ] = shape[ i ];
s *= shape[ i ];
}
}
}
同样，在MATLAB中将线性索引转换为下标
idx = sub2ind( size( A ), i, j )
Julia 还有sub2ind和ind2sub。在 NumPy 中，您可以使用ravel_multi_index和unravel_index。
除了数据局部性之外，跨步数组也很方便，因为它们允许通过操作跨步是负数还是正数来创建数组“视图”。当步幅为负时，我们不是从左到右迭代，而是沿着该维度从右到左迭代。为了支持这种迭代行为，我们需要确定第一个索引元素在底层数据缓冲区中的位置。按照惯例，我们将该索引称为“偏移量”，其计算公式如下
s = ind2sub( size( A ), idx )
一旦我们有了偏移量，我们需要修改我们的get( i, j )函数，如下所示
function strides2offset( shape, strides ) {
var offset = 0;
var i;
for ( i = 0; i < shape.length; i++ ) {
if ( strides[ i ] < 0 ) {
offset -= strides[i] * ( shape[i]-1 ); // increments the offset
}
}
return offset;
}
对于步幅为2x2 的矩阵A2,1，偏移量为0，从而返回上面的原始get函数。当步幅为 时2,-1，偏移量为1；对于-2,1，偏移量为2；对于-2,-1，偏移量为3。因此，我们可以生成以下矩阵视图（假设行优先）
Dims: 2x2
Strides: 2,1
Offset: 0
A = [ 1, 2,
3, 4 ]
Strides: 2,-1
Offset: 1
A = [ 2, 1,
4, 3 ]
Strides: -2,1
Offset: 2
A = [ 3, 4,
1, 2 ]
Strides: -2,-1
Offset: 3
A = [ 4, 3,
2, 1 ]
上述视图突出了跨步数组的优点之一：O(1) 操作。例如，要从左到右翻转矩阵，我们只需要翻转第二维的步幅符号（假设行优先）。为了从上到下翻转，我们翻转第一个维度的步幅符号（假设行优先）。为了从左到右、从上到下翻转，我们翻转两个步幅的符号。所有上述操作均不涉及底层数据缓冲区；我们只需更改跨步数组元数据即可。
次二工业
从下标转换为线性索引很简单，即使考虑负步幅（即步幅数组视图）也是如此。例如，对于任意维度的跨步数组，
function get( i, j ) {
return buffer[ offset + si*i + sj*j ];
}
在这里，我们允许从视图的角度或从底层数据缓冲区的角度返回线性索引。当“偏移量”为 时0，我们假设我们始终将线性索引返回到视图中（这可能与底层数据缓冲区中的线性索引不对应）。换句话说，对于 2x2 矩阵视图(0,0) => 0, (0,1) => 1, (1,0) => 2, (1,1) => 3，总是。从以下角度来看，这是有道理的：在使用视图时，这种映射符合直觉。当我想要时A(0,0)，我希望该元素位于“第一个”线性索引处，即使这不是该元素实际存储在底层数据缓冲区中的位置。
您可以向自己证明，sub2ind当将元素查找扩展到负步长时，对于任何偏移量都会返回相同的索引，如上面详述的那样。
工业2sub
这里提出的问题是如何实现 的相反
sub2ind，并支持负步幅。
对于正步长（以及因此的偏移量0），我们可以使用模算术来恢复下标。例如，考虑求解 NxMxL 跨步数组的线性索引的方程。
idx = offset + si*i + sj*j + sk*k
其中，假设行优先，si = nj*nk, sj = nk, sk = 1和分别ni, nj, nk是维度大小N, M, L。替换值，
idx = 0 + (nj*nk)*i + nk*j + k
可以重新排列
idx = nk*(nj*i + j) + k
如果我们使用 两边取模nk，
idx % nk = k
知道了k，我们重新排列一下初始方程
(idx - k) = nk*(nj*i + j)
(idx - k)/nk = nj*i + j
如果我们使用 两边取模nj，
((idx - k)/nk) % nj = j
知道了j，让我们重新排列初始方程来求解i
(((idx - k)/nk) - j)/nj = i
上述算法可推广到任意数量的维度，并且易于实现（另请参阅Julia和NumPy）。
Dims: 2x2
Strides: 2,1
Offset: 0
A = [ 1, 2,
3, 4 ]
Strides: 2,-1
Offset: 1
A = [ 2, 1,
4, 3 ]
Strides: -2,1
Offset: 2
A = [ 3, 4,
1, 2 ]
Strides: -2,-1
Offset: 3
A = [ 4, 3,
2, 1 ]
然而，上述使用模运算的算法不支持负步幅。如果我们使用上面相同的过程来求解下标i,j,k，我们将从方程开始
idx = offset + nk*(nj*i + j) + k
可以简化为
idx-offset = nk*(nj*i + j) + k
当然，这里的问题是idx-offset可以为负数并且有效地改变了可能值的范围i,j,k（i应该在半开区间 [0,N）上；j在区间 [0,M) 上；以及k区间 [0,L))。
这就提出了一个问题：是否存在将线性索引转换为支持负步长的下标的算法。或者换句话说，是否有一种算法，给定底层数据缓冲区的线性索引，可以返回相应的视图下标？
其他语言/库（例如Julia和NumPy）中的实现似乎只提供对这种offset = 0情况的支持。我正在寻找更通用的东西，它也可以应用于跨步数组视图。
任何指向现有实现/算法的指针将不胜感激。
（编辑 - 我可能正在处理更容易的 nd-index 到 flat 索引情况，而你则专注于相反的情况。现在探索该任务为时已晚 - 我将在早上重新审视这个问题。）
如果偏移量正确，我认为将 nd 索引转换为平面索引的相同公式适用于负步幅和正步幅：
将 (3,4) 数组的索引与其双翻转进行比较：
In [32]: x = np.arange(12).reshape(3,4)
In [33]: y = x[::-1, ::-1]
In [34]: x.strides
Out[34]: (16, 4)
In [35]: y.strides
Out[35]: (-16, -4)
数据缓冲区“开始”可以通过以下方式找到__array_interface__：
In [36]: x.__array_interface__['data']
Out[36]: (166934688, False)
In [37]: y.__array_interface__['data']
Out[37]: (166934732, False)
In [38]: 732-688
Out[38]: 44
缓冲区中有 48 个字节，但偏移量y为 44，即 x[2,3] 元素的开始位置（本例中为 11）。
现在测试元素的平面索引x：
In [39]: x[1,2]
Out[39]: 6 # value
In [40]: x[1:2, 2:3].__array_interface__['data'] # a view
Out[40]: (166934712, False)
In [41]: 688+(1*16)+(2*4) # offset + si*i + sj*j
Out[41]: 712
现在对 执行同样的操作y：
In [42]: y[1:2, 2:3].__array_interface__['data']
Out[42]: (166934708, False)
In [43]: 732+(1*-16)+(2*-4)
Out[43]: 708