numpy - Python, Four-bar linkage angle-time plot

Question

I'm trying to plot the angle vs. time plot for the output angle of a four-bar linkage (angle fi4 in the image below). This angle is calculated using the solution from the https://scholar.cu.edu.eg/?q=anis/files/week04-mdp206-position_analysis-draft.pdf, page 23.

I'm now trying to plot the fi_4(t) plot and am getting some strange results. The diagram displays the input angle fi2 as blue and output angle fi4 as red. Why is the fi2 fluctuating over time? Shouldn't the fi4 have some sort of sine curve?

Am I missing something here?

Four-bar linkage:

The code:

from __future__ import division
import math
import numpy as np
import matplotlib.pyplot as plt

# Input
#lengths of links (tube testing machine actual lengths)
a = 45.5   #mm
b = 250   #mm
c = 140   #mm
d = 244.244  #mm

# Solution for fi2 being a time function, f(time) = angle
f = 16.7/60    #/s
omega = 2 * np.pi * f   #rad/s

t = np.linspace(0, 50, 100)
y = a * np.sin(omega * t)
x = a * np.cos(omega * t)

fi2 = np.arctan(y/x)

# Solution of the vector loop equation
#https://scholar.cu.edu.eg/?q=anis/files/week04-mdp206-position_analysis-draft.pdf
K1 = d/a
K2 = d/c
K3 = (a**2 - b**2 + c**2 + d**2)/(2*a*c)
A = np.cos(fi2) - K1 - K2*np.cos(fi2) + K3
B = -2*np.sin(fi2)
C = K1 - (K2+1)*np.cos(fi2) + K3
fi4_1 = 2*np.arctan((-B+np.sqrt(B**2 - 4*A*C))/(2*A))
fi4_2 = 2*np.arctan((-B-np.sqrt(B**2 - 4*A*C))/(2*A))

# Plot the fi2 time diagram and fi4 time diagram
plt.plot(t, np.degrees(fi2), color = 'blue')
plt.plot(t, np.degrees(fi4_2), color = 'red')
plt.show()

Diagram:

Answer 1

The linespace(0, 50, 100) is too fast. Replacing it with:

t = np.linspace(0, 5, 100)

Second, all the calculations involving the bare np.arctan() are incorrect. You should use np.arctan2(y, x), which determines the correct quadrant (unlike anything based on y/x where the respective signs of x and y are lost). So:

fi2 = np.arctan2(y, x)  # not: np.arctan(y/x)

...

fi4_1 = 2 * np.arctan2(-B + np.sqrt(B**2 - 4*A*C), 2*A)
fi4_2 = 2 * np.arctan2(-B - np.sqrt(B**2 - 4*A*C), 2*A)

Putting some labels on your plots and showing both solutions for θ_4:

plt.plot(t, np.degrees(fi2) % 360, color = 'k', label=r'$θ_2$')
plt.plot(t, np.degrees(fi4_1) % 360, color = 'b', label=r'$θ_{4_1}$')
plt.plot(t, np.degrees(fi4_2) % 360, color = 'r', label=r'$θ_{4_2}$')
plt.xlabel('t [s]')
plt.ylabel('degrees')
plt.legend()
plt.show()

With these mods, we get: