Notebook

Bifurcation Diagram¶

In this notebook, we generate the bifurcation diagram.

In [1]:

%matplotlib qt
import numpy as np
from TriangleFunctions import *
import math
import matplotlib.pyplot as plt

Start by looking at triangle code, and determining method to go about it. I will need to come up with a way to keep iterating until get a point that matches, and then graph all the previous points until the last one that matches.

That is, take all the previous points up until the last one that matches, and compute over them to determine their y value on the graph. The given set of points will then be paired with an alpha value. Can geometrically find their y value useing the method in my notebook.

Will also need to start each one on both sides to ensure get all possible points.

Data Setup¶

In [2]:

def return_iso_points(angle1, right_point):
    """anlge1 is top angle, angle2 is right angle.
    The left point is on origin, the right point is defined by you, the top point is found computationally"""
    
    #find left angle
    angle2 = ((180 - angle1) / 2)
    angle3 = angle2
    
    #define left line, slope and intersection.  Use tan to define slope
    m1 = math.tan(math.radians(angle3))
    b1 = 0
    mb1 = [m1, b1]
    
    #define right line
    m2 = -(math.tan(math.radians(angle2)))
    b2 = newB(m2, right_point)
    mb2 = [m2, b2]
    
    # find out when the line defined by left and right side intersect to get your final point.
    pointT = inter(mb2, mb1)
    
    return (0, 0), pointT, right_point

In [3]:

## want to edit this function.  Add into method the while loop, and iteration step value of the angle.
# need to use return point and feed as input into triangle so as to get initial conditions of the triangle
def orbit_points(startPoint, points):
    """
    startPoint: the point to start the process from within the triangle
    points: the set of points that define the outer triangle
    """
    
    # current point to draw from
    current_point = startPoint
    # track visited points for graphing
    visited_points = []
    indx = 0
    
    # define lines of triangle.  Of form [m, b]
    # bottom line is on x axis so don't have to draw it.
    l1 = findFunc(points[0], points[1])
    l2 = findFunc(points[1], points[2])
    
    # define perpindicular functions.  Returns m'
    lp1 = pFunc(l1[0])
    lp2 = pFunc(l2[0])
    
    # iterate until point is in loop.  When point is in set, get its index
    while True:
        
        if(current_point in visited_points):
            # get indx
            indx = visited_points.index(current_point)
            break
          
        
        visited_points.append(current_point)
        
        # for each perpendicular line, define the exact function that passes through the current point to find next step
        lOne = [lp1, newB(lp1, current_point)]
        inter1 = inter(l1, lOne)
        
        lTwo = [lp2, newB(lp2, current_point)]
        inter2 = inter(l2, lTwo)
        
        # this corresponds to the botton line which lies along the x axis so need a different function rule for it.  
        inter3 = [current_point[0], 0]
        
        dist1 = dist(current_point, inter1)
        dist2 = dist(current_point, inter2)
        #dist3 = dist(current_point, inter3)
        ## NEW CODE
        dist3 = current_point[1]
        
        # use largest to pick the longest distance to travel within triangle
        largest = heapq.nlargest(3, [dist1, dist2, dist3])
        cur_largest = largest[0]
        next_indx = largest.index(cur_largest) + 1
        next_largest = largest[next_indx]
        
        """### NEW CODE UPDATE, AFTER BIF 2, OLD IN RED COMMENTS
        if(cur_largest == dist1):
            current_point = inter1
            
        if(cur_largest == dist2):
            current_point = inter2
            
        if(cur_largest == dist3):
            current_point = inter3"""
        
        
        # prevent program from drawing outside of triangle
        
        while(True):
            if((cur_largest == dist1) and (inter1[1] >= points[1][1])):
                cur_largest = largest[next_indx]
            if((cur_largest == dist2) and (inter2[1] >= points[1][1])):
                cur_largest = largest[next_indx]
            if((cur_largest == dist3) and (inter3[1] >= points[1][1])):
                cur_largest = largest[next_indx]
            else:
                break
        
        # if two distances are the same pick a random one to follow
        if(cur_largest == next_largest):
            # make a random choice between the two
            largest_pick = random.choice([cur_largest, next_largest])

            if(largest_pick == dist1):
                current_point = inter1
                #visited_points.append(current_point)

            if(largest_pick == dist2):
                current_point = inter2
                #visited_points.append(current_point)

            if(largest_pick == dist3):
                current_point = inter3
                #visited_points.append(current_point)
         
        # if all different distances and within rules, move to the new point
        if(cur_largest != next_largest):   
            
            if(cur_largest == dist1):
                current_point = inter1
                #visited_points.append(current_point)
            if(cur_largest == dist2):
                current_point = inter2
                #visited_points.append(current_point)
            if(cur_largest == dist3):
                current_point = inter3
                #visited_points.append(current_point)
                
        
    
    return visited_points, indx

In [ ]:

%%time
# test
triangle = return_iso_points(2, [1, 0])
pnts, indx = orbit_points([0, 0], triangle)
#print(pnts), print(indx)

Generate Data¶

In [15]:

np.random.uniform(low=0, high=(1/3))

Out[15]:

0.07427972123566462

In [17]:

# setup for loop to get the set
def gen_data(angle_start, angle_step_size):
    points = []
    indices = []
    x_angle = []
    
    num_steps = round((angle_start/angle_step_size))
    
    curr_angle = angle_start
    
    for i in range(0, num_steps):
        
        x_angle.append(curr_angle)
        
        # ADD RANDOM POINT
        rand_pnt = np.random.uniform(low=0, high=(1/3))
        
        triangle = return_iso_points(curr_angle, [(1/3),0])
        
        # START IN MIDDLE TO MAKE UPPER LEVELS WORK
        #pnts, indx = orbit_points([((1/6)-.001), 0], triangle)
        
        ## NEW RANDOM POINT TEST
        pnts, indx = orbit_points([rand_pnt, 0], triangle)
        
        points.append(pnts)
        indices.append(indx)
        
        curr_angle = curr_angle - angle_step_size 
    
    return points, indices, x_angle

Clean Data¶

Begin by recovering the wanted data points from each set

In [18]:

# get points we care about
def gen_pure_data(points, indices):
    pure_points = []
    for i in range(0, len(indices)):
        curr_indx = indices[i]
        curr_point = points[i]
        
        pp = curr_point[curr_indx:]
        pure_points.append(pp)
    return pure_points

In [19]:

# method get_point_dist will give the edge value of the given point.
def get_y_dist(angle, point):
    #math.sin(x) takes radians, angle comes in degrees
    # math.radians(x) converts from degrees to radians
    
    # construct angle off of alpha, assume angle is alpha
    angle = (180 - angle) / 2
    
    # first get y value of point, as that is what matters:
    h = point[1]
    angle = math.radians(angle)
    dist = h / (math.sin(angle))
    return dist

In [20]:

# for each pure point get distances, will have to simultaneously determine side and add length
# KEEP IN MIND THE INTERVAL LENGTH WE ARE DOING IS 1, SO AXIS IS 3
def pure_dists(pure_points, angle_set):
    scaled_dists = []
    for i in range(0, len(pure_points)):
        curr_pnt = pure_points[i]
        sub_set = []
        for j in range(0, len(curr_pnt)):
            curr_sub_pnt = curr_pnt[j]
            # determine if x or not
            if(curr_sub_pnt[1] == 0):
                pnt = (1/3) + curr_sub_pnt[0]
                sub_set.append(pnt)
                continue
            
            # determine side of axis of point
            # left axis
            # WAS .5
            if(curr_sub_pnt[0] < (1/6)):
                dist = get_y_dist(angle_set[i], curr_sub_pnt)
                pnt = (1/3) - dist
                sub_set.append(pnt)
            # if on right axis    
            if(curr_sub_pnt[0] > (1/6)):
                dist = get_y_dist(angle_set[i], curr_sub_pnt)
                pnt = (2/3) + dist
                sub_set.append(pnt)
                
        scaled_dists.append(sub_set)
    
    return scaled_dists

In [21]:

# reflect points across axis.  Axis is 1.5.
# make into set at end to ensure unique points
def mirror_points(scaled_dists):
    
    axis = .5
    mirrored_dists = []
    
    for i in scaled_dists:
        sub_mirrors = []
        for j in i:
            # first append point, then look at opposite
            sub_mirrors.append(j)
            
            d = abs((axis - j))
            if (j >= axis):
                nd = axis - d
                sub_mirrors.append(nd)
                """if(nd in sub_mirrors):
                    continue
                else:
                    sub_mirrors.append(nd)"""
            if (j < axis):
                nd = axis + d
                sub_mirrors.append(nd)
                """if(nd in sub_mirrors):
                    continue
                else:
                    sub_mirrors.append(nd)"""
        
        mirrored_dists.append(sub_mirrors)     
        
    return mirrored_dists

