matlab

Matlab中基于共享内存的矩阵 Part 2

目录#

共享内存的操作#

Windows下可以用Windows的API进行操作，而Linux下则可以通过POSIX的API进行操作，虽然过程有点区别，但终究是大同小异。

创建&数据写入#

为了简单测试，这里就直接用Windows版Python自带的mmap当实验：

1
2
3
4
5
6
7
8
9
10
11
12
13
>>> import numpy as np
>>> a = np.random.randn(32*4096,4096)  # 先占用4G内存
>>> a[0,0]
-0.5271477716515403
>>> a[-1,-1]
1.1194140920933267
>>> b = a.tobytes()  # 额外的4G内存，peak 8G
>>> del a  # 释放4G内存
>>> import mmap
>>> m = mmap.mmap(0, 32*4096*4096*8, 'Local\\MyShMem1')
>>> m.write(b)  # 额外的4G内存，peak 8G
>>> del b  # 释放4G内存
>>>

这段代码会常占4G的内存，顶峰的内存会达到数据量的2倍。（在写入共享内存的时候）

注：Linux版Python不支持这种创建共享内存的方式。

数据读取#

数据读取是用matlab的mex写的：

shmem_win_data_read.c：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
#include <mex.h>
#include <stdio.h>
#include <windows.h>
#pragma comment(lib, "user32.lib")

TCHAR szName[] = TEXT("Local\\MyShMem1");
void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) {
    //mxUint64 ptr = mxGetUint64s(prhs[0])[0];
    mxArray* pArr = mxCreateDoubleMatrix(32*4096, 4096, 0);
    HANDLE hMapFile = OpenFileMapping(FILE_MAP_ALL_ACCESS, FALSE, szName); // param: R/W, do not inherit the name, name of mapping object
    if (hMapFile == NULL) {
        printf("Failed to open file mapping: %d\n", GetLastError());
        return;
    }
    double* b = (double*)MapViewOfFile(hMapFile, FILE_MAP_ALL_ACCESS, 0, 0, 32*4096*4096*sizeof(double)); // Last: bufsize
    if (b == NULL) {
        printf("Failed to map view of file: %d\n", GetLastError());
        return;
    }
    for(int i = 0; i < 5; i++)
        printf("%lf\n", b[i]);
    double* originalArr = mxGetPr(pArr);
    mxFree(originalArr);
    printf("Original ptr: %lld\n", originalArr);
    mxSetPr(pArr, b);
    printf("Current ptr: %lld\n", mxGetPr(pArr));
    plhs[0] = pArr;
    if (nlhs > 1) {
        // handle
        mwSize d[] = { 1 };
        pArr = mxCreateNumericArray(1, d, mxUINT64_CLASS, mxREAL);
        *(unsigned long long*)mxGetData(pArr) = (unsigned long long)hMapFile;
        plhs[1] = pArr;
    }
}

然后编译，运行：

1
2
3
4
5
mex 'shmem_win_data_read.c' -g
a = cell(1);
[a{1}, h] = shmem_win_data_read;
b = a{1};
% 用变量b进行正常读写

返回的第一个参数是数据，至于为什么要用cell，就是为了规避matlab自带的GC机制。一旦将数组的指针利用mxSetPr设置为一个非通过调用mxMalloc得到的值（如c语言自带的malloc），GC机制尝试对它进行free的时候，后果自己可以想象。

返回的第二个参数是共享内存的handle，需要在使用完后释放掉。

cell的内存机制后面再说，注意的是，在调用mxSetPr之前，一定要用mxGetPr和mxFree把最初创建返回数组pArr时所申请到的内存释放掉，否则就会发生内存泄漏。

可以看出，数据是完全一致的。

资源释放#

一旦用完，就该妥善处理后事了（笑）。你大可以试试clear直接清掉，看看matlab会不会崩。

shmem_win_data_free.c：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
#include <mex.h>
#include <stdio.h>
#include <windows.h>
#pragma comment(lib, "user32.lib")
#include <matrix.h>
#include <stdlib.h>

