Created
November 5, 2012 14:08
-
-
Save precious/4017348 to your computer and use it in GitHub Desktop.
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
#!/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