3.1.2 Example: Simply supported beam with constant cross-sectional area
Figure 34 illustrates a simply supported beam subjected to an evenly distributed load \( q \) per length unit and a horizontal force \( P \).
Figure 34: Simply supported beam with an evenly distributed load \( q \) per length unit.
The differential equation for the deflection \( U(X) \) is given by: $$ \begin{equation} \frac{d^2U}{dX^2}+\frac{P}{EI}U=-\frac{q}{2EI}(L^2-X^2), \ P>0 \tag{3.26} \end{equation} $$
with the following boundary conditions: $$ \begin{equation} U(-L)=U(L)=0 \tag{3.27} \end{equation} $$ where \( EI \) denotes the flexural rigidity. The equation (3.26) was linearized by assuming small deflections. Alternatively we may assume \( \frac{dU}{dX}(0)=0 \) due to the symmetry.
For convenience we introduce dimensionless variables: $$ \begin{equation} x=\frac{x}{L},\ u=\frac{P}{qL^2}\cdot U,\ \theta^2=\frac{PL^2}{EI} \tag{3.28} \end{equation} $$
which by substitution in (3.26) and (3.27) yield: $$ \begin{equation} \frac{d^2u}{dx^2}+\theta^2 \cdot u=\theta^2\frac{(1-x^2)}{2},\qquad -1 < x < 1 \tag{3.29} \end{equation} $$
with the corresponding boundary conditions: $$ \begin{equation} u(-1)=0,\ \ u(1)=0 \tag{3.30} \end{equation} $$
Equation (3.29) with boundary conditions (3.30) have the analytical solution: $$ \begin{equation} u(x)=\frac{1}{\theta^2}\cdot \left[ \frac{\cos(\theta x)}{\cos(\theta)}-1 \right]-\frac{(1-x^2)}{2} \tag{3.31} \end{equation} $$
The buckling load for this case is given by \( P_k=\dfrac{\pi EI}{4L^2} \) such that $$ \begin{equation} 0\leq \theta \leq \frac{\pi}{2} \tag{3.32} \end{equation} $$
Numerical solution.
We wish to solve the second order linear ODE in (3.29) by means of shooting methods and choose the alternative approach as outlined in equation (3.19).
Two similar sub-problems are then two be solved; they differ only with the source term and the boundary conditions.
Sub-problem 1 $$ \begin{align} u_0''(x)&=-\theta ^2 \, u_0(x) + \theta ^2 \, \frac{(1-x^2)}{2}, & \qquad u_0(-1)=0,\quad \ u_0'(-1)=0 \tag{3.33}\\ \tag{3.34} \end{align} $$
Sub-problem 2 $$ \begin{align} u_1''(x)&=-\theta^2 \, \cdot u_1(x) & \qquad u_1(-1)=0,\quad \ u_1'(-1)=1 \tag{3.35} \\ \tag{3.36} \end{align} $$
From the superposition principle for linear equation given in equation n (3.19) a complete solution is obtained as a combination of the solutions to the two sub-problems: $$ \begin{equation} u(x)=u_0(x)-\frac{u_0(1)}{u_1(1)}\cdot u_1(x) \tag{3.37} \end{equation} $$
Now, to solve the equation numerically we write (3.33) and (3.35) as a system of first order ODEs:
Sub-problem 1 as a system of first order ODEs
For convenience we introduce the notation \( u_0=y_1 \) and \( u_o'=y_2 \): $$ \begin{align} y_1'=&y_2 \tag{3.38} \\ y_2' =&-\theta^2\cdot \left ( y_1+\frac{(1-x^2)}{2} \right ) \nonumber \end{align} $$
with the corresponding initial conditions $$ \begin{equation} y_1(-1)=0,\qquad y_2(-1)=0 \tag{3.39} \end{equation} $$
Sub-problem 2 as a system of first order ODEs
With \( u_1=y_1 \) and \( u_1'=y_2 \): $$ \begin{align} y_1' =& y_2 \tag{3.40} \\ y_2' =& -\theta^2y_1 \nonumber \end{align} $$ and initial conditions $$ \begin{equation} y_1(-1)=0,\qquad y_2(-1)=1 \tag{3.41} \end{equation} $$
Both sub-problems must be integrated from \( x=-1 \) to \( x=1 \) such that
\( u_0(1) \) and \( u_1(1) \) may be found and the solution to the original
problem in equation (3.29) may be constructed according to (3.37). A
python-code beam_deflect_shoot_constant.py for the problem is listed
below. Here the resulting solution is compared with the analytical
solution too.
# src-ch2/beam_deflect_shoot_constant.py;ODEschemes.py @ git@lrhgit/tkt4140/src/src-ch2/ODEschemes.py;
from ODEschemes import euler, heun, rk4
from numpy import cos, sin
import numpy as np
from matplotlib.pyplot import *
# change some default values to make plots more readable
LNWDT=2; FNT=15; FMLY='helvetica';
rcParams['lines.linewidth'] = LNWDT
rcParams['font.size'] = FNT
rcParams['font.family'] = FMLY
def SubProblem1(y, x):
""" system 1 of Governing differential equation on beam bending with constant cross section
Args:
y(array): an array containg y and its derivatives up to second order. (RHS)
x(array): an array
Returns:
yout(array): dydx RHS of reduced system of 1st order equations
"""
yout = np.zeros_like(y)
yout[:] = [y[1],-theta2*(y[0]+0.5*(1-x**2))]
return yout
def SubProblem2(y, x):
""" system 2 of Governing differential equation on beam bending with constant cross section
Args:
y(array): an array containg y and its derivatives up to second order. (RHS)
x(array): an array
Returns:
yout(array): dydx RHS of reduced system of 1st order equations
"""
yout = np.zeros_like(y)
yout[:] = [y[1], -theta2*y[0]]
return yout
# === main program starts here ===
N = 20 # number of elements
L = 1.0 # half the length of the beam
x = np.linspace(-L, L, N + 1) # allocate space
theta = 1 # PL**2/EI
theta2 = theta**2
solverList = [euler, heun, rk4]
solver = solverList[2]
s = [0, 1] # guessed values
# === shoot ===
y0Sys1 = [0, s[0]] # initial values of u and u'
y0Sys2 = [0, s[1]]
u0 = solver(SubProblem1, y0Sys1,x)
u1 = solver(SubProblem2, y0Sys2,x)
u0 = u0[:,0] # extract deflection from solution data
u1 = u1[:,0]
u = u0 -(u0[-1]/u1[-1])*u1 # interpolate to find correct solution
ua = (1/theta2)*(cos(theta*x)/cos(theta) - 1)-(1 - x**2)/2 # analytical solution
legendList=[] # empty list to append legends as plots are generated
plot(x,u,'y',)
plot(x,ua,'r:')
legendList.append('shooting technique')
legendList.append('analytical solution')
## Add the labels
legend(legendList,loc='best',frameon=False)
ylabel('u')
xlabel('x')
grid(b=True, which='both', color='0.65',linestyle='-')
#savefig('../fig-ch2/beam_deflect_shoot_constant.png', transparent=True)
#savefig('../fig-ch2/beam_deflect_shoot_constant.pdf', transparent=True)
#sh
show()
The output of beam_deflect_shoot_constant.py is shown in Figure 35.
Figure 35: The analytical and numerical solutions for simply supported beam with an evenly distributed load \( q \) per length unit. The results are produced with a shooting method implemented in beam_deflect_shoot_constant.py.