void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) {
    mxArray* cell_arr = mxGetCell(prhs[0], 0);
    double* new_ptr = (double *)mxMalloc(sizeof(double) * 5);
    double* original_ptr = mxGetPr(cell_arr);
    new_ptr[0] = 1.2;
    new_ptr[1] = 0.0;
    new_ptr[2] = 0.0;
    new_ptr[3] = 0.0;
    printf("new ptr: %lld\n", new_ptr);
    printf("original ptr: %lld\n", original_ptr);
    mxSetPr(cell_arr, new_ptr);
    printf("current ptr: %lld\n", mxGetPr(cell_arr));
    mxSetM(cell_arr, 2);
    mxSetN(cell_arr, 2);
    UnmapViewOfFile(original_ptr);
    if (nrhs > 1) {
        HANDLE hMapFile = *(HANDLE*)mxGetData(prhs[1]);
        CloseHandle(hMapFile);
    }
}

为了简单实验起见，这里就直接用把内存改回matlab的托管内存了（用mxMalloc申请的内存），然后手动把之前改的内存释放掉。

1
2
mex 'shmem_win_data_free.c' -g
shmem_win_data_free(a, h);  % 这里就用a作为参数传入

合理利用matlab的“特性”，a{1}和b都会同时修改为一个新的数组，接着就可以直接clear掉，不需要再担心什么了。

完整使用#

1
2
3
4
5
6
7
8
9
v = zeros(1, 4096);
parfor x = 1:4096
    a = cell(1);
    [a{1}, h] = shmem_win_data_read;
    b = a{1};
    c = sum(b(:, x));
    v(x) = c;
    shmem_win_data_free(a, h);
end

这时候大可安心看系统的内存使用情况，妈妈再也不用担心内存占用啦。

番外：直接使用mxSetPr的效果#

mxSetPr在mex中一般用法效果如下：

1
2
3
4
5
6
7
8
9
10
11
12
#include <mex.h>
#include <stdio.h>

void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) {
    const mxArray* arg = prhs[0];
    double* ptr = mxGetPr(arg);
    printf("Original dataArray PTR: 0x%p\n", ptr);
    double* new_ptr = (double*)mxMalloc(sizeof(double) * mxGetM(arg) * mxGetN(arg));
    printf("Allocated dataArray PTR: 0x%p\n", new_ptr);
    mxSetPr(arg, new_ptr);
    printf("New dataArray PTR: 0x%p\n", mxGetPr(arg));
}

上面代码的作用：直接申请一块新的内存，将这块内存的地址直接赋值到matlab（作为输入参数传入mex函数）的矩阵数组中。运行结果如下：

可以看到，在mex函数体内对数组进行的任何修改是无法直接影响到输入参数的。

番外2：必须用cell中的元素作为输出变量#

还有一个关键点，就是必须用a{1}当输出变量，如：

1
2
3
4
5
6
7
a{1} = randn(5);
b = a{1};
% 通过b进行数据读取

do_cleanup_within_cell(a);
% 上面的函数能够同时作用于a（输入变量）和b（通过a引用赋值）
% 执行结果：a = 1×1 cell, a{1} = [], b = []

而把一般变量作为输出，再套用cell，则只能影响到cell里面的元素，不会影响到最先赋值的变量

1
2
3
4
5
6
7
8
t = randn(5); % 先赋值到一个新变量
a{1} = t; % 再通过引用赋值给cell
b = a{1};
% 通过t或者b进行数据读取

do_cleanup_within_cell(a);
% 上面的函数仅能作用于a（输入变量）和b（通过a引用赋值），而最原始的t则不受影响
% 执行结果：t = 5×5 double, a = 1×1 cell, a{1} = [], b = []

其中do_cleanup_within_cell可为如下（直接把输入数组清空为[]）：

