Skip to content

Instantly share code, notes, and snippets.

@precious
Created November 5, 2012 14:08
Show Gist options
  • Save precious/4017348 to your computer and use it in GitHub Desktop.
Save precious/4017348 to your computer and use it in GitHub Desktop.
#!/usr/bin/python
# -*- coding: utf-8 -*-
from math import *
import sys
# generator for values [a,b] with step 'step'
def frange(a,b,step):
while a <= b:
yield a
a += step
# return function that calls f n times
def get_f_n(n,f):
def f_n(x,alpha):
res = f(x,alpha)
for i in range(n - 1):
res = f(res,alpha)
return res
return f_n
def get_list_average(_list):
return sum(_list)/float(len(_list))
# return indexes of list element with minimum distance
# seems like this function could be written better:(
def get_elems_with_min_distance(_list):
assert(len(_list) > 1)
min_distance = [0,1]
for i in range(len(_list) - 1):
for j in range(i + 1,len(_list)):
if abs(_list[min_distance[0]] - _list[min_distance[1]]) > abs(_list[i] - _list[j]):
min_distance = [i,j]
return min_distance
def detect_doubling(values,prev_preiod):
groups = [[value] for value in sorted(values)] # now we have len(values) groups
next_period = prev_preiod*2
# marge groups until their number > next_period
while len(groups) > next_period:
i,j = get_elems_with_min_distance(map(get_list_average,groups))
# merge groups i, j
groups[i].extend(groups[j])
del groups[j]
# now it's time to detect whether we actually have next_period values or prev_preiod
critical_value = 0.05
averages = sorted(map(get_list_average,groups))
# print sorted(values)
# print sorted(averages)
for i in range(len(averages))[::2]:
if abs(averages[i] - averages[i+1]) < critical_value:
return False
return True
# initialization
f = lambda x, alpha: x**2 + x - 3*alpha
period = 1
step = 0.005
x_values = []
y_values = []
fixed_points = []
num_of_iterations = 49
for alpha in frange(0.1,5.6,step):
values = []
for i in range(1,num_of_iterations + 1):
f_n = get_f_n(i,f)
try:
res = f_n(0,alpha)
for j in range(100):
res = f_n(res,alpha)
values.append(res)
except OverflowError:
print alpha, '[OverflowError]'
values = None
break
if not values:
break
#print alpha
if detect_doubling(values,period):
period *= 2
print 'doubled in interval (%g,%g], now period is %d' % (alpha - step, alpha, period)
fixed_points.append(alpha - step)
# generate data for plot
x_values += [alpha]*num_of_iterations
y_values += values
#print values
# draw plot
# from matplotlib import pyplot
# if '-f' in sys.argv:
# for fp in fixed_points:
# pyplot.plot([fp,fp],[-2.5,1.5],'r--')
# pyplot.plot(x_values,y_values,',')
# pyplot.ylabel('fixed points')
# pyplot.xlabel('alpha')
# pyplot.show()
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment