Semi-Monocoque Theory: corrective solutions


In [1]:
from pint import UnitRegistry
import sympy
import networkx as nx
import numpy as np
import matplotlib.pyplot as plt
import sys
%matplotlib inline
from IPython.display import display

Import Section class, which contains all calculations


In [2]:
from Section import Section

Initialization of sympy symbolic tool and pint for dimension analysis (not really implemented rn as not directly compatible with sympy)


In [3]:
ureg = UnitRegistry()
sympy.init_printing()

Define sympy parameters used for geometric description of sections


In [4]:
A, A0, t, t0, a, b, h, L, E, G = sympy.symbols('A A_0 t t_0 a b h L E G', positive=True)

We also define numerical values for each symbol in order to plot scaled section and perform calculations


In [5]:
values = [(A, 150 * ureg.millimeter**2),(A0, 250  * ureg.millimeter**2),(a, 80 * ureg.millimeter), \
          (b, 20 * ureg.millimeter),(h, 35 * ureg.millimeter),(L, 2000 * ureg.millimeter), \
          (t, 0.8 *ureg.millimeter),(E, 72e3 * ureg.MPa), (G, 27e3 * ureg.MPa)]
datav = [(v[0],v[1].magnitude) for v in values]

Second example: Rectangular section with 6 nodes and 2 loops

Define graph describing the section:

1) stringers are nodes with parameters:

  • x coordinate
  • y coordinate
  • Area

2) panels are oriented edges with parameters:

  • thickness
  • lenght which is automatically calculated

In [6]:
stringers = {1:[(2*a,h),A],
             2:[(a,h),A],
             3:[(sympy.Integer(0),h),A],
             4:[(sympy.Integer(0),sympy.Integer(0)),A],
             5:[(a,sympy.Integer(0)),A],
             6:[(2*a,sympy.Integer(0)),A]}

panels = {(1,2):t,
          (2,3):t,
          (3,4):t,
          (4,5):t,
          (5,6):t,
          (6,1):t,
          (5,2):t}

Define section and perform first calculations


In [7]:
S1 = Section(stringers, panels)

In [8]:
S1.cycles


Out[8]:
$$\left [ \left [ 2, \quad 5, \quad 6, \quad 1\right ], \quad \left [ 2, \quad 3, \quad 4, \quad 5\right ]\right ]$$

Plot of S1 section in original reference frame

Define a dictionary of coordinates used by Networkx to plot section as a Directed graph. Note that arrows are actually just thicker stubs


In [9]:
start_pos={ii: [float(S1.g.node[ii]['ip'][i].subs(datav)) for i in range(2)] for ii in S1.g.nodes() }

In [10]:
plt.figure(figsize=(12,8),dpi=300)
nx.draw(S1.g,with_labels=True, arrows= True, pos=start_pos)
plt.arrow(0,0,20,0)
plt.arrow(0,0,0,20)
#plt.text(0,0, 'CG', fontsize=24)
plt.axis('equal')
plt.title("Section in starting reference Frame",fontsize=16);


Plot of S1 section in inertial reference Frame

Section is plotted wrt center of gravity and rotated (if necessary) so that x and y are principal axes. Center of Gravity and Shear Center are drawn


In [11]:
positions={ii: [float(S1.g.node[ii]['pos'][i].subs(datav)) for i in range(2)] for ii in S1.g.nodes() }

In [12]:
x_ct, y_ct = S1.ct.subs(datav)

plt.figure(figsize=(12,8),dpi=300)
nx.draw(S1.g,with_labels=True, pos=positions)
plt.plot([0],[0],'o',ms=12,label='CG')
plt.plot([x_ct],[y_ct],'^',ms=12, label='SC')
#plt.text(0,0, 'CG', fontsize=24)
#plt.text(x_ct,y_ct, 'SC', fontsize=24)
plt.legend(loc='lower right', shadow=True)
plt.axis('equal')
plt.title("Section in pricipal reference Frame",fontsize=16);


Expression of inertial properties in principal reference frame


In [13]:
sympy.simplify(S1.Ixx), sympy.simplify(S1.Iyy), sympy.simplify(S1.Ixy), sympy.simplify(S1.θ)


Out[13]:
$$\left ( \frac{3 A}{2} h^{2}, \quad 4 A a^{2}, \quad 0, \quad 0\right )$$

