# Průhyb prutu¶

Vypracoval: Filip Horák, 2017/2018

U prutu dle obrázku sestrojte graf průhybu prutu a nalezněte nejvíce prohnutý průřez.

## Řešení¶

Ze začátku je nutné naimportovat potřebné knihovny.

In [1]:

    %matplotlib inline
    import sympy as sp
    import numpy as np
    import matplotlib.pyplot as plt
    from IPython.core.display import Image
    sp.init_printing()

Vyvoláme obrázek úlohy.

In [2]:

    Image(filename='04-o-pruhyb-1.png', width=500)

Out[2]:

Prut částečně uvolníme a stanovíme deformační podmínku: $$w_B=0$$

In [3]:

    Image(filename='04-o-pruhyb-2.png', width=500)

Out[3]:

Zavedeme potřebné proměnné jako symboly.

In [4]:

    FB,M,L,q=sp.symbols('FB M L q')
    E,J=sp.symbols('E J')
    x1,x2,x3=sp.symbols('x1 x2 x3')
    w1=sp.Function('w1')
    w2=sp.Function('w2')
    w3=sp.Function('w3')
    C1,C2,C3,C4,C5,C6=sp.symbols('C1 C2 C3 C4 C5 C6')

Vyjádříme vnitřní účinky ohybového momentu jednotlivých úseků: $$x_{1} \in \langle 0,8l ;l\rangle,\quad x_{2} \in \langle 0,3l ;0,8l\rangle,\quad x_{3} \in \langle 0 ;0,3l\rangle$$

In [5]:

    Image(filename='04-o-pruhyb-3.png', width=950)

Out[5]:
In [6]:

    M1=FB*(L-x1)-M
    M2=FB*(0.2*L+(0.8*L-x2))-M-q*(0.8*L-x2)**2/2
    M3=FB*(0.2*L+0.5*L+(0.3*L-x3))-M-q*0.5*L*(0.5*L/2+(0.3*L-x3))
    M1,M2,M3

Out[6]:
$$\left ( FB \left(L - x_{1}\right) - M, \quad FB \left(1.0 L - x_{2}\right) - M - \frac{q}{2} \left(0.8 L - x_{2}\right)^{2}, \quad FB \left(1.0 L - x_{3}\right) - 0.5 L q \left(0.55 L - x_{3}\right) - M\right )$$

Napíšeme diferenciální rovnice pro jednotlivé úseky.

In [7]:

    deqn1=w1(x1).diff(x1,x1)+M1/E/J
    deqn2=w2(x2).diff(x2,x2)+M2/E/J
    deqn3=w3(x3).diff(x3,x3)+M3/E/J
    deqn1,deqn2,deqn3

Out[7]:
$$\left ( \frac{d^{2}}{d x_{1}^{2}} \operatorname{w_{1}}{\left (x_{1} \right )} + \frac{1}{E J} \left(FB \left(L - x_{1}\right) - M\right), \quad \frac{d^{2}}{d x_{2}^{2}} \operatorname{w_{2}}{\left (x_{2} \right )} + \frac{1}{E J} \left(FB \left(1.0 L - x_{2}\right) - M - \frac{q}{2} \left(0.8 L - x_{2}\right)^{2}\right), \quad \frac{d^{2}}{d x_{3}^{2}} \operatorname{w_{3}}{\left (x_{3} \right )} + \frac{1}{E J} \left(FB \left(1.0 L - x_{3}\right) - 0.5 L q \left(0.55 L - x_{3}\right) - M\right)\right )$$

Řešení první rovnice.

In [8]:

    dsol1=sp.dsolve(deqn1,w1(x1))
    dsol1

Out[8]:
$$\operatorname{w_{1}}{\left (x_{1} \right )} = C_{1} + C_{2} x_{1} + \frac{FB x_{1}^{3}}{6 E J} + \frac{x_{1}^{2} \left(- FB L + M\right)}{2 E J}$$

Řešení druhé rovnice.

In [9]:

    dsol2=sp.dsolve(deqn2,w2(x2))
    dsol2_modified=dsol2.subs([(C1,C3),(C2,C4)])
    sp.N(dsol2_modified,3)

Out[9]:
$$\operatorname{w_{2}}{\left (x_{2} \right )} = C_{3} + C_{4} x_{2} + \frac{0.0417 q x_{2}^{4}}{E J} + \frac{x_{2}^{3}}{E J} \left(0.167 FB - 0.133 L q\right) + \frac{x_{2}^{2}}{E J} \left(- 0.5 FB L + 0.16 L^{2} q + 0.5 M\right)$$

