Mass-spring system

Mass-spring system#

Introduction#

Wait a few seconds until Python interaction ready! is shown in the top.

After that follow the instructions in the Grasple Exercise.

On this webpage the system of differential equations

\[\begin{split} \frac{d\mathbf{x}}{dt}=A\mathbf{x},\quad \text{where}\quad A=\begin{pmatrix}0&1\\-\frac{k}{M}&-\frac{\gamma }{M}\end{pmatrix} \end{split}\]

is considered, with initial condition

\[\begin{split} \mathbf{x}(0)=\begin{pmatrix}u(0)\\v(0)\end{pmatrix} \end{split}\]

This webpage solves this initial-value problem symbolically, using sympy, and then uses numpy and matplotlib to display the obtained exact solution.

The exact solution is in this case is given by

\[ \mathbf{x}(t) = c_1e^{\lambda_1t}\mathbf{v}_1+c_2e^{\lambda_2t}\mathbf{v}_2. \]

provided \(A\) is (complex) diagonalizable and \(c_1,c_2\) are defined by the initial conditions.

import micropip
await micropip.install("ipympl")
import matplotlib
if not hasattr(matplotlib.rcParams, '_get'):
    matplotlib.rcParams._get = matplotlib.rcParams.get
%matplotlib widget
import sympy as sp
import numpy as np
import matplotlib.pyplot as plt
from IPython.display import HTML, Markdown, clear_output
import matplotlib.animation as animation

def Mass_spring_system(gamma,k,M,u0,v0):

    def solve_ivp(M,gamma,k,u0,v0):
        t = sp.symbols('t')
        u = sp.Function('u')
        v = sp.Function('v')
        ode = [ sp.diff(u(t),t)-v(t) , sp.diff(v(t),t) + k/M*u(t) + gamma/M*v(t) ]
        sol = sp.dsolve(ode)
        constants = sp.solve([sol[0].rhs.subs(t,0)-u0,sol[1].rhs.subs(t,0)-v0],sp.symbols('C1 C2'))
        u = sol[0].rhs.subs(constants)
        v = sol[1].rhs.subs(constants)
        return u,v

    t = sp.symbols('t')

    u,v = solve_ivp(M,gamma,k,u0,v0)

    t_vals = np.linspace(0, 40, 201)
    lam_u = sp.lambdify(t,u,modules=['numpy'])
    u_vals = lam_u(t_vals)
    if isinstance(u_vals, int):
        u_vals = u_vals*np.ones_like(t_vals)
    lam_v = sp.lambdify(t,v,modules=['numpy'])
    v_vals = lam_v(t_vals)
    if isinstance(v_vals, int):
        v_vals = v_vals*np.ones_like(t_vals)

    N = 20
    x_vals = np.linspace(0, 1,N)
    x_vals[0] = 0
    x_vals[1] = 0
    x_vals[2] = 0.05
    for n in np.arange(3,N):
        x_vals[n] = x_vals[n-1]*(-1.0)
    x_vals[-1] = 0
    x_vals[-2] = 0
    x_vals[-3] = 0

    fig = plt.figure()
    ax = fig.add_subplot(111, autoscale_on=False, xlim=(-0.2, 0.2), ylim=(-1.2, 1.2))
    ax.set_aspect('equal')
    ax.grid()
    ax.plot([-0.1,0.1], [1,1], '-', lw=2)
    ax.set_yticks(np.arange(-1.25,1.25,0.25))
    ax.set_yticklabels(-np.arange(-1.25,1.25,0.25))
    ax.set_ylabel('$u(t)$')

    line, = ax.plot([], [], '-', lw=2)
    time_template = '$t$ = %.1fs'
    time_text = ax.text(0, 1.05, '', transform=ax.transAxes)
    mass, = ax.plot([], [], '-', lw=2)
    
    def init():
        line.set_data([], [])
        time_text.set_text('')
        mass.set_data([], [])

        return line, time_text, mass

    def animate(i):
        thisx = x_vals
        thisy = np.linspace(1,-u_vals[i],N).tolist()

        thisx_mass = [-0.05,0.05,0.05,-0.05,-0.05]
        thisy_mass = [0.05,0.05,-0.05,-0.05,0.05]-u_vals[i]

        line.set_data(thisx, thisy)
        time_text.set_text(time_template % (t_vals[i]))
        mass.set_data(thisx_mass, thisy_mass)
        return line, time_text, mass

    ani = animation.FuncAnimation(fig, animate, np.arange(0, len(u_vals)),
                                  interval=25, blit=True, init_func=init,repeat=False)
    
    plt.close(fig);
    
    display(Markdown('# The animation is being rendered, please be patient.'))
    ani_jshtml = ani.to_jshtml()
    ani_html = HTML(ani_jshtml)
    clear_output()
    display(ani_html)
    display(Markdown('# You can slow down the animation (if you want) by using the - button.'))
    
    pass

INPUT CELL#

In the cell below you should put the values of the parameters you see in Grasple.

Please note: Each attempt in Grasple contains different values for the parameters.

After you run the cell the results of the calculations are shown with the set parameter values.

Blow the figure are controls that you can use to manipulate the animation.

''' Set the given values in the following lines '''
gamma = 2.5
k = 5.1
M = 2
u0 = 0.9
v0 = -0.1

''' The next line performs the calculations and shows the results. '''
Mass_spring_system(gamma,k,M,u0,v0)