Linearizing 3 differential equations to make 6 linear equations

Hello I am having trouble trying to find the correct model for this coupled spring system. The scenario is the following we have: Ceiling - Spring - Mass(1) - Spring(2) - Mass(2) - Spring (3) - Mass(3) End.

I came up with the following system of differential equations in the 2nd order to model this problem.

$x_1^{''}=[-k_1x_1-k_2(x_2-x_1)]/m_1$

$x_2^{''}=[k_2(x_2-x_1)-k_3(x_3-x_2)]/m_2$

$x_3^{''}=-k_3(x_5-x_3)/m_3$

Is this the correct model? Afterwards I am trying to linearize these equations into 6 differential equations that I can input in matlab and plot the position of each spring.

So I linearized them and obtained the following:

$y_1^{'}=y_2$

$y_2^{'}=(-k_1y_1-k_2(y_3-y_1)-k_3(y_5-y_3)/m_1$

$y_3^{'}=y_4$

$y_4^{'}=(k_2(y_3-y_1)+k_3(y_5-y_3)/m_2$

$y_5^{'}=y_6$

$y_6^{'}=(-k_3(y_5-y_3))/m_3$

I am not sure if this is correct or not. When I plot them in matlab I get something very strange that doesn't make sense. Does anyone know how the spring mass problems work in Diff EQ? I think your model is off. The force on mass 1 should depend only on the lengths of springs 1 and 2, no? Consider that if you hold mass 2 in place and pull on mass 3, there is no force transferred to mass 1.

Also, in the equation for $x_2''$, I think there is a sign error. The two terms should have opposite signs since they're pulling the mass in opposite directions.

Here's my code in Python. I'm not sure it's right. It's largely ripped off from here: http://wiki.scipy.org/Cookbook/CoupledSpringMassSystem The main/only difference from your problem is that here we include the springs' natural lengths L1, L2, and L3.

from scipy.integrate import odeint
import matplotlib.pyplot as plt

def vectorfield(w,t,p):
y1, y2, y3, y4, y5, y6 = w
m1, m2, m3, k1, k2, k3, L1, L2, L3 = p

# Create y' (the vector)
f = [y2,
( -k1 * (y1 - L1) + k2 * (y3 - y1 - L2)) / m1,
y4,
( -k2 * (y3 - y1 - L2) + k3 * (y5 - y3 - L3) ) / m2,
y6,
( -k3 *(y5 - y3 - L3)  )/m3 ]
return f
# Parameter values
# Masses:
m1 = 10.0
m2 = 15.5
m3 = 12.2
# Spring constants
k1 = 4.0
k2 = 2.0
k3 = 1.0
# Natural lengths
L1 = 0.1
L2 = 2.0
L3 = 3.7

# Initial conditions
y1 = 0.2
y2 = 0.0
y3 = 5.25
y4 = 0.0
y5 = 8.0
y6 = 0.0 # ODE solver parameters
abserr = 1.0e-8
relerr = 1.0e-6
stoptime = 100.0
numpoints = 2500

# Create the time samples for the output of the ODE solver.
# I use a large number of points, only because I want to make
# a plot of the solution that looks nice.
t = [stoptime * float(i) / (numpoints - 1) for i in range(numpoints)]

# Pack up the parameters and initial conditions:
p = [m1, m2, m3, k1, k2, k3, L1, L2, L3]
w0 = [y1, y2, y3, y4, y5, y6]

# Call the ODE solver.
wsol = odeint(vectorfield, w0, t, args=(p,),
atol=abserr, rtol=relerr)

plt.xlabel('t')
plt.grid(True)
plt.hold(True)
lw = 1

plt.plot(t, wsol[:,0], 'b', linewidth=lw)
plt.plot(t, wsol[:,2], 'g', linewidth=lw)
plt.plot(t, wsol[:,4], 'r', linewidth=lw)

plt.legend((r'$y_1$', r'$y_3$', r'$y_5$'))
plt.title('Mass Displacements for the\nCoupled Spring-Mass System')
plt.show()

• I was thinking if you pulled on mass 3 it was a dominoe effect of the whole system. Since you elongate Spring 3 which exerts a force on 2 which exerts a force on 1. – adam Aug 7 '14 at 3:51
• I believe that domino effect would be reflected in the position of $x_2$ if it put force on mass 1. – Matt Adams Aug 7 '14 at 3:57
• I see what you are saying but even when I change the system I still get a messed up plot (See above) – adam Aug 7 '14 at 4:01
• I ran it and my plot seems pretty reasonable. Remove the $k_3$ term from $y_2'$ and fix the sign error in $y_4'$ and I think you'll have it. – Matt Adams Aug 7 '14 at 5:13
• can you post what you have for those equations because I did that and still didnt get it. Also can you post your initial conditions – adam Aug 7 '14 at 5:16

Your second and fourth equations seem to be missing right parentheses. Is this, perhaps, the reason for the strange behavior?

$y_2^{'}=(-k_1y_1-k_2(y_3-y_1)-k_3(y_5-y_3))/m_1$

$y_4^{'}=(-k_2(y_3-y_1)-k_3(y_5-y_3))/m_2$

• I copied and paste the problems from what I had in matlab here and edited them. I probably forgot the ( here. – adam Aug 7 '14 at 3:53