Here are some sympy codes for calculating the element matrices of a beam in finite element method.


In [1]:
from sympy import *

tension_bottom = Symbol('T_b')
K = Symbol('K')
z = Symbol('z')
L = Symbol('L')
E,I = symbols('E I')
tension = (tension_bottom + K*z)
linear_density = Symbol('m')
loading = Symbol('p')

shape_function = Matrix([
1 - 3 * (z/L)**2 + 2 * (z/L)**3,
z * (1 - z/L)**2,
3 * (z/L)**2 - 2 * (z/L)**3,
z**2 / L * (z/L - 1)
])

shape_function_derivation = diff(shape_function, z)

shape_function_derivation_derivation = diff(shape_function, z, z)

M_integral = linear_density * shape_function * shape_function.T

K_integral = shape_function_derivation_derivation * shape_function_derivation_derivation.T *E *I

Kg_integral = shape_function_derivation * shape_function_derivation.T * tension

loading_integral = shape_function * loading

M = integrate(M_integral, (z, 0, L))

K = integrate( K_integral, (z, 0, L))

Kg = integrate(Kg_integral, (z, 0, L))

Loading = integrate(loading_integral, (z, 0, L))

In [2]:
M = expand(M)
K = expand(K)
Kg = expand(Kg)
Loading = expand(Loading)

In [3]:
init_printing()

In [4]:
M


Out[4]:
$$\left[\begin{matrix}\frac{13 L}{35} m & \frac{11 m}{210} L^{2} & \frac{9 L}{70} m & - \frac{13 m}{420} L^{2}\\\frac{11 m}{210} L^{2} & \frac{L^{3} m}{105} & \frac{13 m}{420} L^{2} & - \frac{L^{3} m}{140}\\\frac{9 L}{70} m & \frac{13 m}{420} L^{2} & \frac{13 L}{35} m & - \frac{11 m}{210} L^{2}\\- \frac{13 m}{420} L^{2} & - \frac{L^{3} m}{140} & - \frac{11 m}{210} L^{2} & \frac{L^{3} m}{105}\end{matrix}\right]$$

In [5]:
K


Out[5]:
$$\left[\begin{matrix}\frac{12 E}{L^{3}} I & \frac{6 E}{L^{2}} I & - \frac{12 E}{L^{3}} I & \frac{6 E}{L^{2}} I\\\frac{6 E}{L^{2}} I & \frac{4 E}{L} I & - \frac{6 E}{L^{2}} I & \frac{2 E}{L} I\\- \frac{12 E}{L^{3}} I & - \frac{6 E}{L^{2}} I & \frac{12 E}{L^{3}} I & - \frac{6 E}{L^{2}} I\\\frac{6 E}{L^{2}} I & \frac{2 E}{L} I & - \frac{6 E}{L^{2}} I & \frac{4 E}{L} I\end{matrix}\right]$$

In [6]:
Kg


Out[6]:
$$\left[\begin{matrix}\frac{3 K}{5} + \frac{6 T_{b}}{5 L} & \frac{K L}{10} + \frac{T_{b}}{10} & - \frac{3 K}{5} - \frac{6 T_{b}}{5 L} & \frac{T_{b}}{10}\\\frac{K L}{10} + \frac{T_{b}}{10} & \frac{K L^{2}}{30} + \frac{2 L}{15} T_{b} & - \frac{K L}{10} - \frac{T_{b}}{10} & - \frac{K L^{2}}{60} - \frac{L T_{b}}{30}\\- \frac{3 K}{5} - \frac{6 T_{b}}{5 L} & - \frac{K L}{10} - \frac{T_{b}}{10} & \frac{3 K}{5} + \frac{6 T_{b}}{5 L} & - \frac{T_{b}}{10}\\\frac{T_{b}}{10} & - \frac{K L^{2}}{60} - \frac{L T_{b}}{30} & - \frac{T_{b}}{10} & \frac{K L^{2}}{10} + \frac{2 L}{15} T_{b}\end{matrix}\right]$$

In [7]:
Loading


Out[7]:
$$\left[\begin{matrix}\frac{L p}{2}\\\frac{L^{2} p}{12}\\\frac{L p}{2}\\- \frac{L^{2} p}{12}\end{matrix}\right]$$

In [ ]: