はじめに
今回は, ctypesを使ってPythonから自作のCの関数を呼び出します.
Cの関数の引数としてdouble**型の変数を与え, 結果をdouble**型でPythonに返す方法についての日本語記事がなかったので簡単に調べました.
なので, 正確には多次元配列のやりとりではなくdouble**型変数のやりとりです(小声)
コード
まず以下のような自作の関数を用意します.
今回は引数の行列の全要素にnを足すとします.
libadd.c
#include <stdio.h>
void add_matrix(double **matrix, int row, int col, int n) {
int i, j;
for(i=0; i<row; i++){
for(j=0; j<col; j++){
matrix[i][j] = matrix[i][j] + n;
}
}
}
そして以下のとおりにコンパイル.
$ gcc -cpp -fPIC -shared libadd.c -lm -o libadd.so -O3
-sharedオプションをつけてlibadd.soを作る必要があります.
続いてPython側,
test.py
from ctypes import *
import ctypes.util
from numpy.ctypeslib import ndpointer
import numpy as np
from numpy.random import *
#さっきつくったlibadd.soファイルを指定
lib = np.ctypeslib.load_library("libadd.so",".")
#適当につくります
row = 20
col = 5
n = 5
matrix = rand(row, col)
#doubleのポインタのポインタ型を用意
_DOUBLE_PP = ndpointer(dtype=np.uintp, ndim=1, flags='C')
#add_matrix()関数の引数の型を指定(ctypes)
lib.add_matrix.argtypes = [_DOUBLE_PP, c_int32, c_int32, c_int32]
#add_matrix()関数が返す値の型を指定(今回は返り値なし)
lib.add_matrix.restype = None
#おまじない
tp = np.uintp
mpp = (matrix.__array_interface__['data'][0] + np.arange(matrix.shape[0])*matrix.strides[0]).astype(tp)
#int型もctypeのc_int型へ変換して渡す
n = ctypes.c_int(n)
row = ctypes.c_int(row)
col = ctypes.c_int(col)
print("before:", matrix)
lib.add_matrix(mpp, row, col, n)
print("after:", matrix)
こんな感じで用意してやります.
気をつけること（自分がハマったところ）は以下です.
- lib.add_matrix.argtypes()とlib.add_matrix.restype()のadd_matrix部分は自分が呼び出したい関数名を指定すること.
- integer型の変数もctypes用の型に変換すること.
- Python側では, mppではなくmatrixに結果が格納されること. (Cのadd_matrix()関数の返り値に設定しなくとも良い)
また_DOUBLE_PPはPOINTER(POINTER(c_double))でいけるかと思ったのですが, C側へうまく参照渡しされませんでした.(他に良い方法などあればご教示ください.)
実行結果は以下のようになります.
実行結果
$ python test.py
before: [[ 0.38897722 0.46313741 0.84976397 0.64423638 0.95435434]
[ 0.60333226 0.37742093 0.3465723 0.77567722 0.74132703]
[ 0.87302221 0.43188915 0.1211892 0.03272387 0.16480218]
[ 0.60387532 0.37643299 0.72456328 0.67993591 0.6592758 ]
[ 0.63544435 0.7261329 0.16825874 0.11123031 0.97664729]
[ 0.58328763 0.85650723 0.77429203 0.12418374 0.44881621]
[ 0.25106831 0.22417529 0.22837031 0.17608283 0.49529552]
[ 0.75357385 0.3937937 0.99262338 0.23625638 0.37777575]
[ 0.887013 0.85391559 0.32386967 0.30482811 0.67282479]
[ 0.96701627 0.23598535 0.29455679 0.37472803 0.28614415]
[ 0.82208389 0.87058549 0.21354254 0.18766391 0.13951804]
[ 0.22836767 0.4051488 0.79205244 0.17489158 0.83368298]
[ 0.99464033 0.91165387 0.48920652 0.69251231 0.3369773 ]
[ 0.39567888 0.08322727 0.50138724 0.00972839 0.53206287]
[ 0.4981778 0.81622292 0.73376253 0.19702753 0.56725497]
[ 0.37709128 0.8219189 0.59531858 0.60284025 0.75683106]
[ 0.11898215 0.93335749 0.39432259 0.9733907 0.355181 ]
[ 0.54600745 0.36843876 0.52919308 0.22250149 0.66595212]
[ 0.33897066 0.68039774 0.89937153 0.77007464 0.57667845]
[ 0.08015923 0.42939639 0.82431951 0.05494513 0.89979669]]
after: [[ 5.38897722 5.46313741 5.84976397 5.64423638 5.95435434]
[ 5.60333226 5.37742093 5.3465723 5.77567722 5.74132703]
[ 5.87302221 5.43188915 5.1211892 5.03272387 5.16480218]
[ 5.60387532 5.37643299 5.72456328 5.67993591 5.6592758 ]
[ 5.63544435 5.7261329 5.16825874 5.11123031 5.97664729]
[ 5.58328763 5.85650723 5.77429203 5.12418374 5.44881621]
[ 5.25106831 5.22417529 5.22837031 5.17608283 5.49529552]
[ 5.75357385 5.3937937 5.99262338 5.23625638 5.37777575]
[ 5.887013 5.85391559 5.32386967 5.30482811 5.67282479]
[ 5.96701627 5.23598535 5.29455679 5.37472803 5.28614415]
[ 5.82208389 5.87058549 5.21354254 5.18766391 5.13951804]
[ 5.22836767 5.4051488 5.79205244 5.17489158 5.83368298]
[ 5.99464033 5.91165387 5.48920652 5.69251231 5.3369773 ]
[ 5.39567888 5.08322727 5.50138724 5.00972839 5.53206287]
[ 5.4981778 5.81622292 5.73376253 5.19702753 5.56725497]
[ 5.37709128 5.8219189 5.59531858 5.60284025 5.75683106]
[ 5.11898215 5.93335749 5.39432259 5.9733907 5.355181 ]
[ 5.54600745 5.36843876 5.52919308 5.22250149 5.66595212]
[ 5.33897066 5.68039774 5.89937153 5.77007464 5.57667845]
[ 5.08015923 5.42939639 5.82431951 5.05494513 5.89979669]]
やっためう！すごいめう！
これやるだけでも, 処理の内容によっては相当な高速化になります.
これでGILから解放されるとかどうとか読んだ気がするのでマルチスレッディング使いたいし, もしくはマルチプロセッシングかハイブリッド並列で使いたいのだけど...今のところうまくいってません.