Integral de Duhamel

\begin{equation*} m \frac{d \dot{u}}{d \tau} = F(\tau) \end{equation*}

Reordenando

\begin{equation*} d \dot{u} = \frac{F(\tau)}{m} \, d \tau \end{equation*}

Si $u_{0} = 0$ y $t \geqslant \tau$

\begin{equation*} d u = \frac{d \dot{u}}{\omega} \, \sin \big[ \omega \bigl( t - \tau \bigr) \bigr] \end{equation*}

Integrando

\begin{align*} u = \frac{1}{m \omega} \int_{0}^{t} F(\tau) \sin \big[ \omega \bigl( t - \tau \bigr) \bigr] \, d \tau &= \int_{0}^{t} \frac{F(\tau)}{m \omega} \cos \big( \omega \tau \bigr) \, d \tau \, \sin \big( \omega t \bigr) - \int_{0}^{t} \frac{F(\tau)}{m \omega} \sin \big( \omega \tau \bigr) \, d \tau \, \cos \big( \omega t \bigr) \\ &= A \sin \big( \omega t \bigr) + B \cos \big( \omega t \bigr) \end{align*}

Vibración no amortiguada con carga impulsiva rectangular

Usando los datos de un ejemplo anterior

omega = 15.6604597634
m = 42.0/981
F0 = 2.0
t1 = 1.0

Definiendo el impulso y las constantes

import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate

def rectangular(tau,F0,t1):
    if tau <= t1:
        return F0
    else:
        return 0
    
def dA(tau,omega,F0,t1):
    return rectangular(tau,F0,t1)*np.cos(omega*tau)

def dB(tau,omega,F0,t1):
    return rectangular(tau,F0,t1)*np.sin(omega*tau)

def A(tau,m,omega,F0,t1):
    return integrate.quad(dA, 0, tau, args=(omega,F0, t1))[0]/(m*omega)

def B(tau,m,omega,F0,t1):
    return -integrate.quad(dB, 0, tau, args=(omega,F0, t1))[0]/(m*omega)

t = np.linspace(0,3,3000)
posicion = np.zeros(3000)
for i in range(3000):
    posicion[i] = A(t[i],m,omega,F0,t1)*np.sin(omega*t[i]) + B(t[i],m,omega,F0,t1)*np.cos(omega*t[i]) 
    
plt.figure(figsize=(19,8.5))
plt.plot(t, posicion)
plt.xlabel(r'$t$  seg.')
plt.ylabel(r'$u$  cm.')
plt.grid(True)
plt.show()

Vibración no amortiguada con carga impulsiva senoidal

Usando los datos de un ejemplo anterior

omega = 16.5734727803
m = 25.0/981
F0 = 2.0
t1 = 0.5

Definiendo el impulso y las constantes

import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate

def sinoide(tau,F0,t1):
    if tau <= t1:
        return F0*np.sin(np.pi*tau/t1)
    else:
        return 0
    
def dA(tau,omega,F0,t1):
    return sinoide(tau,F0,t1)*np.cos(omega*tau)

def dB(tau,omega,F0,t1):
    return sinoide(tau,F0,t1)*np.sin(omega*tau)

def A(tau,m,omega,F0,t1):
    return integrate.quad(dA, 0, tau, args=(omega,F0, t1))[0]/(m*omega)

def B(tau,m,omega,F0,t1):
    return -integrate.quad(dB, 0, tau, args=(omega,F0, t1))[0]/(m*omega)

t = np.linspace(0,3,3000)
posicion = np.zeros(3000)
for i in range(3000):
    posicion[i] = A(t[i],m,omega,F0,t1)*np.sin(omega*t[i]) + B(t[i],m,omega,F0,t1)*np.cos(omega*t[i])
    
plt.figure(figsize=(19,8.5))
plt.plot(t, posicion)
plt.xlabel(r'$t$  seg.')
plt.ylabel(r'$u$  cm.')
plt.grid(True)
plt.show()

Vibración no amortiguada con carga impulsiva triangular simétrica

Usando los datos de un ejemplo anterior

omega = 9.90454441153
m = 30.0/981
F0 = 4.5
t1 = 0.3
t2 = 0.6

Definiendo el impulso y las constantes

import numpy as np
import matplotlib.pyplot as plt
from scipy import integrate

def triangular(tau,F0,t1,t2):
    if tau <= t1:
        return (F0*tau)/t1
    elif t1 < tau <= t2:
        return -(F0*(tau - t1))/(t2 - t1) + F0
    else:
        return 0
    
def dA(tau,omega,F0,t1,t2):
    return triangular(tau,F0,t1,t2)*np.cos(omega*tau)

def dB(tau,omega,F0,t1,t2):
    return triangular(tau,F0,t1,t2)*np.sin(omega*tau)

def A(tau,m,omega,F0,t1,t2):
    return integrate.quad(dA, 0, tau, args=(omega,F0,t1,t2))[0]/(m*omega)

def B(tau,m,omega,F0,t1,t2):
    return -integrate.quad(dB, 0, tau, args=(omega,F0,t1,t2))[0]/(m*omega)

t = np.linspace(0,3,3000)
posicion = np.zeros(3000)
for i in range(3000):
    posicion[i] = A(t[i],m,omega,F0,t1,t2)*np.sin(omega*t[i]) + B(t[i],m,omega,F0,t1,t2)*np.cos(omega*t[i]) 
    
plt.figure(figsize=(19,8.5))
plt.plot(t, posicion)
plt.xlabel(r'$t$  seg.')
plt.ylabel(r'$u$  cm.')
plt.grid(True)
plt.show()