In [14]:
S1.symmetry


Out[14]:
[{'edges': [((1, 2), (2, 3)), ((3, 4), (6, 1)), ((4, 5), (5, 6)), (5, 2)],
  'nodes': [(1, 3), (2, 2), (4, 6), (5, 5)]},
 {'edges': [((1, 2), (5, 6)), ((2, 3), (4, 5)), (3, 4), (5, 2), (6, 1)],
  'nodes': [(1, 6), (2, 5), (3, 4)]}]

Compute L matrix: with 6 nodes we expect 3 dofs, two with symmetric load and one with antisymmetric load


In [15]:
S1.compute_L()

In [16]:
S1.L


Out[16]:
$$\left[\begin{matrix}- \frac{1}{2} & 0 & 1\\1 & 0 & 0\\- \frac{1}{2} & 0 & -1\\0 & - \frac{1}{2} & 1\\0 & 1 & 0\\0 & - \frac{1}{2} & -1\end{matrix}\right]$$

Compute H matrix


In [17]:
S1.compute_H()

In [18]:
S1.H.subs(datav)


Out[18]:
$$\left[\begin{matrix}\frac{37}{106} & \frac{8}{53} & - \frac{1}{2}\\- \frac{37}{106} & - \frac{8}{53} & - \frac{1}{2}\\\frac{8}{53} & - \frac{8}{53} & \frac{1}{2}\\\frac{8}{53} & \frac{37}{106} & - \frac{1}{2}\\\frac{16}{53} & - \frac{16}{53} & 0\\- \frac{8}{53} & - \frac{37}{106} & - \frac{1}{2}\\- \frac{8}{53} & \frac{8}{53} & \frac{1}{2}\end{matrix}\right]$$

Compute $\tilde{K}$ and $\tilde{M}$


In [19]:
S1.compute_KM(A,h,t)

In [20]:
S1.Ktilde


Out[20]:
$$\left[\begin{matrix}\frac{3}{2} & 0 & 0\\0 & \frac{3}{2} & 0\\0 & 0 & 4\end{matrix}\right]$$

In [21]:
S1.Mtilde.subs(datav)


Out[21]:
$$\left[\begin{matrix}\frac{296}{371} & \frac{128}{371} & 0\\\frac{128}{371} & \frac{296}{371} & 0\\0 & 0 & \frac{39}{14}\end{matrix}\right]$$

Compute eigenvalues and eigenvectors: results are in the form:

  • eigenvalue
  • multiplicity
  • eigenvector

In [27]:
sol_data = (S1.Ktilde.inv()*(S1.Mtilde.subs(datav))).eigenvects()
sol_data


Out[27]:
$$\left [ \left ( \frac{16}{53}, \quad 1, \quad \left [ \left[\begin{matrix}-1\\1\\0\end{matrix}\right]\right ]\right ), \quad \left ( \frac{39}{56}, \quad 1, \quad \left [ \left[\begin{matrix}0\\0\\1\end{matrix}\right]\right ]\right ), \quad \left ( \frac{16}{21}, \quad 1, \quad \left [ \left[\begin{matrix}1\\1\\0\end{matrix}\right]\right ]\right )\right ]$$

Extract eigenvalues


In [23]:
β2 = [sol[0] for sol in sol_data]
β2


Out[23]:
$$\left [ \frac{16}{53}, \quad \frac{39}{56}, \quad \frac{16}{21}\right ]$$

Extract and normalize eigenvectors


In [26]:
X = [sol[2][0]/sol[2][0].norm() for sol in sol_data]
X


Out[26]:
$$\left [ \left[\begin{matrix}- \frac{\sqrt{2}}{2}\\\frac{\sqrt{2}}{2}\\0\end{matrix}\right], \quad \left[\begin{matrix}0\\0\\1\end{matrix}\right], \quad \left[\begin{matrix}\frac{\sqrt{2}}{2}\\\frac{\sqrt{2}}{2}\\0\end{matrix}\right]\right ]$$

Compute numerical value of $\lambda$


In [25]:
λ = [sympy.N(sympy.sqrt(E*A*h/(G*t)*βi).subs(datav)) for βi in β2]
λ


Out[25]:
$$\left [ 72.6843784311631, \quad 110.397010829098, \quad 115.470053837925\right ]$$

In [ ]: