gpt4 book ai didi

c++ - Python + alglib + NumPy : how to avoid converting arrays to lists?

转载 作者:太空狗 更新时间:2023-10-29 20:21:59 37 4
gpt4 key购买 nike

上下文:我最近发现了 alglib library (用于数值计算),这似乎是我一直在寻找的东西(稳健插值、数据分析......),但在 numpy 或 scipy 中找不到。

但是,我担心的事实是(例如,对于插值)它不接受 numpy 数组作为有效输入格式,而是常规 python 列表对象。

问题:我深入研究了代码和文档,发现(正如预期的那样)这个列表格式只是为了转换,因为库无论如何都会将它转换成 ctypes(cpython 库只是底层 C/C++ 的接口(interface)库)。

这就是我担心的地方:在我的代码中,我正在使用 numpy 数组,因为它大大提高了我在其上执行的科学计算的性能。因此,我担心必须将传递给 alglib 例程的任何数据转换为列表(将转换为 ctypes)会对性能产生巨大影响(我正在使用内部可能有数十万个 float 的数组,以及数千个数组)。

问题:你认为我确实会有性能损失,还是你认为我应该开始修改 alglib 代码(仅 python 接口(interface))以便它可以接受 numpy 数组,并且只进行一次转换(从 numpy 数组到 ctypes)?我什至不知道这是否可行,因为它是一个相当大的图书馆......也许你们有更好的想法或建议(即使是关于相似但不同的图书馆)...


编辑

似乎我的问题没有引起太多兴趣,或者我的问题不明确/不相关。或者也许没有人有解决方案或建议,但我怀疑周围有这么多专家:)不管怎样,我已经写了一个小的、快速的、肮脏的测试代码来说明这个问题……

#!/usr/bin/env python

import xalglib as al
import timeit
import numpy as np

def func(x):
return (3.14 *x**2.3 + x**3 -x**2.34 +x)/(1.+x)**2

def fa(x, y, val=3.14):
s = al.spline1dbuildakima(x, y)
return (al.spline1dcalc(s, val), func(val))

def fb(x, y, val=3.14):
_x = list(x)
_y = list(y)
s = al.spline1dbuildakima(_x, _y)
return (al.spline1dcalc(s, val), func(val))

ntot = 10000
maxi = 100
x = np.random.uniform(high=maxi, size=ntot)
y = func(x)
xl = list(x)
yl = list(y)

print "Test for len(x)=%d, and x between [0 and %.2f):" % (ntot, maxi)
print "Function: (3.14 *x**2.3 + x**3 -x**2.34 +x)/(1.+x)**2"
a, b = fa(xl, yl)
err = np.abs(a-b)/b * 100
print "(x=3.14) interpolated, exact =", (a, b)
print "(x=3.14) relative error should be <= 1e-2: %s (=%.2e)" % ((err <= 1e-2), err)

if __name__ == "__main__":
t = timeit.Timer(stmt="fa(xl, yl)", setup="from __main__ import fa, xl, yl, func")
tt = timeit.Timer(stmt="fb(x, y)", setup="from __main__ import fb, x, y, func")
v = 1000 * t.timeit(number=100)/100
vv = 1000 * tt.timeit(number=100)/100
print "%.2f usec/pass" % v
print "%.2f usec/pass" % vv
print "%.2f %% less performant using numpy arrays" % ((vv-v)/v*100.)

并运行它,我得到:

"""
Test for len(x)=10000, and x between [0 and 100.00):
Function: (3.14 *x**2.3 + x**3 -x**2.34 +x)/(1.+x)**2
(x=3.14) interpolated, exact = (3.686727834705164, 3.6867278531266905)
(x=3.14) relative error should be <= 1e-2: True (=5.00e-07)
25.85 usec/pass
28.46 usec/pass
10.09 % less performant using numpy arrays
"""

性能损失在 8% 到 14% 之间波动,这对我来说是巨大的......

最佳答案

您可以创建自己的 wrap 函数,将 numpy 数组的数据缓冲区直接传递给 vector 的数据指针,这不会复制数据,并且可以大大加快 wrap 函数的速度。以下代码将 x.ctypes.data 传递给 x_vector.ptr.p_ptr,其中 x 是一个 numpy 数组。

当你传递 numpy 数组时,你必须确保数组的元素在连续的内存中。以下代码不检查此项。

import xalglib as al
import numpy as np
import ctypes

def spline1dbuildakima(x, y):
n = len(x)
_error_msg = ctypes.c_char_p(0)
__c = ctypes.c_void_p(0)
__n = al.c_ptrint_t(n)
__x = al.x_vector(cnt=n, datatype=al.DT_REAL, owner=al.OWN_CALLER,
last_action=0,ptr=al.x_multiptr(p_ptr=x.ctypes.data))
__y = al.x_vector(cnt=n, datatype=al.DT_REAL, owner=al.OWN_CALLER,
last_action=0,ptr=al.x_multiptr(p_ptr=y.ctypes.data))

al._lib_alglib.alglib_spline1dbuildakima(
ctypes.byref(_error_msg),
ctypes.byref(__x),
ctypes.byref(__y),
ctypes.byref(__n),
ctypes.byref(__c))

__r__c = al.spline1dinterpolant(__c)
return __r__c

def func(x):
return (3.14 *x**2.3 + x**3 -x**2.34 +x)/(1.+x)**2

def fa(x, y, val=3.14):
s = spline1dbuildakima(x, y)
return al.spline1dcalc(s, val), func(val)

def fb(x, y, val=3.14):
s = al.spline1dbuildakima(x, y)
return al.spline1dcalc(s, val), func(val)

ntot = 10000
maxi = 100
x = np.random.uniform(high=maxi, size=ntot)
y = func(x)
xl = list(x)
yl = list(y)

import time
start = time.clock()
for i in xrange(100):
a, b = fa(x, y)
print time.clock()-start
err = np.abs(a-b)/b * 100
print a, b, err

start = time.clock()
for i in xrange(100):
a, b = fb(xl, yl)
print time.clock()-start
err = np.abs(a-b)/b * 100
print a, b, err

输出是:

0.722314760822 <- seconds of numpy array version
3.68672728107 3.68672785313 1.55166878281e-05
3.22011891502 <- seconds of list version
3.68672728107 3.68672785313 1.55166878281e-05

关于c++ - Python + alglib + NumPy : how to avoid converting arrays to lists?,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/9448915/

37 4 0
Copyright 2021 - 2024 cfsdn All Rights Reserved 蜀ICP备2022000587号
广告合作:1813099741@qq.com 6ren.com