cython でコードをセットアップして実行する方法がわかりません。
cdef、doubleなどをコードの関連部分に追加しました。
setup.py もちろん、hello という名前は使用されていません。 cython doc
from distutils.core import setup
from distutils.extension import Extension
from Cython.Distutils import build_ext
ext_modules = [Extension("hello", ["hello.pyx"])]
setup(
name = 'Hello world app',
cmdclass = {'build_ext': build_ext},
ext_modules = ext_modules
)
編集
質問を再定義して、より具体的な例にします。
この ODE のシステムを Cython で実行したいと考えています。が最良の選択かどうかはわかりませんdoubleが、数値は 10 から 10 のオーダーであり、整数ではありません。
musしたがって、モジュールを呼び出すスクリプトで定義された変数である、より大きなコードでこれを呼び出したいと思います。pyxsetup.pyファイルをコンパイルするためのファイルがあります。私は何をする必要があるのか よくわからないので、今この歌を呼ぶことができます. モジュールに 3bodyproblem という名前を付けたとします。次に、scipt でそれを呼び出して、import 3bodyproblem「3bodyproblem.3bodyproblem(この入力は何でしょうか)」を実行します。
私は頌歌のためにcythonの紹介を読んでいますが、私の例でそれらの例を使用する方法がわかりません。また、rk 形式にする必要がある場合は、最初のコードの下にあるコードを参照してください。
コード 1
cdef deriv(double u, dt):
cdef double u[0], u[1], u[2]
return [u[3],
u[4],
u[5],
-mus * u[0] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[1] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[2] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5]
dt = np.linspace(0.0, t12sec, t12sec)
u = odeint(deriv, u0, dt, atol = 1e-13, rtol = 1e-13)
x, y, z, vx, vy, vz = u.T
コード 2
cdef deriv(dt, double u):
cdef double u[0], u[1], u[2], u[3], u[4], u[5]
return [u[3],
u[4],
u[5],
-mus * u[0] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[1] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[2] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5]
solver = ode(deriv).set_integrator('dop853')
solver.set_initial_value(u0)
z = np.zeros(300000)
for ii in range(300000):
z[ii] = solver.integrate(dt[ii])[0]
x, y, z, x2, y2, z2 = z.T
コード 1 をコンパイルしようとすると、cython はより大きな.pyファイルで定義された変数を認識しません。再コンパイルせずに他の場所で使用できるように、cython コードで変数を設定することを避けたいです。エラーは次のとおりです。
Error compiling Cython file:
------------------------------------------------------------
...
import numpy as np
from scipy.integrate import odeint
cdef deriv(double u, dt):
cdef double u[0], u[1], u[2]
^
------------------------------------------------------------
ODEcython.pyx:7:17: 'u' redeclared
Error compiling Cython file:
------------------------------------------------------------
...
import numpy as np
from scipy.integrate import odeint
cdef deriv(double u, dt):
cdef double u[0], u[1], u[2]
^
------------------------------------------------------------
ODEcython.pyx:7:23: 'u' redeclared
Error compiling Cython file:
------------------------------------------------------------
...
import numpy as np
from scipy.integrate import odeint
cdef deriv(double u, dt):
cdef double u[0], u[1], u[2]
^
------------------------------------------------------------
ODEcython.pyx:7:29: 'u' redeclared
Error compiling Cython file:
------------------------------------------------------------
...
-mus * u[0] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[1] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[2] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5]
dt = np.linspace(0.0, t12sec, t12sec)
^
------------------------------------------------------------
ODEcython.pyx:16:28: undeclared name not builtin: t12sec
Error compiling Cython file:
------------------------------------------------------------
...
-mus * u[1] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[2] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5]
dt = np.linspace(0.0, t12sec, t12sec)
u = odeint(deriv, u0, dt, atol = 1e-13, rtol = 1e-13)
^
------------------------------------------------------------
ODEcython.pyx:17:16: Cannot convert 'object (double, object)' to Python object
Error compiling Cython file:
------------------------------------------------------------
...
-mus * u[1] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
-mus * u[2] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5]
dt = np.linspace(0.0, t12sec, t12sec)
u = odeint(deriv, u0, dt, atol = 1e-13, rtol = 1e-13)
^
------------------------------------------------------------
ODEcython.pyx:17:20: undeclared name not builtin: u0
Error compiling Cython file:
------------------------------------------------------------
...
cdef deriv(double u, dt):
cdef double u[0], u[1], u[2]
return [u[3],
u[4],
u[5],
-mus * u[0] / (u[0] ** 2 + u[1] ** 2 + u[2] ** 2) ** 1.5,
^
------------------------------------------------------------
ODEcython.pyx:11:17: undeclared name not builtin: mus
building 'ODEcython' extension
creating build
creating build/temp.linux-x86_64-2.7
x86_64-linux-gnu-gcc -pthread -fno-strict-aliasing -DNDEBUG -g -fwrapv -O2 -Wall -Wstrict-prototypes -fPIC -I/usr/include/python2.7 -c ODEcython.c -o build/temp.linux-x86_64-2.7/ODEcython.o
ODEcython.c:1:2: error: #error Do not use this file, it is the result of a failed Cython compilation.
error: command 'x86_64-linux-gnu-gcc' failed with exit status 1
多分これが役立つでしょう
したがって、次のようなファイルでモジュールを使用したい: (異なる deriv 関数ですが、アイデアが同じであることは無視してください)
import numpy as np
from scipy.integrate import ode
import pylab
from mpl_toolkits.mplot3d import Axes3D
from scipy.optimize import brentq
from scipy.optimize import fsolve
me = 5.974 * 10 ** 24 # mass of the earth
mm = 7.348 * 10 ** 22 # mass of the moon
G = 6.67259 * 10 ** -20 # gravitational parameter
re = 6378.0 # radius of the earth in km
rm = 1737.0 # radius of the moon in km
r12 = 384400.0 # distance between the CoM of the earth and moon
rs = 66100.0 # distance to the moon SOI
Lambda = np.pi / 6 # angle at arrival to SOI
M = me + mm
d = 300 # distance the spacecraft is above the Earth
pi1 = me / M
pi2 = mm / M
mue = 398600.0 # gravitational parameter of earth km^3/sec^2
mum = G * mm # grav param of the moon
mu = mue + mum
omega = np.sqrt(mu / r12 ** 3)
# distance from the earth to Lambda on the SOI
r1 = np.sqrt(r12 ** 2 + rs ** 2 - 2 * r12 * rs * np.cos(Lambda))
vbo = 10.85 # velocity at burnout
h = (re + d) * vbo # angular momentum
energy = vbo ** 2 / 2 - mue / (re + d) # energy
v1 = np.sqrt(2.0 * (energy + mue / r1)) # refer to the close up of moon diagram
# refer to diagram for angles
theta1 = np.arccos(h / (r1 * v1))
phi1 = np.arcsin(rs * np.sin(Lambda) / r1)
#
p = h ** 2 / mue # semi-latus rectum
a = -mue / (2 * energy) # semi-major axis
eccen = np.sqrt(1 - p / a) # eccentricity
nu0 = 0
nu1 = np.arccos((p - r1) / (eccen * r1))
# Solving for the eccentric anomaly
def f(E0):
return np.tan(E0 / 2) - np.sqrt((1 - eccen) / (1 + eccen)) * np.tan(0)
E0 = brentq(f, 0, 5)
def g(E1):
return np.tan(E1 / 2) - np.sqrt((1 - eccen) / (1 + eccen)) * np.tan(nu1 / 2)
E1 = fsolve(g, 0)
# Time of flight from r0 to SOI
deltat = (np.sqrt(a ** 3 / mue) * (E1 - eccen * np.sin(E1)
- (E0 - eccen * np.sin(E0))))
# Solve for the initial phase angle
def s(phi0):
return phi0 + deltat * 2 * np.pi / (27.32 * 86400) + phi1 - nu1
phi0 = fsolve(s, 0)
nu = -phi0
gamma = 0 * np.pi / 180 # angle in radians of the flight path
vx = vbo * (np.sin(gamma) * np.cos(nu) - np.cos(gamma) * np.sin(nu))
# velocity of the bo in the x direction
vy = vbo * (np.sin(gamma) * np.sin(nu) + np.cos(gamma) * np.cos(nu))
# velocity of the bo in the y direction
xrel = (re + 300.0) * np.cos(nu) - pi2 * r12
yrel = (re + 300.0) * np.sin(nu)
u0 = [xrel, yrel, 0, vx, vy, 0]
def deriv(u, dt):
return [u[3], # dotu[0] = u[3]
u[4], # dotu[1] = u[4]
u[5], # dotu[2] = u[5]
(2 * omega * u[4] + omega ** 2 * u[0] - mue * (u[0] + pi2 * r12) /
np.sqrt(((u[0] + pi2 * r12) ** 2 + u[1] ** 2) ** 3) - mum *
(u[0] - pi1 * r12) /
np.sqrt(((u[0] - pi1 * r12) ** 2 + u[1] ** 2) ** 3)),
# dotu[3] = that
(-2 * omega * u[3] + omega ** 2 * u[1] - mue * u[1] /
np.sqrt(((u[0] + pi2 * r12) ** 2 + u[1] ** 2) ** 3) - mum * u[1] /
np.sqrt(((u[0] - pi1 * r12) ** 2 + u[1] ** 2) ** 3)),
# dotu[4] = that
0] # dotu[5] = 0
dt = np.linspace(0.0, 259200.0, 259200.0) # secs to run the simulation
u = odeint(deriv, u0, dt)
x, y, z, x2, y2, z2 = u.T
fig = pylab.figure()
ax = fig.add_subplot(111, projection='3d')
ax.plot(x, y, z, color = 'r')
# adding the moon
phi = np.linspace(0, 2 * np.pi, 100)
theta = np.linspace(0, np.pi, 100)
xm = rm * np.outer(np.cos(phi), np.sin(theta)) + r12 - pi2 * r12
ym = rm * np.outer(np.sin(phi), np.sin(theta))
zm = rm * np.outer(np.ones(np.size(phi)), np.cos(theta))
ax.plot_surface(xm, ym, zm, color = '#696969', linewidth = 0)
ax.auto_scale_xyz([-8000, 385000], [-8000, 385000], [-8000, 385000])
# adding the earth
xe = re * np.outer(np.cos(phi), np.sin(theta)) - pi2 * r12
ye = re * np.outer(np.sin(phi), np.sin(theta))
ze = re * np.outer(np.ones(np.size(phi)), np.cos(theta))
ax.plot_surface(xe, ye, ze, color = '#4169E1', linewidth = 0)
ax.auto_scale_xyz([-8000, 385000], [-8000, 385000], [-8000, 385000])
pylab.show()