Řešení třetí rovnice.

In [10]:

    dsol3=sp.dsolve(deqn3,w3(x3))
    dsol3_modified=dsol3.subs([(C1,C5),(C2,C6)])
    sp.N(dsol3_modified,3)

Out[10]:
$$\operatorname{w_{3}}{\left (x_{3} \right )} = C_{5} + C_{6} x_{3} + \frac{x_{3}^{3}}{E J} \left(0.167 FB - 0.0833 L q\right) + \frac{x_{3}^{2}}{E J} \left(- 0.5 FB L + 0.138 L^{2} q + 0.5 M\right)$$

Abychom mohli s rovnicemi dále pracovat musíme nejprve vyřešit průhyby $w_1, w_2,w_3$ pomocí příkazu solve.

In [11]:

    w1_=sp.solve(dsol1,w1(x1))
    w2_=sp.solve(dsol2_modified,w2(x2))
    w3_=sp.solve(dsol3_modified,w3(x3))
    w1,w2,w3=w1_[0],w2_[0],w3_[0]
    w1,sp.N(w2,3),sp.N(w3,3)

Out[11]:
$$\left ( \frac{1}{E J} \left(E J \left(C_{1} + C_{2} x_{1}\right) - \frac{FB L}{2} x_{1}^{2} + \frac{FB x_{1}^{3}}{6} + \frac{M x_{1}^{2}}{2}\right), \quad \frac{0.00167}{E J} \left(600.0 E J \left(C_{3} + C_{4} x_{2}\right) - 300.0 FB L x_{2}^{2} + 100.0 FB x_{2}^{3} + 96.0 L^{2} q x_{2}^{2} - 80.0 L q x_{2}^{3} + 300.0 M x_{2}^{2} + 25.0 q x_{2}^{4}\right), \quad \frac{1}{E J} \left(E J \left(C_{5} + C_{6} x_{3}\right) - 0.5 FB L x_{3}^{2} + 0.167 FB x_{3}^{3} + 0.138 L^{2} q x_{3}^{2} - 0.0833 L q x_{3}^{3} + 0.5 M x_{3}^{2}\right)\right )$$

Vyjádříme derivace $w'_1,w'_2,w'_3$.

In [12]:

    dw1=w1.diff(x1)
    dw2=w2.diff(x2)
    dw3=w3.diff(x3)
    dw1,sp.N(dw2,3),dw3

Out[12]:
$$\left ( \frac{1}{E J} \left(C_{2} E J - FB L x_{1} + \frac{FB x_{1}^{2}}{2} + M x_{1}\right), \quad \frac{0.00167}{E J} \left(600.0 C_{4} E J - 600.0 FB L x_{2} + 300.0 FB x_{2}^{2} + 192.0 L^{2} q x_{2} - 240.0 L q x_{2}^{2} + 600.0 M x_{2} + 100.0 q x_{2}^{3}\right), \quad \frac{1}{E J} \left(C_{6} E J - 1.0 FB L x_{3} + 0.5 FB x_{3}^{2} + 0.275 L^{2} q x_{3} - 0.25 L q x_{3}^{2} + 1.0 M x_{3}\right)\right )$$

Vyjádříme energii napjatosti W.

In [13]:

    dW1=M1**2/2/E/J
    dW2=M2**2/2/E/J
    dW3=M3**2/2/E/J
    W=sp.Integral(dW1,[x1,0.8*L,L])+sp.Integral(dW2,[x2,0.3*L,0.8*L])+sp.Integral(dW3,[x3,0,0.3*L])
    W

Out[13]:
$$\int_{0.8 L}^{L} \frac{\left(FB \left(L - x_{1}\right) - M\right)^{2}}{2 E J}\, dx_{1} + \int_{0.3 L}^{0.8 L} \frac{1}{2 E J} \left(FB \left(1.0 L - x_{2}\right) - M - \frac{q}{2} \left(0.8 L - x_{2}\right)^{2}\right)^{2}\, dx_{2} + \int_{0}^{0.3 L} \frac{1}{2 E J} \left(FB \left(1.0 L - x_{3}\right) - 0.5 L q \left(0.55 L - x_{3}\right) - M\right)^{2}\, dx_{3}$$

