I am a student of mathematics and for my modelling seminar, I need to dip a little bit into engineering in order to model the movements of a suspension bridge under influence of wind. This will be modeled as a one-dimensional vibrating beam and so I am trying to get a grasp on the following partial differential equation:
$$ \rho \frac{\partial^2}{\partial t^2} v(x,t) = -EI \frac{\partial^4}{\partial x^4} v(x,t) -\delta \frac{\partial}{\partial t}v(x,t) +q(x,t)$$
with boundary condition
$$ v(0,t) = v(L,t) = \frac{\partial^2}{\partial x^2} v(0,t) = \frac{\partial^2}{\partial x^2} v(L,t) = 0 $$
Here, $v(x,t)$ is supposed to give the downward deflection at time $t$ of the point at length $x$ in the beam, subject to the load $q(x,t)$. $E$ is Young's modulus, $I$ the second moment of inertia of the beam cross-section, $\rho$ is the mass per unit length and $\delta$ some small dampening coefficient.
I have tried to approximate a solution to this in MATLAB by discreticizing the spatial variable (representing the fourth derivative with finite differences) and introducing a dummy variable for the first derivative in time, thereby transforming the whole thing into an ordinary differential equation in time.
However, my problem is that solutions $v$ seem to exhibit extremely high-frequency oscillations that make the ode45 function take forever to finish and let the resultant graph look almost like a solid block of color even on a time scale of [0 10]. Is a Euler beam supposed to oscillate quite as much? I'm really having trouble even obtaining any sort of useable result.
This is especially the case when trying to insert realistic values for $E$, i.e. at least $10^9$. (I notice though that many numerical simulations across the internet use very small values (<=100) -what might be the justification for that?)
P.S. Apologies if any of this is badly or insufficiently phrased - again, as a mathematician I am a little out of my depth here.
EDIT:
Note that for calculation I have actually divided the equation by $\rho$ on both sides, so we get
$$ \frac{\partial^2}{\partial t^2} v(x,t) = -\frac{EI}{\rho} \frac{\partial^4}{\partial x^4}v(x,t) - \delta' \frac{\partial}{\partial t} v(x,t) + g(x,t) $$
where $g$ is now the downward acceleration (only due to gravity in this case). For actual values, I want to take:
- $E = 2 \cdot 10^{11} N/m^2$
- $I = 20 m^4$. This comes from the formula $ I = \frac{bh^3}{12} $ where $b = 30m$ is the width and $h = 2m$ the height of the beam cross-section.
- $\rho = 2 \cdot 10^5 kg/m$
- $\delta' = 0.01 kg/ms$
- $g(x,t) = 10 m/s^2$ constant, obviously a rough approximation of acceleration due to gravity
- $L = 2000 m$
However, the combined parameter $a := \frac{EI}{\rho}$ which would thus end up as $a = 2 \cdot 10^7 m^4/s^2$ and I am not getting the ODE solver to terminate unless I take something as low as $a = 1 m^4/s^2$ (still takes several minutes).
