matlab 将Fortran数组转换为C数组(列为主到行为主)

fbcarpbf  于 7个月前  发布在  Matlab
关注(0)|答案(2)|浏览(89)

给定一个任意形状的x,如何使o = x(idx)按照C而不是Fortran进行索引?也就是说,以行为主而不是列为主遍历展开的数组。
一个动机是,当遍历row-major时,某些条件应该导致一个排序数组:isequal(x(x > 0), sort(x))

说明

  • "* o = x(idx)indexes per C” 表示匹配Python NumPy的o = x[idx]。设xre = reorder(x),左= MATLAB,右= Python:
  • xre(:) == x.flatten()
  • xre(i1:i2, i3:i4, i5:i6) == x[j1:j2, j3:j4, j5:j6] ; j = i - 1

尝试直播
Definitions--“row-major”一般化为“outer-first”,“column-major”一般化为“inner-first”,指的是扁平数组的遍历,permute单独 * 无法 * 实现(理论上也不能,因为形状不能改变)。

演示

MATLAB:

s = @(x, i0, i1, i2)reshape(x(i0, i1, i2), 1, []);

xF = reshape(1:2*3*4, [2 3 4]);
{s(xF, 1:2, 1, 1), s(xF, 1, 1:2, 1), s(xF, 1, 1, 1:2)}
{[1 2]}    {[1 3]}    {[1 7]}

Python:

import numpy as np
xC = np.arange(1, 2*3*4+1).reshape(2, 3, 4)
print(xC[0:2, 0, 0], xC[0, 0:2, 0], xC[0, 0, 0:2])
[ 1 13] [1 5] [1 2]

MATLAB permute

xp = permute(xF, [3 2 1]);
{s(xp, 1:2, 1, 1), s(xp, 1, 1:2, 1), s(xp, 1, 1, 1:2)}
{[1 7]}    {[1 3]}    {[1 2]}

permute本身并不能实现我们所要求的。即使它以某种方式这样做了,它也保证是越界的,这取决于i0, i1, i2

nfzehxib

nfzehxib1#

根据你的要求,这里是列为主与行为主的顺序对我来说意味着什么。MATLAB和Fortran对数组使用列优先的存储顺序。C/C++对数组使用行优先的内存顺序。在MATLAB中:
创建一个2x 3x 4的3D数组,其中数字1到24按内存顺序依次排列:

>> x = reshape(1:24,2,3,4)
x(:,:,1) =
     1     3     5
     2     4     6
x(:,:,2) =
     7     9    11
     8    10    12
x(:,:,3) =
    13    15    17
    14    16    18
x(:,:,4) =
    19    21    23
    20    22    24

演示内存顺序:

>> x(:)
ans =
     1
     2
     3
     4
     5
     6
     7
     8
     9
    10
    11
    12
    13
    14
    15
    16
    17
    18
    19
    20
    21
    22
    23
    24

现在将其传递到C-mex例程中进行演示:
mex代码(bare bones no argument checking):

// File array3d_naive.c

/* Includes ----------------------------------------------------------- */

#include "mex.h"
#include <string.h> // for memcpy( )

/* Gateway ------------------------------------------------------------ */

void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
    // A 3D array of shape 2x3x4
    double array[2][3][4];
    void *ptarget, *psource;
    psource = mxGetData(prhs[0]);
    ptarget = array;
    // Copy the data into the array in MATLAB memory order
    memcpy(ptarget,psource,2*3*4*sizeof(double));
    // Print out the (2,2,2) element, [1][1][1] in C
    mexPrintf("element [1][1][1] = %f\n",array[1][1][1]);
}

编译、运行和比较:

>> mex array3D_naive.c
Building with 'MinGW64 Compiler (C)'.
MEX completed successfully.
>> x(2,2,2)
ans =
    10
>> array3D_naive(x)
element [1][1][1] = 18.000000

不符合预期
现在在MATLAB端使用permute,通过反转维度将数据存储顺序更改为行优先:

>> y = permute(x,[3,2,1])
y(:,:,1) =
     1     3     5
     7     9    11
    13    15    17
    19    21    23
y(:,:,2) =
     2     4     6
     8    10    12
    14    16    18
    20    22    24

y是一个4x 3x 2形状的数组,并且数据的存储器顺序已经改变(它 * 不 * 只是整形)。例如,数据的新存储顺序:

>> y(:)
ans =
     1
     7
    13
    19
     3
     9
    15
    21
     5
    11
    17
    23
     2
     8
    14
    20
     4
    10
    16
    22
     6
    12
    18
    24

将其传递到我们的C-mex例程中:

>> array3D_naive(y)
element [1][1][1] = 10.000000

现在它匹配。MATLAB permute()函数通过反转任意形状数组的维度,将列为主的内存顺序更改为行为主的内存顺序(反之亦然)。在MATLAB方面,数组的形状当然会改变,因为我们颠倒了维度。但是MATLAB端的原始2x 3x 4形状(column-major)匹配C端2x 3x 4上所需的相同形状,以便您可以以类似的方式进行索引(基于1与基于0的除外),我们只需首先使用permute()函数重新排列内存。无需额外的压平或整形。
这就是从列为主切换到行为主(反之亦然)对我来说意味着什么。如果你不想重新排列内存,那么你通常会在C中将数组设置为反向维度,例如。在这种情况下,这将是数组[4][3][2]。
我对Python和NumPy不够熟悉,无法进行这种类比。

fsi0uk1n

fsi0uk1n2#

xre = permute(reshape(x, flip(size(x))), flip(1:ndims(x)));

“行优先”是“外部优先”,而“列优先”是“内部优先”,因此要将一个N维Tensor转换为另一个,我们首先展平,然后反向填充,然后翻转维度的顺序。

逐行+示例

xF = reshape(1:2*3*4*5, [2 3 4 5]);

temp = xF(:);  % flatten (already done by `reshape`, can skip)
temp = reshape(temp, flip(size(xF)));     % fill `5` first, then `4`, ...
xC   = permute(temp, flip(1:ndims(xF)));  % shape is `[5 4 3 2]`, need `[2 3 4 5]`

生成NumPy(Python,C order)引用:

import numpy as np
xN = np.arange(1, 2*3*4*5+1).reshape(2, 3, 4, 5)

比较:

{xC(1:2,1,1,1), xC(1,1:2,1,1), xC(1,1,1:2,1), xC(1,1,1,1:2)}
%     {[1 61]}    {[1 21]}    {[1 6]}    {[1 2]}
(xN[0:2,0,0,0], xN[0,0:2,0,0], xN[0,0,0:2,0], xN[0,0,0,0:2])
# (array([ 1, 61]), array([ 1, 21]), array([1, 6]), array([1, 2]))

(我在MAT中使用了额外的代码来避免显示问题,例如:1x1x1x2 double

C转Fortran

xre = reshape(permute(x, flip(1:ndims(x))), size(x));

相关问题