
In this lecture we will estimate the ground velocity and displacement by numerically integrating the acceleration record. We'll use the **trapezoidal rule** to perform the integration. We'll write a function to implement the calculation, but first let's explore the basic idea. Recall that an integration is simply the area beneath a line, defined by a function e.g. $f(t)$, between two limits, $t=a$ and $t=b$.

<Image
	src="/images/courses/multi-degree-of-freedom-dynamics-modal-analysis-and-seismic-response-simulation-in-python/integrating-ground-acceleration/img1.jpg"
	alt="Function f(t) with the shaded area enclosed by the x-axis, f(t) and the vertical boundaries identified by the limits t=a and t=b. | EngineeringSkills.com"
	width="300"
	height="199"
/>
<Caption>
	Fig 1. Function $f(t)$ with the shaded area enclosed by the x-axis, $f(t)$ and
	the vertical boundaries identified by the limits $t=a$ and $t=b$.
</Caption>

The shaded area is given by,

$$
\text{Area} = \int_a^b f(t)\:\mathrm{d}t
\tag{1}
$$

In our case, the function is the acceleration record. To calculate the area, we implement the trapezoidal rule. Instead of a line defined by an analytical function of $t$, consider the same line but now defined by a series of data points.

<Image
	src="/images/courses/multi-degree-of-freedom-dynamics-modal-analysis-and-seismic-response-simulation-in-python/integrating-ground-acceleration/img2.jpg"
	alt="Function defined by a series of data points, f(t_0), f(t_1), f(t_n) at equally spaced intervals t_0, t_1, t_n. | EngineeringSkills.com"
	width="300"
	height="207"
/>
<Caption>
	Fig 2. Function defined by a series of data points, $f(t_0), f(t_1), f(t_n)$ at
	equally spaced intervals $t_0, t_1, t_n$.
</Caption>

Using the trapezoidal rule we can say that the total area under the line defined by the data points is,

$$
\text{Area} = \frac{\Delta_t}{2}\left[f(t_0) + 2f(t_1)+\cdots + 2f(t_{n-1}) + f(t_n) \right]
\tag{2}
$$

where $\Delta_t$ is the time between data points. From this statement of the trapezoidal rule, we can state that the value of velocity at some point $v_n$ is given by,

$$
v_n = v_{n-1} + \frac{\Delta_t}{2}[a_{n-1}+a_n]
\tag{3}
$$

where $a_n$ and $a_{n-1}$ are the current and previous acceleration values. Once we have the velocity values we can repeat the numerical integration to determine the displacement values.

$$
x_n = x_{n-1} + \frac{\Delta_t}{2}[v_{n-1}+v_n]
\tag{4}
$$

Next we'll write a function to implement this calculation for the complete acceleration record. The function will take as arguments, the acceleration record to integrate, the associated time vector, the initial velocity and the initial displacement. The initial velocity and displacement will always be zero in our context, but we pass them as arguments for the sake of maintaining a general reusable function.

```python
def calcVelDisp(A, Time, vel_0, disp_0):
delT = Time[1]-Time[0] #Define the timestep

#Re-assign input arguments to variables used inside the function
vel_im1 = vel_0 #Initial velocity
disp_im1 = disp_0 #Initial displacement

#Initialise containers to hold Velocity and Displacement
Velocity = np.zeros(len(Time))
Velocity[0] = vel_im1
Displacement = np.zeros(len(Time))
Displacement[0] = disp_im1

#Loop through data and integrate (Trapezoidal rule)
for i, m in enumerate(A):
    if i>0:
        A_im1 = A[i-1]
        A_i = A[i]

        #Integrate acceleration values to get velocities
        vel_i = vel_im1 + 0.5*delT*(A_i + A_im1)
        #Integrate velocity values to get displacements
        disp_i = disp_im1 + 0.5*delT*(vel_i+vel_im1)

        #Store data
        Velocity[i] = vel_i
        Displacement[i] = disp_i

        #Update values for next loop iteration
        vel_im1 = vel_i
        disp_im1 = disp_i

return Velocity, Displacement
```

Now we can call our function and assign the arrays returned to `vel` and `disp` variables.

```python
vel, disp = calcVelDisp(acc, time, 0, 0)
```

Now we can visualise the output.

```python
#Plotting
fig, axes = plt.subplots(figsize=(12,14),nrows=4,ncols=1)
axes[0].plot(time,acc,'-',label='Ground acceleration')
axes[0].set_xlim([0,30])
axes[0].set_ylabel('Acceleration ($m/s^2$)')
axes[0].set_title('Ground acceleration')
axes[0].grid()

axes[1].plot(time,vel,'-',label='Ground velocity')
axes[1].set_xlim([0,30])
axes[1].set_ylabel('velocity (m/s)')
axes[1].set_title('Ground velocity')
axes[1].grid()

axes[2].plot(time,disp,'-',label='Ground displacement')
axes[2].set_xlim([0,30])
axes[2].set_ylabel('displacement (m)')
axes[2].set_title('Ground displacement')
axes[2].grid()

axes[3].plot(time,pos_EC,'-',label='Relative structural displacement')
axes[3].set_xlim([0,30])
axes[3].set_xlabel('Time (sec)')
axes[3].set_ylabel('Displacement (m)')
axes[3].set_title('Relative structural response')
axes[3].grid()
plt.show()
```

<Image
	src="/images/courses/multi-degree-of-freedom-dynamics-modal-analysis-and-seismic-response-simulation-in-python/integrating-ground-acceleration/img3.jpg"
	alt="Ground acceleration, velocity, displacement and relative structural displacement (top to bottom). | EngineeringSkills.com"
	width="1000"
	height="1109"
	caption="Fig 3. Ground acceleration, velocity, displacement and relative structural displacement (top to bottom)."
/>

### Plot ground and absolute structural displacement

It's also interesting to plot the ground and absolute structural displacement, overlaid on top of each other.

```python
fig = plt.figure()
axes = fig.add_axes([0.1,0.1,2,2])
axes.plot(time,disp,'-',label='Ground displacement')
axes.plot(time,disp+pos_EC,'-',label='Absolute structural displacement')
axes.set_xlim([0,30])
axes.set_xlabel('Time (sec)')
axes.set_ylabel('Displacement (m)')
axes.set_title('Ground and structure displacement')
axes.grid()
axes.legend(loc='lower right')
plt.show()
```

<Image
	src="/images/courses/multi-degree-of-freedom-dynamics-modal-analysis-and-seismic-response-simulation-in-python/integrating-ground-acceleration/img4.jpg"
	alt="Ground displacement and absolute structural displacement. | EngineeringSkills.com"
	width="1000"
	height="681"
	caption="Fig 4. Ground displacement and absolute structural displacement."
/>


Please log in or enroll to continue

Please log in or enroll to continue