Sestavíme deformační podmínku a následně energii napjatosti zderivujeme podle síly $F_B$. $$\frac{\partial W}{\partial F_B} = 0$$

In [14]:

    dp=W.diff(FB)
    dp_=dp.doit()
    dp_

Out[14]:
$$\frac{0.162666666666667 FB L^{3}}{E J} + \frac{0.0501875 L^{4} q}{E J} + \frac{0.161666666666667 L^{3}}{E J} \left(1.0 FB - 1.3 L q\right) + \frac{0.009 L^{3}}{E J} \left(1.0 FB - 0.5 L q\right) + \frac{0.18 L^{2}}{E J} \left(- 2 FB L + M\right) + \frac{0.045 L^{2}}{E J} \left(- 2.0 FB L + 0.775 L^{2} q + 1.0 M\right) + \frac{0.275 L^{2}}{E J} \left(- 2.0 FB L + 1.12 L^{2} q + 1.0 M\right) + \frac{0.2 L}{E J} \left(FB L^{2} - L M\right) + \frac{0.5 L}{E J} \left(1.0 FB L^{2} - 0.32 L^{3} q - 1.0 L M\right) + \frac{0.3 L}{E J} \left(1.0 FB L^{2} - 0.275 L^{3} q - 1.0 L M\right)$$

Sestavíme okrajové podmínky hladkosti a spojitosti: $$\begin{align} w_1(0,8L) & = w_2(0,8L) \\ w'_1(0,8L) & = w'_2(0,8L) \\ w_2(0,3L) & = w_3(0,3L) \\ w'_2(0,3L) & = w'_3(0,3L) \\ w_3(0) & = 0 \\ w'_3(0) & = 0 \end{align}$$

In [15]:

    bc1=w1.subs(x1,0.8*L)-w2.subs(x2,0.8*L)
    bc2=dw1.subs(x1,0.8*L)-dw2.subs(x2,0.8*L)
    bc3=w2.subs(x2,0.3*L)-w3.subs(x3,0.3*L)
    bc4=dw2.subs(x2,0.3*L)-dw3.subs(x3,0.3*L)
    bc5=w3.subs(x3,0)
    bc6=dw3.subs(x3,0)

Vyřešíme soustavu rovnic sestavenou z okrajových podmínek, čímž získáme konstanty $C_1,C_2,C_3,C_4,C_5$,$C_6$ a sílu $F_B$.

In [16]:

    sol_bc=sp.linsolve([sp.expand(bc1),sp.expand(bc2),sp.expand(bc3),\
                        sp.expand(bc4),sp.expand(bc5),sp.expand(bc6),\
                        sp.expand(dp_)],[C1,C2,C3,C4,C5,C6,FB])
    C1_,C2_,C3_,C4_,C5_,C6_,FB_=next(iter(sol_bc))
    sp.N(C1_,3),sp.N(C2_,3),sp.N(C3_,3),\
    sp.N(C4_,3),sp.N(C5_,3),sp.N(C6_,3),sp.N(FB_,3)

Out[16]:
$$\left ( \frac{L^{2}}{E J} \left(- 0.0167 L^{2} q + 1.18 \cdot 10^{-16} M\right), \quad \frac{L}{E J} \left(0.0808 L^{2} q - 1.39 \cdot 10^{-16} M\right), \quad \frac{L^{2}}{E J} \left(0.000337 L^{2} q + 9.71 \cdot 10^{-18} M\right), \quad - \frac{L}{E J} \left(0.0045 L^{2} q + 5.55 \cdot 10^{-17} M\right), \quad 0, \quad 0, \quad 0.192 L q + \frac{1.5 M}{L}\right )$$

Konstanty zpětně dosadíme rovnic průhybů a natočení.

In [17]:

    w1_solved=w1.subs({C1:C1_,C2:C2_,FB:FB_})
    dw1_solved=dw1.subs({C1:C1_,C2:C2_,FB:FB_})
    w2_solved=w2.subs({C3:C3_,C4:C4_,FB:FB_})
    dw2_solved=dw2.subs({C3:C3_,C4:C4_,FB:FB_})
    w3_solved=w3.subs({C5:C5_,C6:C6_,FB:FB_})
    dw3_solved=dw3.subs({C5:C5_,C6:C6_,FB:FB_})