Data Demonstration¶

Now we make everything

In [22]:

%%time
# enter in angle to start at, and step size
points, indices, x_angle = gen_data(150, 0.05)

CPU times: user 1min 11s, sys: 73.9 ms, total: 1min 12s
Wall time: 1min 12s

In [23]:

pure_points = gen_pure_data(points, indices)

In [24]:

dists = pure_dists(pure_points, x_angle)

In [25]:

final_point_set = mirror_points(dists)

For both possible orientations, mirrored on one another. To get 180 mark at end, append empy point to final set, and x_angle

In [21]:

final_point_set.append([])
x_angle.append(180)

In [28]:

# angles is x, final_point_set is y
for xe, ye in zip(x_angle, final_point_set):
    plt.xticks(np.arange(min(x_angle)-.02, max(x_angle)+1, 30))
    plt.yticks(np.arange(0, 1.001, .5))
    plt.scatter([xe] * len(ye), ye, s=.04, color='black')

Both orientations, no mirroring

In [27]:

for xe, ye in zip(x_angle, dists):
    plt.xticks(np.arange(min(x_angle)-.01, max(x_angle)+1, 90))
    plt.yticks(np.arange(0, 1.001, .5))
    plt.scatter([xe] * len(ye), ye, s=.09, color='black')

Note that copied bars are forming at 0 1, 2, 3.

The thing giving the jumps is when the orientation is flipped.

flipping¶

Mirror image just gives you flipped view, not actual diagram though All kinks worked out, now acts as expected although the same pretty mirror imagery stuff is not going on as in the one above :)

In [14]:

def gen_data2(angle_start, angle_step_size):
    points = []
    indices = []
    points2 = []
    indices2 = []
    x_angle = []
    
    num_steps = round((angle_start/angle_step_size))
    
    curr_angle = angle_start
    
    for i in range(0, num_steps):
        
        x_angle.append(curr_angle)
        
        triangle = return_iso_points(curr_angle, [4,0])
        
        pnts, indx = orbit_points([(1/12), 0], triangle)
        pnts2, indx2 = orbit_points([(3/12), 0], triangle)
        
        points.append(pnts)
        indices.append(indx)
        points2.append(pnts)
        indices2.append(indx)
        
        curr_angle = curr_angle - angle_step_size 
    
    return points, indices, points2, indices2, x_angle

In [15]:

#COPY
def gen_pure_data2(points, indices, points2, indices2):
    pure_points = []
    for i in range(0, len(indices)):
        curr_indx = indices[i]
        curr_point = points[i]
        
        pp = curr_point[curr_indx:]
        
        curr_indx2 = indices2[i]
        curr_point2 = points2[i]
        
        pp2 = curr_point[curr_indx2:]
        
        ppp = pp + pp2
        
        pure_points.append(ppp)
    return pure_points

For just one orientation throughout

In [16]:

%%time
# enter in angle to start at, and step size
points1, indices1, points2, indices2, x_angle = gen_data2(150, 0.1)

---------------------------------------------------------------------------
KeyboardInterrupt                         Traceback (most recent call last)
<timed exec> in <module>

<ipython-input-14-78116d1de604> in gen_data2(angle_start, angle_step_size)
     17 
     18         pnts, indx = orbit_points([(1/12), 0], triangle)
---> 19         pnts2, indx2 = orbit_points([(3/12), 0], triangle)
     20 
     21         points.append(pnts)

<ipython-input-4-b25bddba3670> in orbit_points(startPoint, points)
     25     while True:
     26 
---> 27         if(current_point in visited_points):
     28             # get indx
     29             indx = visited_points.index(current_point)

KeyboardInterrupt:

In [17]:

pure_points2 = gen_pure_data2(points1, indices1, points2, indices2)

---------------------------------------------------------------------------
NameError                                 Traceback (most recent call last)
<ipython-input-17-c59fd22056d7> in <module>
----> 1 pure_points2 = gen_pure_data2(points1, indices1, points2, indices2)

NameError: name 'points1' is not defined

In [ ]:

dists2 = pure_dists(pure_points2, x_angle)

In [29]:

for xe, ye in zip(x_angle, dists2):
    plt.xticks(np.arange(min(x_angle)-.1, max(x_angle)+1, 90))
    plt.yticks(np.arange(0, 1.001, .5))
    plt.scatter([xe] * len(ye), ye, s=.1, color='black')

In [ ]: