-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathhenon_reversed_stability.py
122 lines (90 loc) · 2.77 KB
/
henon_reversed_stability.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
import numpy as np
from matplotlib import pyplot as plt
plt.style.use('dark_background')
import copy
def henon_map(x, y, a, b):
"""
"""
x_next = 1 - a*x**2 + y
y_next = b*x
return x_next, y_next
def reverse_henon_stability(max_iterations, a, b, x_range, y_range):
"""
Find the reverse Henon map stability.
Args:
max_iterations: int, number of iterations of the henon equation
a: float, constant
b: float, constant
x_range: int, number of pixels on the x-axis
y_range: int, number of pixels on the y-axis
Returns:
"""
xl, xr = -3, 3
yl, yr = 0.8, -0.8
x_list = np.arange(xl, xr, (xr - xl)/x_range)
y_list = np.arange(yl, yr, -(yl - yr)/y_range)
array = np.meshgrid(x_list[:x_range], y_list[:y_range])
x2 = np.zeros(x_range)
y2 = np.zeros(y_range)
iterations_until_divergence = np.meshgrid(x2, y2)
for i in iterations_until_divergence:
for j in i:
j += max_iterations
not_already_diverged = np.ones(np.shape(iterations_until_divergence))
not_already_diverged = not_already_diverged[0] < 1000
for k in range(max_iterations):
x_array_copied = copy.deepcopy(array[0]) # copy array to prevent premature modification of x array
# henon map applied to array
array[0] = array[1] / b
array[1] = a * (array[1] / b)**2 + x_array_copied - 1
r = (array[0]**2 + array[1]**2)**0.5
diverging = (r > 100) & not_already_diverged
not_already_diverged = np.invert(diverging) & not_already_diverged
iterations_until_divergence[0][diverging] = k
return iterations_until_divergence[0]
def assemble_plot():
"""
Assemble a plot of the Henon boundary and attractor for varying constants
Args:
None
Returns:
None (saves .png images)
"""
x, y = 1, 1
a, b = 1.4, 0.3
ls = [[x, y]]
steps = 100000
X = [0 for i in range(steps)]
Y = [0 for i in range(steps)]
X[0], Y[0] = 0, 0 # initial point
for i in range(steps-1):
if abs(X[i] + Y[i])**2 < 1000:
X[i+1] = henon_map(X[i], Y[i], a, b)[0]
Y[i+1] = henon_map(X[i], Y[i], a, b)[1]
for t in range(0, 1000):
a = 1 + t/2000
b = 0.3
steps = 10000
X = [0 for i in range(steps)]
Y = [0 for i in range(steps)]
X[0], Y[0] = 0, 0 # initial point
for i in range(steps-1):
if abs(X[i] + Y[i])**2 < 1000:
X[i+1] = henon_map(X[i], Y[i], a, b)[0]
Y[i+1] = henon_map(X[i], Y[i], a, b)[1]
plt.plot(X, Y, ',', color='white', alpha = 0.5, markersize = 0.2)
plt.imshow(reverse_henon_stability(20, a, b,
x_range=2000,
y_range=1400),
extent=[-3, 3, -0.8, 0.8],
vmin=1,
vmax=14,
aspect=2.2,
cmap='inferno',
alpha=1)
plt.axis('off')
plt.tick_params(labelsize=6)
# plt.show()
plt.savefig('henon_reversed_{0:03d}.png'.format(t), dpi=450, bbox_inches='tight', pad_inches=0)
plt.close()
assemble_plot()