Metadata-Version: 2.4
Name: beamfeapy
Version: 0.4.0
Summary: 3D FEM beam solver: Euler-Bernoulli & Timoshenko beams, tapered sections, releases, thermal/prestress loads, modal & buckling analysis, section groups, Excel & HDF5 I/O, Plotly plots
Author: Domenico Gaudioso
License-Expression: MIT
Project-URL: Homepage, https://github.com/DomenicoGaudioso/beamfeapy
Project-URL: Repository, https://github.com/DomenicoGaudioso/beamfeapy
Project-URL: Documentation, https://github.com/DomenicoGaudioso/beamfeapy#readme
Project-URL: Issues, https://github.com/DomenicoGaudioso/beamfeapy/issues
Keywords: FEM,beam,frame,structural,Euler-Bernoulli,Timoshenko,prestress,tapered,finite-element
Classifier: Development Status :: 4 - Beta
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.9
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Programming Language :: Python :: 3.13
Classifier: Operating System :: OS Independent
Classifier: Topic :: Scientific/Engineering
Classifier: Topic :: Scientific/Engineering :: Mathematics
Classifier: Topic :: Scientific/Engineering :: Physics
Classifier: Intended Audience :: Science/Research
Classifier: Intended Audience :: Developers
Classifier: Typing :: Typed
Requires-Python: >=3.9
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy>=1.24
Requires-Dist: scipy>=1.10
Provides-Extra: plot
Requires-Dist: plotly>=6.1.1; extra == "plot"
Requires-Dist: kaleido>=1.0; extra == "plot"
Provides-Extra: matplotlib
Requires-Dist: matplotlib>=3.7; extra == "matplotlib"
Provides-Extra: report
Requires-Dist: matplotlib>=3.7; extra == "report"
Requires-Dist: python-docx>=1.1; extra == "report"
Provides-Extra: excel
Requires-Dist: pandas>=1.5; extra == "excel"
Requires-Dist: openpyxl>=3.1; extra == "excel"
Provides-Extra: hdf5
Requires-Dist: h5py>=3.0; extra == "hdf5"
Provides-Extra: fast
Requires-Dist: pypardiso>=0.4; extra == "fast"
Provides-Extra: all
Requires-Dist: plotly>=6.1.1; extra == "all"
Requires-Dist: kaleido>=1.0; extra == "all"
Requires-Dist: matplotlib>=3.7; extra == "all"
Requires-Dist: python-docx>=1.1; extra == "all"
Requires-Dist: pandas>=1.5; extra == "all"
Requires-Dist: openpyxl>=3.1; extra == "all"
Requires-Dist: h5py>=3.0; extra == "all"
Provides-Extra: dev
Requires-Dist: pytest>=7.0; extra == "dev"
Requires-Dist: plotly>=6.1.1; extra == "dev"
Requires-Dist: kaleido>=1.0; extra == "dev"
Requires-Dist: matplotlib>=3.7; extra == "dev"
Requires-Dist: python-docx>=1.1; extra == "dev"
Requires-Dist: pandas>=1.5; extra == "dev"
Requires-Dist: openpyxl>=3.1; extra == "dev"
Requires-Dist: h5py>=3.0; extra == "dev"
Dynamic: license-file

# beamfeapy

<div align="center">
  <img src="img/Design%206.png" alt="beamfeapy Logo" width="216">
</div>

A Python finite-element solver for the **static, modal and buckling analysis** of **3D frame structures** composed of Euler-Bernoulli and Timoshenko beams — including tapered (non-prismatic) elements, end releases, thermal loads, prestress, nodal settlements, load cases, Excel/HDF5 I/O, export to external solvers, Plotly visualization and a Streamlit web UI.

## 📖 Documentation

Full documentation is published as a **bilingual website** (same topics in both languages):

- 🌐 **Site:** <https://domenicogaudioso.github.io/beamfeapy/>
- 🇬🇧 **English:** <https://domenicogaudioso.github.io/beamfeapy/en.html>
- 🇮🇹 **Italiano:** <https://domenicogaudioso.github.io/beamfeapy/it.html>

The source pages live in [`docs/`](docs/) (Jekyll + just-the-docs).

## Features

- **3D Euler-Bernoulli beam element** (12 DOFs per element: axial + bi-axial bending + torsion)
- **Timoshenko beam element** (shear deformability), `shear=True`
- **Tapered beam element** (non-prismatic / variable section) — exact force-based stiffness, one element with no mesh required
- **Nodal loads** (forces and moments)
- **Concentrated loads in span** (force/moment at internal point ξ∈[0,1])
- **Distributed loads**: uniform, partial, trapezoidal (forces and distributed moments)
- **Thermal loads**: uniform, linear gradient, and **generic nonlinear profiles** along the section height (eigenstress, EN 1991-1-5)
- **Nodal settlements** (imposed displacements/rotations)
- **Prestress**: equivalent-load method for parabolic/straight/eccentric cables, and **3D cable geometry** (anchor + deviation forces from polyline)
- **End releases** (hinges / internal releases via static condensation)
- **Support reactions** and **internal forces** (N, Vy, Vz, T, My, Mz)
- **Load cases**: assign loads to cases; solve single cases or combinations with **multiplicative coefficients** (`solve(cases={"G": 1.35, "Q": 1.5})`)
- **Modal analysis**: masses derived from chosen load cases (distributed + concentrated) with coefficients; natural frequencies, periods, mode shapes, mass participation (validated vs OpenSees)
- **Excel I/O** (Node/Material/Section/Element/Support/Load sheets, WOBridge-style)
- **Results I/O**: save & read results in **Excel** (`.xlsx`, static) or **HDF5** (`.h5`, big-data format for static + modal + buckling in one file)
- **Model export** to external solvers: **OpenSees** (TCL & OpenSeesPy), **SAP2000** (`.s2k`), **MIDAS** (`.mct`), **Robot** (`.str`), **Straus7** (text) — with exact local-axis mapping
- **Springs & cables**: axial springs (bilateral / tension-only / compression-only), ground elastic supports, and **cable elements** (tension-only co-rotational bar with Ernst modulus, and exact elastic catenary) solved with a non-linear Newton scheme
- **Modal & buckling analysis**: natural frequencies/periods/mode shapes and critical load multipliers (λ), with optional per-case section groups — **also with springs and cables present** (cables linearized about their non-linear equilibrium, so the analysis always runs). Cable self-weight is included in the mass; the suspension-bridge example matches the classical beam-under-tension theory to ≈2%
- **Plotly plots**: loads (per load case), internal-force diagrams (European convention), deformed shape, reactions
- **Sparse solver** (COO→CSR assembly + scipy SuperLU, optional MKL Pardiso via `pip install beamfeapy[fast]`) for large models, with vectorized batched assembly, closed-form distributed loads and sparse ARPACK eigensolvers for modal & buckling — within ~2× of OpenSees on large models; see [Performance & Benchmarks](https://domenicogaudioso.github.io/beamfeapy/en-29-performance.html)
- **`solve_many`**: solve many load combinations reusing a single factorization
- **Web UI (Streamlit)** (`app.py`): import a model from Excel, edit it in tables, run analyses and visualize 3D results in the browser — see the [Web UI guide](https://domenicogaudioso.github.io/beamfeapy/en-28-web-ui.html)
- **Command-line interface** (`beamfeapy` / `bfp`): run `solve`, `modal`, `buckling`, `export` and `report` directly from an Excel model, with a colored ASCII logo — the python snake as a simply-supported beam on two tree trunks, body shaded like its bending-moment diagram (`beamfeapy logo`) — see [Command-Line Interface](#command-line-interface-cli)

## Installation

```bash
# From source (development), with all extras (plot + Excel):
pip install -e ".[all]"

# Base only (numpy, scipy):
pip install -e .
```

Optional extras: `plot` (Plotly + kaleido), `excel` (pandas + openpyxl), `all`, `dev`.

Once published on PyPI: `pip install beamfeapy[all]`

**Requirements:** Python ≥ 3.9, numpy ≥ 1.24, scipy ≥ 1.10

## Quick Start

```python
from beamfeapy import Model, Material, Section

m = Model()
m.add_node(1, 0, 0, 0)
m.add_node(2, 4, 0, 0)

mat = Material(E=210e9, nu=0.3, alpha=1.2e-5)
sec = Section(A=1e-2, Iy=2e-5, Iz=3e-5, J=1e-5)
m.add_beam(1, 1, 2, mat, sec)

m.fix(1)                        # fixed support at node 1
m.add_nodal_load(2, Fy=-10000)  # vertical force at tip

res = m.solve()
print(res.displacements(2))     # [ux, uy, uz, rx, ry, rz]
print(res.reactions(1))         # [Fx, Fy, Fz, Mx, My, Mz]
```

## Command-Line Interface (CLI)

Installing the package adds a stylized `beamfeapy` command (alias `bfp`) so you can run analyses straight from an Excel model — no Python code required:

<div align="center">
  <img src="img/snake_anim.gif" alt="beamfeapy animated CLI logo" width="430">
  <img src="img/shot_solve.png" alt="beamfeapy CLI solve" width="430">
</div>

```bash
beamfeapy logo                                 # show the colored app logo
beamfeapy new model.xlsx                        # build a model with a guided wizard
beamfeapy template model.xlsx                 # generate an input template
beamfeapy info model.xlsx                      # model summary (nodes, loads, cases)
beamfeapy solve model.xlsx -o out.xlsx --cases G=1.35 Q=1.5
beamfeapy modal model.xlsx -n 12 --mass-cases G=1.0 Q=0.3 -o modal.h5
beamfeapy buckling model.xlsx --cases P
beamfeapy plot model.xlsx --what deformed --cases G=1.35 Q=1.5 --open   # interactive HTML, opens in browser
beamfeapy plot model.xlsx --what forces --component Mz --png            # static PNG instead
beamfeapy export model.xlsx out.tcl            # → OpenSees / SAP2000 / MIDAS / Robot / Straus7
beamfeapy report model.xlsx -o report.docx     # Word calc report
```

**Build a model** from scratch without leaving the terminal with `beamfeapy new model.xlsx`: a guided wizard asks for nodes, material, section, elements, supports and nodal loads, then writes the Excel model — ready to `solve` or `plot`.

**Visualize** the model and results with `beamfeapy plot`. `--what` selects `model` (geometry + supports), `loads`, `deformed` (deflected shape, auto-scaled), `forces` (internal-force diagrams N/V/T/M, pick one with `--component Mz`), `reactions`, or `all` (writes every figure to a folder). By default it writes an **interactive 3D Plotly HTML** and prints a clickable `file://` link to open in the browser (add `--open` to launch it automatically); pass `--png` for a static Matplotlib image instead:

<div align="center">
  <img src="img/shot_plot.png" alt="beamfeapy CLI plot" width="560">
</div>

<div align="center">
  <img src="img/plot_model_example.png" alt="model plot" width="270">
  <img src="img/plot_deformed_example.png" alt="deformed plot" width="270">
  <img src="img/plot_forces_example.png" alt="internal forces plot" width="270">
</div>

Run `beamfeapy logo` to play the **animated** logo: a little rainbow python that slithers across the terminal **from right to left**, with the `beamfeapy` wordmark shimmering below; it animates for a few seconds in an interactive terminal — use `--loop` to keep it going, `--still` for the static logo. The static logo is the **python snake drawn as a simply-supported beam resting on two cut-trunk supports**. Its body carries the FEM stress-contour gradient — blue at the supports, red at midspan — so it doubles as the bending-moment diagram of a simply-supported beam, while the `beamfeapy` wordmark uses the same blue→red gradient. A compact "mini snake" header precedes the operational commands. `beamfeapy logo --image` shows instead the photo-realistic logo rendered in truecolor half-blocks (from `img/2.png`, embedded compressed via `scripts/_gen_cli_logo.py`). `beamfeapy <command> --help` shows per-command options.

Load combinations are passed as `CASE=COEF` pairs (`G=1.35 Q=1.5` → dict with multiplicative coefficients) or plain names (`G Q`). Output is colorized automatically (24-bit ANSI, with Virtual-Terminal support on Windows); use `--color never` or set `NO_COLOR` to disable (falls back to a stylized box banner), and `--no-banner` to hide the header. The CLI needs no image library at runtime.

## Web UI (Streamlit)

Prefer a graphical interface? A Streamlit app (`app.py`) lets you import a model from Excel, edit it in editable tables, run analyses (static / modal / buckling) and visualize 3D results in the browser:

```bash
pip install beamfeapy streamlit plotly openpyxl
streamlit run app.py
```

See the [Web UI guide](https://domenicogaudioso.github.io/beamfeapy/en-28-web-ui.html) ([italiano](https://domenicogaudioso.github.io/beamfeapy/it-28-web-ui.html)) for the full walkthrough with screenshots.

## API Reference

### Model Construction

| Method | Description |
|--------|-------------|
| `Model()` | Create an empty model |
| `m.add_node(id, x, y, z)` | Add a node with ID and coordinates |
| `m.add_beam(id, ni, nj, mat, sec, ...)` | Add a 3D beam element |
| `m.add_tapered_beam(id, ni, nj, mat, ...)` | Add a tapered (variable-section) beam element |
| `m.add_section(id, A=..., Iy=..., ...)` | Register a section by ID for reuse |

### Materials & Sections

```python
mat = Material(E=210e9, nu=0.3, alpha=1.2e-5)   # steel, SI
mat = Material(E=30e9, nu=0.2, alpha=1.0e-5)     # concrete

sec = Section(A=1e-2, Iy=2e-5, Iz=3e-5, J=1e-5)                # basic (EB)
sec = Section(A=0.18, Iy=5.4e-3, Iz=1.35e-3, J=2e-3,
              Asy=5/6*0.18, Asz=5/6*0.18)                       # Timoshenko
sec = Section(A=..., Iy=..., Iz=..., J=..., h_y=0.6, h_z=0.3)  # with heights for thermal
```

**Section parameters:**

| Parameter | Description |
|-----------|-------------|
| `A` | Cross-sectional area |
| `Iy` | Moment of inertia about local y-axis (bending in x-z plane) |
| `Iz` | Moment of inertia about local z-axis (bending in x-y plane) |
| `J` | Torsional constant |
| `Asy`, `Asz` | Effective shear areas (Timoshenko only) |
| `h_y`, `h_z` | Section heights along local y/z (for thermal gradients) |

### Supports

```python
m.fix(1)                                    # fixed (all 6 DOFs restrained)
m.pin(2)                                    # pin (3 translations restrained)
m.support(3, ux=True, uy=True, rz=True)    # custom: restrain specific DOFs
```

### Loads

| Method | Description |
|--------|-------------|
| `m.add_nodal_load(node, Fx=..., case=...)` | Force/moment at a node (global axes) |
| `m.add_distributed_load(elem, comp, qi, [qj], [a], [b], [frame], [case])` | Distributed load (uniform/partial/trapezoidal) |
| `m.add_concentrated_load(elem, xi, Fx=..., frame=..., case=...)` | Concentrated load at internal point ξ∈[0,1] |
| `m.add_thermal_load(elem, dT_axial=..., dT_grad_y=..., dT_grad_z=..., case=...)` | Thermal load (uniform + linear gradient) |
| `m.add_thermal_profile(elem, profile, axis=..., width=..., case=...)` | Nonlinear thermal profile along section height |
| `m.add_prestress(elem, P, e_i=..., e_j=..., plane=..., sag=..., case=...)` | Prestress via equivalent loads |
| `m.add_cable_prestress(P, points, [elements], case=...)` | Prestress from 3D cable geometry |
| `m.add_settlement(node, dof, value)` | Imposed displacement at a node |

### Solution & Results

```python
res = m.solve()                    # dense solver (default)
res = m.solve(sparse=True)         # sparse solver (large models)
res = m.solve(cases=["G", "Q"])    # load-case combination
res = m.solve(cases="G")           # single load case

res.displacements(node)            # array [ux, uy, uz, rx, ry, rz]
res.displacement(node, "uy")       # single DOF value
res.reactions(node)                 # array [Fx, Fy, Fz, Mx, My, Mz]
res.element_forces[elem_id]        # end forces in local coords (12,)
```

### Post-processing

```python
from beamfeapy import postprocess

di = postprocess.internal_forces(res, elem_id, n=101)
# Returns dict with keys: x, N, Vy, Vz, T, My, Mz

dd = postprocess.element_displacements(res, elem_id, n=51)
# Returns dict with keys: x, u_local (n×6 array)

pts = postprocess.deformed_shape_global(res, elem_id, n=51, scale=100)
# Returns n×3 array of global coordinates (deformed, scaled)
```

### Plotting (requires `pip install beamfeapy[plot]`)

```python
from beamfeapy.plotting import (
    plot_model, plot_loads, plot_diagram,
    plot_deformed, plot_reactions, plot_internal_forces,
)

plot_loads(m, case="G").show()        # structure + loads for load case "G"
plot_diagram(res, "Mz").show()        # Mz diagram (European convention)
plot_deformed(res, scale=200).show()   # deformed shape
plot_reactions(res).show()             # support reactions
plot_internal_forces(res, elem_id=1)   # all 6 diagrams for one element
```

- **`plot_diagram(result, component)`** — component is one of: `N`, `Vy`, `Vz`, `T`, `My`, `Mz`
- **European convention**: negative moment drawn at the extrados (top)

### Excel I/O (requires `pip install beamfeapy[excel]`)

```python
from beamfeapy import Model, read_excel
from beamfeapy.io_excel import write_template

write_template("input.xlsx")             # generate a fillable template
m = read_excel("input.xlsx")             # = Model.from_excel("input.xlsx")
res = m.solve()
res.to_excel("results.xlsx", n_diagram=21)  # displacements, reactions, forces, diagrams
```

## Detailed Feature Guide

### Distributed Loads

```python
# Uniform load over the entire element (local y, -3 kN/m)
m.add_distributed_load(1, "fy", -3000)

# Partial load from x=1m to x=3m (local y, -5 kN/m)
m.add_distributed_load(1, "fy", -5000, a=1.0, b=3.0)

# Trapezoidal load: from 0 to -8 kN/m over the full element
m.add_distributed_load(1, "fy", 0, -8000)

# Trapezoidal partial load from x=2m to x=5m
m.add_distributed_load(1, "fy", -2000, -8000, a=2.0, b=5.0)

# Distributed moments (torsion, bending)
m.add_distributed_load(1, "mx", 500)                    # torsional moment
m.add_distributed_load(1, "mz", 0, -1200, a=2, b=5)     # partial trapezoidal bending moment

# Global-frame distributed load
m.add_distributed_load(1, "fy", -3000, frame="global")
```

**Component** ∈ {`fx`, `fy`, `fz`, `mx`, `my`, `mz`}.  
**Frame**: `"local"` (default) or `"global"`.

### Concentrated Loads in Span

```python
# Force at 35% of element length
m.add_concentrated_load(1, 0.35, Fy=-50000)

# Moment at 70% of element length
m.add_concentrated_load(1, 0.70, Mz=80000)

# Force in global coordinates
m.add_concentrated_load(1, 0.5, Fz=-20000, frame="global")
```

### Thermal Loads

```python
# Uniform temperature increase (+20°C)
m.add_thermal_load(1, dT_axial=20.0)

# Temperature gradient along z (requires section.h_z)
m.add_thermal_load(1, dT_axial=20.0, dT_grad_z=15.0)

# Nonlinear thermal profile (EN 1991-1-5 style)
m.add_thermal_profile(1, lambda s: 15 * (0.5 + s / 0.3), axis="z")
# Or from discrete points:
m.add_thermal_profile(1, [(-0.15, 0), (0.05, 2.5), (0.15, 15)], axis="z", width=0.30)
```

### Settlements

```python
# Vertical settlement of 5 mm at node 2
m.add_settlement(2, "uz", -0.005)
```

### Prestress

```python
# Parabolic cable: P=2 MN, sag=0.35 m, eccentricity zero at ends
m.add_prestress(1, P=2.0e6, sag=0.35)

# Straight eccentric cable: e=0.20 m (constant)
m.add_prestress(2, P=1.5e6, e_i=0.20, e_j=0.20, plane="z")

# Generic eccentricity profile e(ξ):
m.add_prestress(3, P=1.0e6, profile=lambda xi: 0.3 * (1 - (2*xi - 1)**2))

# 3D cable geometry
pts = [(0, 0.3, 0), (5, -0.05, 0), (10, 0.3, 0)]
m.add_cable_prestress(P=3.0e6, points=pts)
```

### End Releases (Hinges)

```python
# Moment-free (hinge) at end j (Mz = 0)
m.add_beam(1, 1, 2, mat, sec, releases_j=["rz"])

# Releases at both ends
m.add_beam(2, 2, 3, mat, sec, releases_i=["rz"], releases_j=["ry", "rz"])
```

Allowed release names: `ux`, `uy`, `uz`, `rx`, `ry`, `rz` (local DOFs).

### Timoshenko Beam

```python
# Shear-deformable beam: provide effective shear areas
sec = Section(A=0.18, Iy=5e-3, Iz=2e-3, J=3e-3,
              Asy=5/6*0.18, Asz=5/6*0.18)
m.add_beam(1, 1, 2, mat, sec, shear=True)
```

### Tapered Beam (Variable Section)

```python
from beamfeapy import VariableSection

# Method 1: continuous function via VariableSection
vs = VariableSection.rectangular(b=0.30, h=lambda xi: 0.70 * (1 - 0.6 * xi))
m.add_tapered_beam(1, 1, 2, mat, vs)

# Method 2: sections at i and j ends (linear interpolation)
m.add_section(1, A=1.2e-2, Iy=4e-5, Iz=6e-5, J=3e-5)
m.add_section(2, A=0.6e-2, Iy=1e-5, Iz=1.5e-5, J=0.8e-5)
m.add_tapered_beam(1, 1, 2, mat, section_i=1, section_j=2)

# Method 3: sections at multiple stations
m.add_section(3, A=1.0e-2, Iy=3e-5, Iz=4e-5, J=2e-5)
m.add_tapered_beam(2, 2, 3, mat, stations={0.0: 1, 0.4: 3, 1.0: 2})
```

### Load Cases & Combinations

```python
m.add_distributed_load(2, "fy", -20e3, case="G")   # permanent loads
m.add_nodal_load(2, Fx=30e3, case="Q")              # variable loads
m.add_thermal_load(1, dT_axial=15, case="T")         # thermal

m.load_cases()                    # → ['G', 'Q', 'T']
res = m.solve(cases=["G", "Q"])   # combine G + Q
res = m.solve(cases="G")          # single case
res = m.solve()                    # all cases combined
```

### Section Orientation (ref_vector and roll)

By default, the local y-axis is determined automatically. For 3D frames (e.g. portal in the x-y plane), use `ref_vector` or `roll`:

```python
# Portal frame in x-y plane: local z = global Z
m.add_beam(2, 2, 3, mat, sec, ref_vector=(0, 0, 1))

# Roll angle (radians) rotates the section about x-local
m.add_beam(1, 1, 2, mat, sec, roll=0.15)
```

**Convention**: `e_y = ref_vector × e_x` (right-hand rule). The `roll` angle rotates the section about the local x-axis after the default orientation.

### Sparse Solver

For models with many DOFs, use the sparse solver for significant speed and memory savings:

```python
res = m.solve(sparse=True)

# many combinations, a single factorization
results = m.solve_many({"ULS": {"G": 1.35, "Q": 1.5}, "SLS": {"G": 1, "Q": 1}})

# sparse eigensolvers (default for large models)
modal = m.modal(n_modes=6, mass_source={"G": 1.0}, sparse=True)
buck  = m.buckling(n_modes=4, cases="G", sparse=True)
```

Install the optional MKL Pardiso backend for an extra speed-up on large systems:

```bash
pip install beamfeapy[fast]   # uses pypardiso automatically when available
```

See [Performance & Benchmarks](https://domenicogaudioso.github.io/beamfeapy/en-29-performance.html) for the full comparison against OpenSees, PyNite and anastruct.

## Conventions

- **Nodal DOFs**: `[ux, uy, uz, rx, ry, rz]` (global system)
- **Local axes**: x from node i to node j; y, z are principal axes of the section.  
  `Iz` → bending in x-y plane; `Iy` → bending in x-z plane.
- **Normal force**: positive in tension.
- **Units**: user's choice (consistent, e.g. SI: N, m, Pa).

## Project Structure

```
FEM/
├── beamfeapy/                  # core library
│   ├── material.py              # Material (E, nu, alpha, G)
│   ├── section.py               # Section (A, Iy, Iz, J, Asy, Asz, h_y, h_z)
│   ├── node.py                  # Node (id, x, y, z)
│   ├── element.py               # BeamElement3D (stiffness, transformation, releases)
│   ├── tapered.py               # VariableSection, TaperedBeamElement3D
│   ├── loads.py                 # NodalLoad, DistributedLoad, ConcentratedLoad, Thermal, Prestress, Settlement
│   ├── integration.py           # Gauss-Legendre quadrature
│   ├── model.py                 # Model: assembly, constraints, load cases, solution, Result
│   ├── postprocess.py           # internal forces, deformed shape
│   ├── io_excel.py              # Excel I/O (read_model, write_results, write_template)
│   └── plotting/                # Plotly visualizations
├── examples/                   # basic examples (save HTML to output/)
├── usage_examples/              # comprehensive usage examples (see below)
├── tests/                      # pytest tests vs analytical solutions
├── validation/                  # validation scripts vs OpenSees / analytical
├── benchmark/                   # benchmarks vs OpenSees / PyNite / CALFEM / pystran / anastruct (static, modal, buckling)
└── pyproject.toml               # packaging (pip install -e ".[all]")
```

## Usage Examples

The `usage_examples/` directory contains self-contained scripts covering every feature:

| # | File | Description |
|---|------|-------------|
| 01 | `01_cantilever_nodal_load.py` | Cantilever beam with tip force — the "hello world" |
| 02 | `02_simply_supported_beam.py` | Simply supported beam with uniform distributed load |
| 03 | `03_distributed_loads.py` | Partial and trapezoidal distributed loads |
| 04 | `04_concentrated_loads_in_span.py` | Concentrated forces and moments at internal points |
| 05 | `05_thermal_loads.py` | Uniform and gradient thermal loads |
| 06 | `06_thermal_profile.py` | Nonlinear thermal profile (EN 1991-1-5, eigenstress) |
| 07 | `07_settlements.py` | Imposed displacements (settlements) |
| 08 | `08_prestress_parabolic.py` | Prestress with parabolic cable (equivalent loads) |
| 09 | `09_prestress_cable_geometry.py` | Prestress from 3D cable polyline geometry |
| 10 | `10_timoshenko.py` | Timoshenko vs Euler-Bernoulli comparison |
| 11 | `11_end_releases.py` | End releases (hinges, pins, partial restraints) |
| 12 | `12_tapered_beam.py` | Tapered (variable-section) beam, single element |
| 13 | `13_tapered_beam_stations.py` | Tapered beam with section IDs and multiple stations |
| 14 | `14_3d_portal_frame.py` | 3D portal frame: columns + beam, multiple load types |
| 15 | `15_load_cases.py` | Load cases and combinations |
| 16 | `16_ref_vector_and_roll.py` | Section orientation with ref_vector and roll angle |
| 17 | `17_support_types.py` | Fix, pin, roller, and custom support conditions |
| 18 | `18_internal_forces.py` | Post-processing internal forces (N, V, M, T diagrams) |
| 19 | `19_plotting.py` | All Plotly visualization functions |
| 20 | `20_excel_io.py` | Excel input/output workflow |
| 21 | `21_sparse_solver.py` | Sparse solver for large models |
| 22 | `22_continuous_beam.py` | Multi-span continuous beam with various load patterns |
| 23 | `23_frame_with_hinges.py` | Frame with internal hinges and distributed loads |
| 24 | `24_prestress_secondary_moments.py` | Prestress in a hyperstatic structure (secondary moments) |

Run any example:

```bash
cd FEM
python usage_examples/01_cantilever_nodal_load.py
```

## Testing

```bash
pip install -e ".[dev]"
python -m pytest tests -q
```

Tests verify results against known analytical solutions: cantilever deflection (`PL³/3EI`), simply supported beam with uniform load (`5qL⁴/384EI`), thermal expansion, Timoshenko shear deformability, prestress, and mode shapes (validated vs OpenSees).

## Validation

See `validation/` for cross-validation against OpenSees and analytical benchmarks (errors at machine precision, ≈1e-10%).

## License

MIT — see `LICENSE`.
