前回のこのブログの続きをやる。Cのソースコードがpython2用に書かれているので、このサイトを参考にしてpython3用に書き直した。
スポンサーリンク
Cコードのpython2から3への書き換え¶
書き換え前のこのサイトから拝借したソースコード。
#include "Python.h"
/* make this static if you don't want other code to call this function */
/* I don't make it static because want to access this via ctypes */
/* static */
int iterate_point(double x0, double y0, int max_iterations) {
int iteration = 0;
double x=x0, y=y0, x2=x*x, y2=y*y;
while (x2+y2<4 && iteration<max_iterations) {
y = 2*x*y + y0;
x = x2-y2 + x0;
x2 = x*x;
y2 = y*y;
iteration++;
}
return iteration;
}
/* The module doc string */
PyDoc_STRVAR(mandelbrot__doc__,
"Mandelbrot point evalutation kernel");
/* The function doc string */
PyDoc_STRVAR(iterate_point__doc__,
"x,y,max_iterations -> iteration count at that point, up to max_iterations");
/* The wrapper to the underlying C function */
static PyObject *
py_iterate_point(PyObject *self, PyObject *args)
{
double x=0, y=0;
int iterations, max_iterations=1000;
/* "args" must have two doubles and may have an integer */
/* If not specified, "max_iterations" remains unchanged; defaults to 1000 */
/* The ':iterate_point' is for error messages */
if (!PyArg_ParseTuple(args, "dd|i:iterate_point", &x, &y, &max_iterations))
return NULL;
/* Verify the parameters are correct */
if (max_iterations < 0) max_iterations = 0;
/* Call the C function */
iterations = iterate_point(x, y, max_iterations);
/* Convert from a C integer value to a Python integer instance */
return PyInt_FromLong((long) iterations);
}
/* A list of all the methods defined by this module. */
/* "iterate_point" is the name seen inside of Python */
/* "py_iterate_point" is the name of the C function handling the Python call */
/* "METH_VARGS" tells Python how to call the handler */
/* The {NULL, NULL} entry indicates the end of the method definitions */
static PyMethodDef mandelbrot_methods[] = {
{"iterate_point", py_iterate_point, METH_VARARGS, iterate_point__doc__},
{NULL, NULL} /* sentinel */
};
/* When Python imports a C module named 'X' it loads the module */
/* then looks for a method named "init"+X and calls it. Hence */
/* for the module "mandelbrot" the initialization function is */
/* "initmandelbrot". The PyMODINIT_FUNC helps with portability */
/* across operating systems and between C and C++ compilers */
PyMODINIT_FUNC
initmandelbrot(void)
{
/* There have been several InitModule functions over time */
Py_InitModule3("mandelbrot", mandelbrot_methods,
mandelbrot__doc__);
}
上のコードはpython2用に書かれているので、これをpython3用のコードに書き換える。
#include "Python.h"
/* make this static if you don't want other code to call this function */
/* I don't make it static because want to access this via ctypes */
/* static */
int iterate_point(double x0, double y0, int max_iterations) {
int iteration = 0;
double x=x0, y=y0, x2=x*x, y2=y*y;
while (x2+y2<4 && iteration<max_iterations) {
y = 2*x*y + y0;
x = x2-y2 + x0;
x2 = x*x;
y2 = y*y;
iteration++;
}
return iteration;
}
/* The module doc string */
PyDoc_STRVAR(mandelbrot__doc__,
"Mandelbrot point evalutation kernel");
/* The function doc string */
PyDoc_STRVAR(iterate_point__doc__,
"x,y,max_iterations -> iteration count at that point, up to max_iterations");
/* The wrapper to the underlying C function */
static PyObject *
py_iterate_point(PyObject *self, PyObject *args)
{
double x=0, y=0;
int iterations, max_iterations=1000;
/* "args" must have two doubles and may have an integer */
/* If not specified, "max_iterations" remains unchanged; defaults to 1000 */
/* The ':iterate_point' is for error messages */
if (!PyArg_ParseTuple(args, "dd|i:iterate_point", &x, &y, &max_iterations))
return NULL;
/* Verify the parameters are correct */
if (max_iterations < 0) max_iterations = 0;
/* Call the C function */
iterations = iterate_point(x, y, max_iterations);
/* Convert from a C integer value to a Python integer instance */
return PyLong_FromLong((long) iterations);
}
/* A list of all the methods defined by this module. */
/* "iterate_point" is the name seen inside of Python */
/* "py_iterate_point" is the name of the C function handling the Python call */
/* "METH_VARGS" tells Python how to call the handler */
/* The {NULL, NULL} entry indicates the end of the method definitions */
static PyMethodDef mandelbrot_methods[] = {
{"iterate_point", py_iterate_point, METH_VARARGS, iterate_point__doc__},
{NULL, NULL} /* sentinel */
};
static struct PyModuleDef mandelbrot_definition = {
PyModuleDef_HEAD_INIT,
"mandelbrot",
"A Python module that plots mandelbrot from C code.",
-1,
mandelbrot_methods
};
/* When Python imports a C module named 'X' it loads the module */
/* then looks for a method named "init"+X and calls it. Hence */
/* for the module "mandelbrot" the initialization function is */
/* "initmandelbrot". The PyMODINIT_FUNC helps with portability */
/* across operating systems and between C and C++ compilers */
PyMODINIT_FUNC
PyInit_mandelbrot(void)
{
/* There have been several InitModule functions over time */
Py_Initialize();
return PyModule_Create(&mandelbrot_definition);
}
スポンサーリンク
ソースファイルのコンパイル¶
!python3 setup.py build_ext --inplace
スポンサーリンク
ライブラリをインポート¶
import ctypes
import mandelbrot
ctypes_iterate_point = ctypes.CDLL\
("./mandelbrot.cpython-36m-x86_64-linux-gnu.so").iterate_point
ctypes_iterate_point.restype = ctypes.c_int
ctypes_iterate_point.argtypes = \
[ctypes.c_double, ctypes.c_double, ctypes.c_int]
x = -2
y = -1
w = 2.5
h = 2.0
NY = 40
NX = 70
RANGE_Y = range(NY)
RANGE_X = range(NX)
def render(iterate_point):
chars = []
append = chars.append
for j in RANGE_Y:
for i in RANGE_X:
it = iterate_point(x+w/NX*i, y+h/NY*j, 1000)
if it == 1000:
append("*")
elif it > 5:
append(",")
elif it > 2:
append(".")
else:
append(" ")
append("\n")
return "".join(chars)
import time
t1 = time.time()
s1 = render(mandelbrot.iterate_point)
t2 = time.time()
s2 = render(ctypes_iterate_point)
t3 = time.time()
assert s1 == s2
print (s1)
print ("as C extension", t2-t1)
print ("with ctypes", t3-t2)
このエラーはPyInt_FromLongをPyLong_FromLongに変えればいいらしい。
import ctypes
import mandelbrot
ctypes_iterate_point = ctypes.CDLL\
("./mandelbrot.cpython-36m-x86_64-linux-gnu.so").iterate_point
ctypes_iterate_point.restype = ctypes.c_int
ctypes_iterate_point.argtypes = \
[ctypes.c_double, ctypes.c_double, ctypes.c_int]
x = -2
y = -1
w = 2.5
h = 2.0
NY = 40
NX = 70
RANGE_Y = range(NY)
RANGE_X = range(NX)
def render(iterate_point):
chars = []
append = chars.append
for j in RANGE_Y:
for i in RANGE_X:
it = iterate_point(x+w/NX*i, y+h/NY*j, 1000)
if it == 1000:
append("*")
elif it > 5:
append(",")
elif it > 2:
append(".")
else:
append(" ")
append("\n")
return "".join(chars)
import time
t1 = time.time()
s1 = render(mandelbrot.iterate_point)
t2 = time.time()
s2 = render(ctypes_iterate_point)
t3 = time.time()
assert s1 == s2
print (s1)
print ("as C extension", t2-t1)
print ("with ctypes", t3-t2)
0.006439924240112305/0.005004167556762695
C extensionの方がctypesよりも1.29倍高速ではあったが、速度差はほとんど消失していると言ってもいいだろう。
スポンサーリンク
スポンサーリンク