Python による Fortran の手続きの呼び出し2016年8月31日 | |
はじめにPython による Fortran の手続きの呼び出しについて。 環境
方針Fortran の手続きを呼び出すため、C でラッパー関数を書いて、CFFI で呼び出す。 コードサンプルコードは以下の通り。 test.f subroutine test(n, a) implicit none integer :: n double precision :: dx double precision, parameter::pi=3.1415926535897932 double precision, dimension(2, n) :: a double precision :: x, f integer i dx = 2*pi/(n-1) do i = 1, n x = i*dx f = sin(x) a(1, i) = x a(2, i) = f enddo end test_f.h #ifndef TEST_H #define TEST_H void test_(int *n, double a[][2]); #endif /* TEST_H */ test.c #include "test_f.h" void test(int n, double a[][2]) { test_(&n, a); } Makefile (Windows) all: gfortran -Wall -O2 -o test_f.o -c test.f gcc -Wall -O2 -c test.c gcc -shared -o test.so test_f.o test.o clean: rm -f *.o Makefile (Linux) all: gfortran -fPIC -Wall -O2 -o test_f.o -c test.f gcc -fPIC -Wall -O2 -c test.c gcc -shared -o test.so test_f.o test.o clean: rm -f *.o x = [0, 2π] の範囲で sin 関数を計算した結果をリストで返すものである。 注意点
Python から C の関数を呼び出すPython から CFFI を用いて C の関数を呼び出す。 call_c_func_cffi.py from matplotlib import pyplot as plt import numpy as np from cffi import FFI ffi = FFI() ffi.cdef(r''' void test(int n, double a[][2]); ''') lib = ffi.dlopen('./test.so') def test(n): _r = ffi.new('double[%d][2]' % n) lib.test(n, _r) r = np.array([[elem[i] for i in (0, 1)] for elem in _r]) return r r = test(20) for i in range(20): print('%g %g' % (r[i, 0], r[i, 1])) plt.plot(r[:, 0], r[:, 1]) plt.show() Python から Fortran の手続きを呼び出す実は CFFI で直接 Fortran の手続きを呼び出すこともできる (Fortran のコードを C でラッピングできているのだから、当然)。 Makefile を次のように書き換える。 Makefile (Windows) all: gfortran -Wall -O2 -o test_f.o -c test.f gcc -shared -o test.so test_f.o clean: rm -f *.o Makefile (Linux) all: gfortran -fPIC -Wall -O2 -o test_f.o -c test.f gcc -shared -o test.so test_f.o clean: rm -f *.o Fortran の手続きを呼び出すコードは以下の通り。 call_f_func_cffi.py from matplotlib import pyplot as plt import numpy as np from cffi import FFI ffi = FFI() ffi.cdef(r''' void test_(int *n, double a[][2]); ''') lib = ffi.dlopen('./test.so') def test(n): _n = ffi.new('int[]', [n]) _r = ffi.new('double[%d][2]' % n) lib.test_(_n, _r) r = np.array([[elem[i] for i in (0, 1)] for elem in _r]) return r r = test(20) for i in range(20): print('%g %g' % (r[i, 0], r[i, 1])) plt.plot(r[:, 0], r[:, 1]) plt.show() 関数名にアンダースコアを付ける必要があること、引数をすべて配列の形にしなければならないところが少しややこしいかもしれない。 | |
PENGUINITIS |