|
1 |
| -# Plot of the Lorenz Attractor based on Edward Lorenz's 1963 "Deterministic |
2 |
| -# Nonperiodic Flow" publication. |
3 |
| -# http://journals.ametsoc.org/doi/abs/10.1175/1520-0469%281963%29020%3C0130%3ADNF%3E2.0.CO%3B2 |
4 |
| -# |
5 |
| -# Note: Because this is a simple non-linear ODE, it would be more easily |
6 |
| -# done using SciPy's ode solver, but this approach depends only |
7 |
| -# upon NumPy. |
| 1 | +''' |
| 2 | +Plot of the Lorenz Attractor based on Edward Lorenz's 1963 "Deterministic |
| 3 | +Nonperiodic Flow" publication. |
| 4 | +http://journals.ametsoc.org/doi/abs/10.1175/1520-0469%281963%29020%3C0130%3ADNF%3E2.0.CO%3B2 |
| 5 | +
|
| 6 | +Note: Because this is a simple non-linear ODE, it would be more easily |
| 7 | + done using SciPy's ode solver, but this approach depends only |
| 8 | + upon NumPy. |
| 9 | +''' |
8 | 10 |
|
9 | 11 | import numpy as np
|
10 | 12 | import matplotlib.pyplot as plt
|
11 | 13 | from mpl_toolkits.mplot3d import Axes3D
|
12 | 14 |
|
13 | 15 |
|
14 | 16 | def lorenz(x, y, z, s=10, r=28, b=2.667):
|
| 17 | + ''' |
| 18 | + Given: |
| 19 | + x, y, z: a point of interest in three dimensional space |
| 20 | + s, r, b: parameters defining the lorenz attractor |
| 21 | + Returns: |
| 22 | + x_dot, y_dot, z_dot: values of the lorenz attractor's partial |
| 23 | + derivatives at the point x, y, z |
| 24 | + ''' |
15 | 25 | x_dot = s*(y - x)
|
16 | 26 | y_dot = r*x - y - x*z
|
17 | 27 | z_dot = x*y - b*z
|
18 | 28 | return x_dot, y_dot, z_dot
|
19 | 29 |
|
20 | 30 |
|
21 | 31 | dt = 0.01
|
22 |
| -stepCnt = 10000 |
| 32 | +num_steps = 10000 |
23 | 33 |
|
24 | 34 | # Need one more for the initial values
|
25 |
| -xs = np.empty((stepCnt + 1,)) |
26 |
| -ys = np.empty((stepCnt + 1,)) |
27 |
| -zs = np.empty((stepCnt + 1,)) |
| 35 | +xs = np.empty((num_steps + 1,)) |
| 36 | +ys = np.empty((num_steps + 1,)) |
| 37 | +zs = np.empty((num_steps + 1,)) |
28 | 38 |
|
29 |
| -# Setting initial values |
| 39 | +# Set initial values |
30 | 40 | xs[0], ys[0], zs[0] = (0., 1., 1.05)
|
31 | 41 |
|
32 |
| -# Stepping through "time". |
33 |
| -for i in range(stepCnt): |
34 |
| - # Derivatives of the X, Y, Z state |
| 42 | +# Step through "time", calculating the partial derivatives at the current point |
| 43 | +# and using them to estimate the next point<
7585
/div> |
| 44 | +for i in range(num_steps): |
35 | 45 | x_dot, y_dot, z_dot = lorenz(xs[i], ys[i], zs[i])
|
36 | 46 | xs[i + 1] = xs[i] + (x_dot * dt)
|
37 | 47 | ys[i + 1] = ys[i] + (y_dot * dt)
|
38 | 48 | zs[i + 1] = zs[i] + (z_dot * dt)
|
39 | 49 |
|
| 50 | + |
| 51 | +# Plot |
40 | 52 | fig = plt.figure()
|
41 | 53 | ax = fig.gca(projection='3d')
|
42 | 54 |
|
|
0 commit comments