import matplotlib
if not hasattr(matplotlib.RcParams, "_get"):
matplotlib.RcParams._get = dict.get
Rotation#
Click –> Live Code on the top right corner of this screen to investigate some influence lines!
import micropip
await micropip.install("ipympl")
import numpy as np
import sympy as sym
import matplotlib.pylab as plt
from ipywidgets import widgets, interact
from matplotlib.patches import Arc,StepPatch
%matplotlib inline
EI, x = sym.symbols('EI, x',real=True)
Av, Bv = sym.symbols('Av, Bv',real=True)
a, L= sym.symbols('a, L ',positive=True,real=True)
C1, C2, C3, C4 = sym.symbols('C1, C2, C3, C4')
q = -Bv * sym.SingularityFunction(x,L,-1) - Av * sym.SingularityFunction(x,0,-1) + 1 * sym.SingularityFunction(x,a,-1)
V = -sym.integrate(q,x)+C1
M = sym.integrate(V,x)+C2
kappa = M / EI
phi = sym.integrate(kappa,x)+C3
w = - sym.integrate(phi,x)+C4
eq1 = sym.Eq(w.subs(x,0),0)
eq2 = sym.Eq(w.subs(x,L),0)
eq3 = sym.Eq(M.subs(x,0),0)
eq4 = sym.Eq(M.subs(x,L),0)
eq5 = sym.Eq(V.subs(x,-1),0)
eq6 = sym.Eq(V.subs(x,L+1),0)
sol = sym.solve([eq1,eq2,eq3,eq4,eq5,eq6],[C1,C2,C3,C4,Av,Bv])
V_sol = V.subs(sol)
V_subs = V_sol.subs([(L,5),(a,2),(EI,5000)])
V_numpy = sym.lambdify([L,a,EI,x],V_sol.rewrite(sym.Piecewise))
M_sol = M.subs(sol)
M_subs = M_sol.subs([(L,5),(a,2),(EI,5000)])
M_numpy = sym.lambdify([L,a,EI,x],M_sol.rewrite(sym.Piecewise))
Av_sol = sol[Av].rewrite(sym.Piecewise)
Av_numpy = sym.lambdify([L,a,EI],Av_sol)
Bv_sol = sol[Bv].rewrite(sym.Piecewise)
Bv_numpy = sym.lambdify([L,a,EI],Bv_sol)
Find rotation line \(\phi\)#
q = -Bv * sym.SingularityFunction(x,L,-1) - Av * sym.SingularityFunction(x,0,-1) + 1 * sym.SingularityFunction(x,a,-1)
V = -sym.integrate(q,x)+C1
M = sym.integrate(V,x)+C2
kappa = M / EI
phi = sym.integrate(kappa,x)+C3
w = -sym.integrate(phi,x)+C4
eq1 = sym.Eq(w.subs(x,0),0)
eq2 = sym.Eq(w.subs(x,L),0)
eq3 = sym.Eq(M.subs(x,0),0)
eq4 = sym.Eq(M.subs(x,L),0)
eq5 = sym.Eq(V.subs(x,-1),0)
eq6 = sym.Eq(V.subs(x,L+1),0)
sol = sym.solve([eq1,eq2,eq3,eq4,eq5,eq6],[C1,C2,C3,C4,Av,Bv])
w_sol = w.subs(sol)
w_subs = w_sol.subs([(L,5),(a,2),(EI,5000)])
w_numpy = sym.lambdify([L,a,EI,x],w_sol.rewrite(sym.Piecewise))
phi_sol = phi.subs(sol)
phi_subs = phi_sol.subs([(L,5),(a,2),(EI,5000)])
phi_numpy = sym.lambdify([L,a,EI,x],phi_sol.rewrite(sym.Piecewise))
display(sym.simplify(phi_subs.rewrite(sym.Piecewise)))
x_plot = np.linspace(0,5,100)
plt.figure()
plt.plot([0,5],[0,0],color='black',linewidth=2)
plt.plot(x_plot,phi_numpy(L=5,a=2,EI=5000,x=x_plot),color="blue")
plt.xlabel('$x$')
plt.ylabel('rotation')
axs = plt.gca()
axs.grid()
axs.invert_yaxis()
title0 = 'rotation for force at $x_F = 2 $'
axs.set_title(title0);
axs.spines['right'].set_color('none')
axs.spines['top'].set_color('none')
axs.spines['bottom'].set_position('zero')
axs.spines['left'].set_position('zero')
plt.show()
\[\begin{split}\displaystyle \begin{cases} \frac{7}{25000} & \text{for}\: x > 5 \\\frac{- x^{2} + 10 x - 18}{25000} & \text{for}\: x > 2 \\\frac{3 x^{2} - 16}{50000} & \text{for}\: x > 0 \\- \frac{1}{3125} & \text{otherwise} \end{cases}\end{split}\]
Find influence line \(\phi\)#
q = -Bv * sym.SingularityFunction(x,L,-1) - Av * sym.SingularityFunction(x,0,-1) + 1 * sym.SingularityFunction(x,a,-2)
V = -sym.integrate(q,x)+C1
M = sym.integrate(V,x)+C2
kappa = M / EI
phi = sym.integrate(kappa,x)+C3
w = -sym.integrate(phi,x)+C4
eq1 = sym.Eq(w.subs(x,0),0)
eq2 = sym.Eq(w.subs(x,L),0)
eq3 = sym.Eq(M.subs(x,0),0)
eq4 = sym.Eq(M.subs(x,L),0)
eq5 = sym.Eq(V.subs(x,-1),0)
eq6 = sym.Eq(V.subs(x,L+1),0)
sol = sym.solve([eq1,eq2,eq3,eq4,eq5,eq6],[C1,C2,C3,C4,Av,Bv])
w_sol = w.subs(sol)
w_subs = w_sol.subs([(L,5),(a,2),(EI,5000)])
w_inf_numpy = sym.lambdify([L,a,EI,x],w_sol.rewrite(sym.Piecewise))
display(sym.simplify(w_subs.rewrite(sym.Piecewise)))
plt.figure()
x_plot = np.linspace(0,5,100)
plt.plot([0,5],[0,0],color='black',linewidth=2)
plt.plot(x_plot,w_inf_numpy(L=5,a=2,EI=5000,x=x_plot),color="blue")
plt.xlabel('$x_F$')
plt.ylabel('Influence factor rotation')
axs = plt.gca()
axs.grid()
axs.invert_yaxis()
title0 = 'Influence line for rotation at $x_\phi = 2 $'
axs.spines['right'].set_color('none')
axs.spines['top'].set_color('none')
axs.spines['bottom'].set_position('zero')
axs.spines['left'].set_position('zero')
axs.set_title(title0);
plt.show()
\[\begin{split}\displaystyle \begin{cases} \frac{13 \left(x - 5\right)}{150000} & \text{for}\: x > 5 \\\frac{- x^{3} - 2 x + 15 \left(x - 2\right)^{2}}{150000} & \text{for}\: x > 2 \\- \frac{x \left(x^{2} + 2\right)}{150000} & \text{for}\: x > 0 \\- \frac{x}{75000} & \text{otherwise} \end{cases}\end{split}\]
Comparison displacement with \(\phi\) and influence line for \(\phi\) at \(x = ...\)#
%matplotlib widget
from matplotlib.patches import Polygon
x = np.linspace(0,5.0001,100)
fig = plt.figure(figsize=(12, 6))
fig.canvas.toolbar_visible = False
gs = fig.add_gridspec(2,2)
axs0 = fig.add_subplot(gs[0, 1])
axs1 = fig.add_subplot(gs[0, 0])
axs2 = fig.add_subplot(gs[1, 0])
axs3 = fig.add_subplot(gs[1, 1])
def func(a,b):
axs0.clear()
axs1.clear() # Clear the existing plot
axs2.clear()
axs3.clear()
axs0.plot([0,5],[0,0],color='black',linewidth=2)
L_i = 0
y = np.array([[-1/6+L_i,1/4], [L_i,0], [1/6+L_i,1/4]])
p = Polygon(y, facecolor = '#00CC99',edgecolor='black')
axs0.plot([-0.25+L_i,0.25+L_i],[0.25,0.25],color='black')
axs0.add_patch(p)
L_i = 5
y = np.array([[-1/6+L_i,1/4], [L_i,0], [1/6+L_i,1/4]])
p = Polygon(y, facecolor = '#00CC99',edgecolor='black')
axs0.plot([-0.25+L_i,0.25+L_i],[0.25,0.25],color='black')
axs0.plot([-0.25+L_i,0.25+L_i],[0.3,0.3],color='black')
axs0.add_patch(p)
axs0.annotate('1 kN',xy = [b+0.1,-1],ha='left',va='top')
axs0.annotate(text='', xy=(b,0), xytext=(b,-1), arrowprops=dict(arrowstyle='simple'))
axs0.set_ylim([-1,1])
axs0.set_aspect('equal')
axs0.invert_yaxis()
axs0.axis('off')
axs3.plot([0,5],[0,0],color='black',linewidth=2)
L_i = 0
y = np.array([[-1/6+L_i,1/4], [L_i,0], [1/6+L_i,1/4]])
p = Polygon(y, facecolor = '#00CC99',edgecolor='black')
axs3.plot([-0.25+L_i,0.25+L_i],[0.25,0.25],color='black')
axs3.add_patch(p)
L_i = 5
y = np.array([[-1/6+L_i,1/4], [L_i,0], [1/6+L_i,1/4]])
p = Polygon(y, facecolor = '#00CC99',edgecolor='black')
axs3.plot([-0.25+L_i,0.25+L_i],[0.25,0.25],color='black')
axs3.plot([-0.25+L_i,0.25+L_i],[0.3,0.3],color='black')
axs3.add_patch(p)
axs3.annotate('1 kNm',xy = [a,-1],ha='center',va='top')
axs3.annotate(text='', xy=(a-0.4,0), xytext=(a+0.4,0), arrowprops=dict(arrowstyle='simple',connectionstyle='arc3,rad=0.8'))
axs3.set_ylim([-1,1])
axs3.set_aspect('equal')
axs3.invert_yaxis()
axs3.axis('off')
axs1.plot([0,5],[0,0],color='black',linewidth=2)
axs1.grid()
axs1.plot(x,w_numpy(5,b,5000,x),color='blue')
axs1.plot([a-0.4,a,a+0.4],[w_numpy(5,b,5000,a)+phi_numpy(5,b,5000,a)*0.4,w_numpy(5,b,5000,a),w_numpy(5,b,5000,a)-phi_numpy(5,b,5000,a)*0.4],linewidth=3)
axs1.plot([a-0.4,a,a+0.4],[w_numpy(5,b,5000,a),w_numpy(5,b,5000,a),w_numpy(5,b,5000,a)],color='black')
if phi_numpy(5,b,5000,a) > 0:
axs1.annotate(text='', xytext=(a+0.3,w_numpy(5,b,5000,a)), xy=(a+0.3,w_numpy(5,b,5000,a)-phi_numpy(5,b,5000,a)*0.3), arrowprops=dict(arrowstyle='simple',connectionstyle="arc3,rad=.5"))
else:
axs1.annotate(text='', xytext=(a+0.3,w_numpy(5,b,5000,a)), xy=(a+0.3,w_numpy(5,b,5000,a)-phi_numpy(5,b,5000,a)*0.3), arrowprops=dict(arrowstyle='simple',connectionstyle="arc3,rad=-.5"))
axs1.annotate('%.5f rad' % phi_numpy(5,b,5000,a),xy = [a+0.4,(w_numpy(5,b,5000,a)+w_numpy(5,b,5000,a)-phi_numpy(5,b,5000,a)*0.3)/2])
#axs1.annotate(text='', xy=(b,0), xytext=(b,-0.0001), arrowprops=dict(arrowstyle='simple'))
axs1.set_ylim([-0.00005,0.0006])
axs1.invert_yaxis()
title0 = 'Rotation for force at at $x_F = '+str(b)+'$, showing rotation at $x_w = '+str(a)+'$'
axs1.set_title(title0)
axs2.plot([0,5],[0,0],color='black',linewidth=2)
axs2.grid()
axs2.plot(x,w_inf_numpy(5,a,5000,x),color='blue')
axs2.plot(b,w_inf_numpy(5,a,5000,b),marker='o')
axs2.annotate('%.5f rad' % np.round(w_inf_numpy(5,a,5000,b),5),xy = [b,w_inf_numpy(5,a,5000,b)-0.00005])
axs2.set_ylim([-0.00045,0.00045])
axs2.invert_yaxis()
title1 = 'Influence line for displacement at $x_w = '+str(a)+'$'
axs2.set_title(title1)
axs1.spines['right'].set_color('none')
axs1.spines['top'].set_color('none')
axs1.spines['bottom'].set_position('zero')
axs1.spines['left'].set_position('zero')
axs2.spines['right'].set_color('none')
axs2.spines['top'].set_color('none')
axs2.spines['bottom'].set_position('zero')
axs2.spines['left'].set_position('zero')
plt.draw()
interact(func, a = widgets.FloatSlider(min=0, max=5, value=2, step=0.1, description="Location rotation x_phi = ... (m)",readout_format='.1f',style= {'description_width': '250px'},layout = {'width': '500px'}),
b = widgets.FloatSlider(min=0, max=5, value=2, step=0.1, description="Location force x_F = ... (m)",readout_format='.1f',style= {'description_width': '250px'},layout = {'width': '500px'}) );