pyStruct

Python library for discrete structural analysis, featuring 2D Euler-Bernoulli beam frame solvers, nonlinear continuum FEM, and semi-monocoque solvers.

pyStruct is an object-oriented Python library developed for discrete and continuum structural mechanics. Originating from the Aerospace Structures curriculum during my MSc in Aeronautical Engineering at Politecnico di Milano, it provides modular, educational, and high-performance solvers for aerospace structural sizing, frame deformation, and nonlinear elastoplastic solid mechanics.

Architectural Overview

The library is organized into three core structural analysis modules alongside constitutive material modeling tools:

  1. Beam Solver (BEAM): 2D Euler-Bernoulli finite element solver for frame structures subjected to nodal and distributed mechanical loads.
  2. Nonlinear FEM (FEM): 2D quadrilateral continuum finite element formulation with full Newton-Raphson nonlinear iterative solver, capable of modeling elastoplasticity and large deformations.
  3. Semi-Monocoque Analysis (MONOCOQUE): Discrete shear-web and axial-flange solver tailored to thin-walled aerospace structures (wing boxes and fuselage sections).
  4. Constitutive Modeling (material.py): Support for linear elastic and bilinear elasto-plastic constitutive material models with hardening.
pyStruct/
├── BEAM/
│   └── BeamSolver.py      # Euler-Bernoulli beam finite element assembly
├── FEM/
│   ├── fem.py             # Quadrilateral element stiffness & internal force routines
│   ├── mesh.py            # Structured 2D quad mesh generator
│   └── newtonSolver.py    # Damped Newton-Raphson nonlinear solver
├── MONOCOQUE/
│   └── MonocoqueSolver.py # Discrete thin-walled semi-monocoque solver
└── material.py            # Linear and bilinear constitutive laws

1. Frame & Beam Finite Element Solver

The BEAM module models arbitrary planar frame networks using 2-node beam elements. Each node incorporates 3 degrees of freedom: horizontal displacement \(u\), vertical deflection \(w\), and in-plane rotation \(\theta = \frac{dw}{dx}\).

Formulation

For an element of length \(L\), cross-sectional area \(A\), second moment of area \(I\), and Young’s modulus \(E\), the local element stiffness matrix \(\mathbf{k}_e\) relates nodal generalized displacements \(\mathbf{u}_e = [u_1, w_1, \theta_1, u_2, w_2, \theta_2]^T\) to generalized forces \(\mathbf{f}_e\):

\[\mathbf{k}_e = \begin{bmatrix} \frac{EA}{L} & 0 & 0 & -\frac{EA}{L} & 0 & 0 \\ 0 & \frac{12EI}{L^3} & \frac{6EI}{L^2} & 0 & -\frac{12EI}{L^3} & \frac{6EI}{L^2} \\ 0 & \frac{6EI}{L^2} & \frac{4EI}{L} & 0 & -\frac{6EI}{L^2} & \frac{2EI}{L} \\ -\frac{EA}{L} & 0 & 0 & \frac{EA}{L} & 0 & 0 \\ 0 & -\frac{12EI}{L^3} & -\frac{6EI}{L^2} & 0 & \frac{12EI}{L^3} & -\frac{6EI}{L^2} \\ 0 & \frac{6EI}{L^2} & \frac{2EI}{L} & 0 & -\frac{6EI}{L^2} & \frac{4EI}{L} \end{bmatrix}\]

Arbitrary orientations in 2D space are assembled into the global system \(\mathbf{K} \mathbf{U} = \mathbf{F}\) via coordinate transformation matrices \(\mathbf{T}_e\), incorporating self-weight distributed loads \(\mathbf{q} = [0, -\rho A g]^T\) and point loads.

2. 2D Continuum FEM & Nonlinear Newton-Raphson Solver

For continuum 2D plane stress mechanics, pyStruct implements 4-node isoparametric bilinear quadrilateral elements (Q4).

Nonlinear Elasto-Plastic Equilibrium

When material yields or large deflections occur, internal resistance forces \(\mathbf{F}_{\text{int}}(\mathbf{u})\) depend nonlinearly on the displacement vector \(\mathbf{u}\). The equilibrium residual is given by:

\[\mathbf{R}(\mathbf{u}) = \mathbf{F}_{\text{ext}} - \mathbf{F}_{\text{int}}(\mathbf{u}) = \mathbf{0}\]

where internal nodal forces are integrated numerically via Gauss quadrature over the element volume \(\Omega_e\):

\[\mathbf{F}_{\text{int}}^e = \int_{\Omega_e} \mathbf{B}^T \boldsymbol{\sigma}(\boldsymbol{\epsilon}) \, d\Omega\]

To resolve the system under bilinear plasticity with yield stress \(\sigma_y\) and tangent modulus \(E_t\), FEM/newtonSolver.py employs a damped Newton-Raphson scheme:

\[\mathbf{K}_T^{(k)} \Delta \mathbf{u}^{(k)} = \mathbf{F}_{\text{ext}} - \mathbf{F}_{\text{int}}\left(\mathbf{u}^{(k)}\right)\] \[\mathbf{u}^{(k+1)} = \mathbf{u}^{(k)} + \omega \, \Delta \mathbf{u}^{(k)}\]

where \(\mathbf{K}_T = \frac{\partial \mathbf{F}_{\text{int}}}{\partial \mathbf{u}}\) is the tangent stiffness matrix and \(\omega \in (0, 1]\) is a relaxation factor ensuring robust convergence.

3. Aerospace Semi-Monocoque Analysis

For aircraft structures, the MONOCOQUE module implements the classical discrete idealization:

  • Longerons / Stringers: Carry purely normal (axial) direct stresses \(\sigma = N / A\).
  • Skin Panels: Thin panels carrying purely shear flows \(q = \tau t\) (\(\text{N/m}\)).

Given global cross-sectional shear forces \(T_x, T_y\) and twisting moment \(M_z\), the solver determines the shear center, torsional constant, and closed-cell shear flow distribution by enforcing kinematic compatibility and moment equilibrium across multi-cell assemblies.