Vertical Drop

Simulating a vertical drop using the spring-mass-damper model.

Project Overview

This project implements a computational model of a spring-mass-damper system, a fundamental system in mechanical engineering and physics used to understand oscillatory behavior and damped motion.

The system models the dynamics of a human being falling from a certain initial height and landing on the ground. This system is modeled using the spring-mass-damper system's second order ODE:

$$m \frac{d^2y}{dt^2} = -k \Delta y - c \frac{dy}{dt} - mg$$

We use the spring-mass-damper ODE to simulate three phenomena that occur in our scenario:

  1. The ground as a damping force, "stopping" the person from falling further, denoted by the term $$-c\frac{dy}{dt}$$
  2. The person re-positioning themselves to a standing position after landing as a spring force, "springing" themselves back up, denoted by the term $$-k \Delta y$$
  3. The person falling downward due to the force of gravity, denoted by the term, $$-mg$$

Where the sum of these forces equals the mass of the person times their acceleration. We denote the force of gravity as being negative since it acts downwards, the spring force as negative as it acts opposite of the direction in which the person compresses when they land, and the damping force as negative because it acts opposite of the direction of their velocity when they land. The program has fixed k and c values for each terrain type:

  1. Rocky: k = 80000 and c = 2000
  2. Grassy: k = 30000 and c = 1200
  3. Muddy: k = 8000, c = 3500

Note that when the person is falling downward and still in the air, the only force that is present is the force of gravity and the spring and damping forces are both zero. Upon landing, both the spring and damping forces come into play. This is an example of contact non-linearity, which is common in Impact mechanics, Robotics (foot–ground contact), Vehicle suspension models, and Biomechanics.

The code iterates through each time step and uses Euler's method to calculate the Ground Reaction Force (GRF), the total mechanical energy, and the total dissipated energy. For Euler's method to make sense, we must convert our ODE to a first-order ODE: \[ \dot{y} = v \] \[ \dot{v} = \frac{1}{m}\left(-mg + F_{\text{spring}} + F_{\text{damp}}\right) \] Our Ground Reaction Force (GRF) is simply the sum of the spring and the damping forces: $$F_{\text{GRF}} = F_{\text{spring}} + F_{\text{damp}}$$ We assume that acceleration for a single moment each time step, and calculate it as follows: $$a = \frac{F_{\text{spring}} + F_{\text{damp}} + F_{\text{gravity}}}{m}$$ We then calculate the new velocity and position based off of this acceleration: $$ v_{n+1} = v_n + a\,\Delta t $$ $$ y_{n+1} = y_n + v_{n+1}\,\Delta t $$ The new velocity and position values are used to calculate the Kinetic, Gravitational Potential, and Spring Elastic Potential energies, where the total mechanical energy is the sum of all three. The new velocity value is also used to calculate the amount of energy dissipated from the damping force. To further clarify, the total energy is dissipated through the ground upon landing, which is why our damping force is the only force used to calculate our dissipated energy. Here is the derivation of our dissipated energy: $$F_{\text{damp}} = -cv$$ $$dW_{\text{damp}} = (-cv)dy $$ $$dW_{\text{damp}} = -cv * vdt = -cv^2 dt$$ where the velocity is the derivative of the person's position y with respect to time. $$dE_{\text{dissipated}} = -dW_{\text{damp}} = cv^2 dt$$ where the negative sign denotes the energy leaving the system. The dissipated energy is accumulated over each time step. The total energy is calculated at each time step as the sum of the mechanical and accumulated dissipated energies.

This recursion continues until we reach the total simulation time. We store the GRF values for each time step and the energy values for each time step to result in two charts: GRF over time and Energy over time.

In the future, I would like to try iterating through each time step with different methods such as Runge-Kutta, and see how it compares with the current Euler's method. In particular, I want to see if the Energy vs Time Graph would be eliminated of its little humps, which are currently caused by the use of Euler's Method.