Zavedeme číselné hodnoty veličin. $$\begin{align} L & = 1000\, \mathrm{mm} \\ d & = 30\, \mathrm{mm} \\ M & = 500\, \mathrm{Nmm} \\ E & = 2.1\cdot 10^5\, \mathrm{MPa} \\ J & = \frac{\pi d^4}{64} \approx 39760\, \mathrm{mm^4} \end{align}$$

In [18]:

    L_in=1000
    d=30
    M_in=500
    E_in=2.1e5
    J_in=sp.pi*d**4/64

Vytvoříme souřadnice průhybové čáry pro různá $q$.

In [19]:

    xI=np.linspace(0.8*L_in,L_in,20)
    wI10=[-w1_solved.evalf(subs={x1:ii,q:10,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xI]
    wI20=[-w1_solved.evalf(subs={x1:ii,q:20,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xI]
    wI30=[-w1_solved.evalf(subs={x1:ii,q:30,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xI]

    xII=np.linspace(0.3*L_in,0.8*L_in,20)
    wII10=[-w2_solved.evalf(subs={x2:ii,q:10,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xII]
    wII20=[-w2_solved.evalf(subs={x2:ii,q:20,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xII]
    wII30=[-w2_solved.evalf(subs={x2:ii,q:30,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xII]

    xIII=np.linspace(0,0.3*L_in,20)
    wIII10=[-w3_solved.evalf(subs={x3:ii,q:10,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xIII]
    wIII20=[-w3_solved.evalf(subs={x3:ii,q:20,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xIII]
    wIII30=[-w3_solved.evalf(subs={x3:ii,q:30,L:L_in,E:E_in,J:J_in,M:M_in}) for ii in xIII]

Vykreslíme graf průhybu pro různá $q$.

In [20]:

    x=[]

    for ii in xIII:
        x.append(ii)

    for ii in xII:
        x.append(ii)

    for ii in xI:
        x.append(ii)

    w10=wIII10+wII10+wI10
    w20=wIII20+wII20+wI20
    w30=wIII30+wII30+wI30
    fig,ax=plt.subplots()
    plt.ylabel(r'$w\ \mathrm{[mm]}$').set_fontsize(16)
    plt.xlabel(r'$x\ \mathrm{[mm]}$').set_fontsize(16)
    ax.plot(x,w10,label=r'$q=10\ \mathrm{Nmm}$',lw=2)
    ax.plot(x,w20,label=r'$q=20\ \mathrm{Nmm}$',lw=2)
    ax.plot(x,w30,label=r'$q=20\ \mathrm{Nmm}$',lw=2)
    ax.legend(loc='lower left')
    ax.grid(True)

Spočteme bod největšího průhybu při různých hodnotách $q$ pomocí vztahu: \begin{equation} \frac{\mathrm{d}w_2(x_2)}{\mathrm{d}x_2}=0\Rightarrow x_{w_{max}} \end{equation}

In [21]:

    xI10_sol=sp.solve(dw2_solved.subs({q:10,L:L_in,E:E_in,J:J_in,M:M_in}),x2)
    x10=abs(xI10_sol[1])
    wmax10=w2_solved.subs({q:10,L:L_in,E:E_in,J:J_in,M:M_in,x2:xI10_sol[1]})
    xI20_sol=sp.solve(dw2_solved.subs({q:20,L:L_in,E:E_in,J:J_in,M:M_in}),x2)
    x20=abs(xI20_sol[1])
    wmax20=w2_solved.subs({q:20,L:L_in,E:E_in,J:J_in,M:M_in,x2:xI20_sol[1]})
    xI30_sol=sp.solve(dw2_solved.subs({q:30,L:L_in,E:E_in,J:J_in,M:M_in}),x2)
    x30=abs(xI30_sol[1])
    wmax30=w2_solved.subs({q:30,L:L_in,E:E_in,J:J_in,M:M_in,x2:xI30_sol[1]})
    print ('x10={}; wmax={} mm\nx20={}; wmax={} mm\nx30={};\
     wmax={} mm'.format(sp.N(x10,6),-abs(sp.N(wmax10,4)),\
                        sp.N(x20,6),-abs(sp.N(wmax20,4)),\
                        sp.N(x30,6),-abs(sp.N(wmax30,5))))

    x10=577.567; wmax=-4.979 mm
    x20=577.584; wmax=-9.959 mm
    x30=577.590; wmax=-14.940 mm