1
2
3
4
5
6
7
8
void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) {
    mxArray* cell_arr = mxGetCell(prhs[0], 0);
    double* original_ptr = mxGetPr(cell_arr);
    printf("original ptr: 0x%p\n", original_ptr);
    mxSetM(cell_arr, 0);
    mxSetN(cell_arr, 0);
    mxSetPr(cell_arr, NULL);
}

反正个人猜测这种现象肯定跟matlab的GC是脱不了钩的，由于没有官方的内存介绍文档，所以也不深究了。写代码讲究一件事，那就是能用就行。

数据类型#

除了double、float以及(u)int8/16/32/64外，还有complex、sparse、cell、struct和object等数据类型，一个完整的check如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
#include <mex.h>
#include <matrix.h>
#include <stdio.h>

#define CHECK_FUNC(x) printf(#x": %d\n", x(prhs[0]))

void mexFunction(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]) {
    // 1. mxArray attributes
    CHECK_FUNC(mxIsNumeric);
    CHECK_FUNC(mxIsComplex);
    // 2. create, query, and access data types
    // 2.2 noncomplex float
    CHECK_FUNC(mxIsScalar);
    CHECK_FUNC(mxIsDouble);
    CHECK_FUNC(mxIsSingle);
    // 2.3 noncomplex integer
    CHECK_FUNC(mxIsInt8);
    CHECK_FUNC(mxIsUint8);
    CHECK_FUNC(mxIsInt16);
    CHECK_FUNC(mxIsUint16);
    CHECK_FUNC(mxIsInt32);
    CHECK_FUNC(mxIsUint32);
    CHECK_FUNC(mxIsInt64);
    CHECK_FUNC(mxIsUint64);
    // 2.6 sparse
    CHECK_FUNC(mxIsSparse);
    // 2.8 character
    CHECK_FUNC(mxIsChar);
    // 2.9 logical
    CHECK_FUNC(mxIsLogical);
    // 2.10 object
    // mxIsClass(mxArray* arr, char* className)
    // 2.11 structure
    CHECK_FUNC(mxIsStruct);
    // 2.12 cell
    CHECK_FUNC(mxIsCell);
}

接下来简单地验证完Linux的共享内存之后，就可以开始造（已经重复）的轮子了。

Linux版的碎碎念#

手头上只有Win版的R2018b和Linux版的R2017b服务器，想着版本差距不算大，直接改改代码，把WinAPI的Named shared memory改成POSIX的shm_open大概就行了。

把代码扔过去，编译，运行，哦豁完蛋，报了个在Windows下未曾出现过的崩溃错误：

突然告诉我还有个16Byte的header？喂不是吧大哥，我只把数据扔共享内存上了你告诉我header？

好吧，既然你想读，那就让你读呗。我把共享内存创建大一点，把数据扔后面一点，前面留16Byte的零，总算可以吧。然后：

你塔喵还要查header是吧？好好，我把mxGetPr往前16个字节的数据一并复制过去总行了吧。

嗯，这次好像真行了，而且能够正常退出没崩。

至于header是什么，说实话自己也没弄懂。随便弄了几个矩阵，往前挪16个字节print一下，结果如下：

131072*4096大小的double矩阵的header：
00 00 00 00 01 00 00 00 CE FA ED FE 20 00 10 00

7*65537大小的double矩阵：
40 00 38 00 00 00 00 00 CE FA ED FE 20 00 20 00

7*65538大小的double矩阵：
70 00 38 00 00 00 00 00 CE FA ED FE 20 00 10 00

7*65538大小的int32/uint32/single矩阵：
40 00 1C 00 00 00 00 00 CE FA ED FE 20 00 10 00

前面8字节明显跟数组的大小相关，那么后面8字节呢，跟type也不相关啊……特别是对sparse、complex这种数据类型也不变，而且倒数第二个0x10偶尔会变成0x20，也毫无规律可言。

目前只能推测跟matlab的内部GC有关（因为改一下清空变量的时候就会崩掉）。