Note that there are some explanatory texts on larger screens.

plurals
  1. POgetting started with cython by using cdef
    primarykey
    data
    text
    <p>I dont understand how to setup and run code with cython.</p> <p>I added <code>cdef</code>, <code>double</code>, etc to pertinent pieces of my code.</p> <p>setup.py of course the name hello isn't being used. <a href="http://docs.cython.org/src/quickstart/build.html" rel="nofollow">cython doc</a></p> <pre><code>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 ) </code></pre> <p><strong>Edit</strong></p> <p>Redefining the question to make it more of a concrete example.</p> <p>I want to run this system of ODEs through Cython. I am not sure if <code>double</code> is the best choice but the numbers are on the order of magnitudes of 10 to the 10 and are non-integer.</p> <p>So I want to call this in a bigger code where <code>mus</code> is a variable defined in the script calling the module. I have my <code>setup.py</code> file to compile the pyx file. I am not sure with what I need to do so I can call this ode now. Say I name the module 3bodyproblem. I would then call it in the scipt as <code>import 3bodyproblem</code> and then do `3bodyproblem.3bodyproblem(what would this input be)'</p> <p>I have be reading <a href="http://hplgit.github.io/teamods/cyode/cyode-sphinx/main_cyode.html" rel="nofollow">intro to cython for odes</a> but I am not sure how to use their example with mine. Also, if it needs to be in rk format, see the code below the first code.</p> <p><strong>Code 1</strong></p> <pre><code>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 </code></pre> <p><strong>Code 2</strong></p> <pre><code>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 </code></pre> <p>If I just try to compile code 1, cython doesn't know about the variables defined the larger <code>.py</code> file. I want to avoid setting variables in the cython code so it can be used else where with out re compiling. Here are the errors:</p> <pre><code>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 </code></pre> <p><strong>Maybe this will help</strong></p> <p>So I want to use the module in a file like this: (different deriv function but just ignore that the idea is the same)</p> <pre><code>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() </code></pre>
    singulars
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    plurals
    1. This table or related slice is empty.
    1. This table or related slice is empty.
    1. This table or related slice is empty.
 

Querying!

 
Guidance

SQuiL has stopped working due to an internal error.

If you are curious you may find further information in the browser console, which is accessible through the devtools (F12).

Reload