10

简而言之,我试图从 python 调用一个共享库,更具体地说,是从 numpy 调用。共享库是在 C 中使用 sse2 指令实现的。启用优化,即使用 -O2 或 -O1 构建库时,我在通过 ctypes 调用共享库时遇到了奇怪的段错误。禁用优化 (-O0),一切都按预期进行,就像将库直接链接到 c 程序时一样(优化与否)。附上你会发现一个片段,它展示了我系统上的描述行为。启用优化后,gdb 会在 emmintrin.h:113 的 __builtin_ia32_loadupd (__P) 中报告段错误。__P 的值被报告为优化。

测试.c:

#include <emmintrin.h>
#include <complex.h>
void test(const int m, const double* x, double complex* y) {

    int i;
    __m128d _f, _x, _b;
    double complex f __attribute__( (aligned(16)) );
    double complex b __attribute__( (aligned(16)) );
    __m128d* _p;

    b = 1;
    _b = _mm_loadu_pd( (double *) &b );

    _p = (__m128d*) y;

    for(i=0; i<m; ++i) {
        f = cexp(-I*x[i]);
        _f = _mm_loadu_pd( (double *) &f );
        _x = _mm_loadu_pd( (double *) &x[i] );      
        _f = _mm_shuffle_pd(_f, _f, 1);
        *_p = _mm_add_pd(*_p, _f);
        *_p = _mm_add_pd(*_p, _x); 
        *_p = _mm_mul_pd(*_p,_b);
        _p++;
    }
    return;
}

编译器标志: gcc -o libtest.so -shared -std=c99 -msse2 -fPIC -O2 -g -lm test.c

测试.py:

import numpy as np
import os

def zerovec_aligned(nr, dtype=np.float64, boundary=16):
    '''Create an aligned array of zeros.
    '''
    size = nr * np.dtype(dtype).itemsize
    tmp = np.zeros(size + boundary, dtype=np.uint8)
    address = tmp.__array_interface__['data'][0]
    offset = boundary - address % boundary
    return tmp[offset:offset + size].view(dtype=dtype)


lib = np.ctypeslib.load_library('libtest', '.' )
lib.test.restype = None
lib.test.argtypes = [np.ctypeslib.ctypes.c_int,
                     np.ctypeslib.ndpointer(np.float64, flags=('C', 'A') ),
                     np.ctypeslib.ndpointer(np.complex128, flags=('C', 'A', 'W') )]


n = 13
y = zerovec_aligned(n, dtype=np.complex128)
x = np.ones(n, dtype=np.float64)
# x = zerovec_aligned(n, dtype=np.float64)
# x[:] = 1.

lib.test(n,x,y)

从 C 调用 test 按预期工作:

call_from_c.c:

#include <stdio.h>
#include <complex.h>
#include <stdlib.h>
#include <emmintrin.h>

void test(const int m, const double* x, double complex* y);

int main() {

    int i; 
    const int n = 11;
    double complex *y = (double complex*) _mm_malloc(n*sizeof(double complex), 16);
    double *x = (double *) malloc(n*sizeof(double));
    for(i=0; i<n; ++i) {
        x[i] = 1;
        y[i] = 0;
    }

    test(n, x, y);
    for(i=0; i<n; ++i)
            printf("[%f %f]\n", creal(y[i]), cimag(y[i]));

    return 1;

}

编译调用:
gcc -std=c99 -otestc -msse2 -L。-ltest call_from_c.c
导出 LD_LIBRARY_PATH=${LD_LIBRARY_PATH}:.
./testc
... 有效。

我的系统:

  • Ubuntu Linux i686 2.6.31-22-generic
  • 编译器:gcc (Ubuntu 4.4.1-4ubuntu9)
  • Python:Python 2.6.4(r264:75706,2009 年 12 月 7 日,18:45:15)[GCC 4.4.1]
  • 麻木:1.4.0

我已经采取了规定(参见python代码)y是对齐的并且x的对齐应该无关紧要(我认为;虽然明确对齐x并不能解决问题)。

另请注意,我在加载 b 和 f 时使用 _mm_loadu_pd 而不是 _mm_load_pd。对于 C-only 版本 _mm_load_pd 工作(如预期)。但是,当使用 _mm_load_pd 通过 ctypes 调用函数时,总是会出现段错误(与优化无关)。

我已经尝试了几天来解决这个问题但没有成功......而且我正处于将我的显示器殴打致死的边缘。欢迎任何输入。丹尼尔

4

3 回答 3

2

我只是被这个试图从 python 调用一些 SSE 代码所困扰,问题似乎是 GCC 想要假设堆栈在 16 字节边界上对齐(架构上最大的本机类型,即 SSE 类型),并使用该假设计算所有偏移量。当该假设为假时,SSE 指令将陷入困境。

答案似乎是编译

gcc -mstackrealign
这将函数序言更改为始终将堆栈对齐到 16 个字节。

于 2011-06-20T22:19:04.463 回答
1

尝试使用 numpy 构建系统构建您的扩展,以消除潜在的 cflags/ldflags 差异: http ://projects.scipy.org/numpy/wiki/NumpySconsExtExamples

于 2010-08-28T20:32:44.050 回答
-1

您是否尝试升级到 Numpy 1.5.0b2。只需运行以下命令(但要小心它可能会破坏其他东西(您必须重新编译所有 pyrex):

sudo easy_install -U numpy

当我尝试使用 H5PY 时,我遇到了与 ctypes 类似的问题(我必须重新编译 .deb 才能使用最新版本的 numpy),并且最新升级修复了 weave 的主要问题。

于 2010-08-28T03:22:01.670 回答