In [50]:

    

import numpy as np
import matplotlib.pyplot as plt

g = 9.7864
m1 = 1
m2 = 0.05
l1 = 0.5
l2 = 0.2


    

In [51]:

    

def der_ang(p1, p2, c1, c2):
    return p1 * c1 + p2 * c2

def der_p(ang, c1, c2):
    return c1 * np.sin(ang) + c2

def fixang(ang):
    return np.arctan2(np.sin(ang),np.cos(ang))


    

In [96]:

    

dt = 0.01;
t_max = 100;
count = t_max/dt;
x=0

ang1 = np.zeros((count, 9));
ang2 = np.zeros((count, 9));
p1 = np.zeros((count, 9));
p2 = np.zeros((count, 9));

#constantes de uso comum
c1 = 1 / (m1 * l1 ** 2);
c2 = 1 / (m2 * l2 ** 2);
c3 = -1 * m1 * g * l1;
c4 = -1 * m2 * g * l2;
cp = 1 / (m1 * l1 * l2);

out_p1 = [];
out_ang1 =[];

for j in range(1, 9):
  #condicoes iniciais envolvendo energia
  E_0 = 0.2 + (j-1) * 0.2;
  ang1[1, j] = 0;
  ang2[1, j] = np.pi/2;
  p1[1, j] = np.sqrt(2 * m1 * l1 ** 2 *(E_0 - m2 * l2 * g));
  p2[1, j] = 0;  
  print(E_0);
  
  for i in range(1, int(count-1)):
    cv1 = -np.cos(ang1[i, j] - ang2[i, j]) * cp;
    cv2 = -np.sin(ang1[i, j] - ang2[i, j]) * cp * p1[i, j] * p2[i, j];

    k_ang1_1 = der_ang(p1[i, j], p2[i, j], c1, cv1) * dt;
    k_ang2_1 = der_ang(p2[i, j], p1[i, j], c2, cv1) * dt;
    k_p1_1 = der_p(ang1[i, j], c3, cv2) * dt;
    k_p2_1 = der_p(ang2[i, j], c4, -cv2) * dt;

    xang1 = ang1[i, j] + k_ang1_1;
    xang2 = ang2[i, j] + k_ang2_1;
    xp1 = p1[i, j] + k_p1_1;
    xp2 = p2[i, j] + k_p2_1;

    k_ang1_2 = der_ang(xp1, xp2, c1, cv1) * dt;
    k_ang2_2 = der_ang(xp2, xp1, c2, cv1) * dt;
    k_p1_2 = der_p(xang1, c3, cv2) * dt;
    k_p2_2 = der_p(xang2, c4, -cv2) * dt;

    ang1[i+1, j] = fixang(ang1[i, j] + k_ang1_2);
    ang2[i+1, j] = fixang(ang2[i, j] + k_ang2_2);
    if( np.sign(ang1[i+1, j]) != np.sign(ang1[i, j])):
      if( p2[i, j] > 0):
          x = x + 1;
          out_ang1.append(ang2[i, j])
          out_p1.append(p2[i, j])

    p1[i+1, j] = p1[i, j] + k_p1_2;
    p2[i+1, j] = p2[i, j] + k_p2_2;


    

    




/home/danilo/.local/lib/python3.4/site-packages/ipykernel/__main__.py:6: DeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
/home/danilo/.local/lib/python3.4/site-packages/ipykernel/__main__.py:7: DeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
/home/danilo/.local/lib/python3.4/site-packages/ipykernel/__main__.py:8: DeprecationWarning: using a non-integer number instead of an integer will result in an error in the future
/home/danilo/.local/lib/python3.4/site-packages/ipykernel/__main__.py:9: DeprecationWarning: using a non-integer number instead of an integer will result in an error in the future






    




0.2
0.4
0.6000000000000001
0.8
1.0
1.2
1.4000000000000001
1.6

In [97]:

    

plt.plot(np.array(out_ang1) / np.pi, out_p1, '.k')
plt.show()


    

In [92]:

    

plt.plot(ang2, p2, '.k')
plt.show()


    

In [ ]: