Metadata-Version: 2.4
Name: wrinklefe
Version: 1.4.0
Summary: Finite element analysis of wrinkled composite laminates with advanced failure theories
Author: Rani F. Elhajjar
License: MIT
Project-URL: Documentation, https://wrinklefe.readthedocs.io
Project-URL: Repository, https://github.com/elhajjar1/wrinkleFE
Project-URL: Issues, https://github.com/elhajjar1/wrinkleFE/issues
Project-URL: Changelog, https://github.com/elhajjar1/wrinkleFE/blob/main/CHANGELOG.md
Keywords: composites,finite-element,wrinkle,failure,laminate
Classifier: Development Status :: 3 - Alpha
Classifier: Intended Audience :: Science/Research
Classifier: Topic :: Scientific/Engineering
Classifier: Programming Language :: Python :: 3
Classifier: Programming Language :: Python :: 3.10
Classifier: Programming Language :: Python :: 3.11
Classifier: Programming Language :: Python :: 3.12
Classifier: Programming Language :: Python :: 3.13
Requires-Python: >=3.10
Description-Content-Type: text/markdown
License-File: LICENSE
Requires-Dist: numpy>=1.24
Requires-Dist: scipy>=1.10
Requires-Dist: matplotlib>=3.7
Provides-Extra: streamlit
Requires-Dist: streamlit>=1.57; extra == "streamlit"
Requires-Dist: plotly>=6.7; extra == "streamlit"
Provides-Extra: plotly
Requires-Dist: plotly>=6.7; extra == "plotly"
Provides-Extra: dev
Requires-Dist: pytest>=7.0; extra == "dev"
Requires-Dist: pytest-cov>=4.0; extra == "dev"
Requires-Dist: pytest-benchmark>=4.0; extra == "dev"
Requires-Dist: ruff>=0.1; extra == "dev"
Requires-Dist: mypy>=1.5; extra == "dev"
Requires-Dist: scipy-stubs>=1.14; extra == "dev"
Provides-Extra: docs
Requires-Dist: sphinx>=8.0; extra == "docs"
Requires-Dist: myst-parser>=4.0; extra == "docs"
Requires-Dist: furo>=2024.8; extra == "docs"
Requires-Dist: plotly>=6.7; extra == "docs"
Provides-Extra: vtk
Requires-Dist: pyvista>=0.42; extra == "vtk"
Provides-Extra: all
Requires-Dist: wrinklefe[dev,docs,streamlit,vtk]; extra == "all"
Dynamic: license-file

# WrinkleFE

[![Open in Streamlit](https://static.streamlit.io/badges/streamlit_badge_black_white.svg)](https://wrinklefe.streamlit.app/)
[![License: MIT](https://img.shields.io/badge/License-MIT-yellow.svg)](LICENSE)
[![Software DOI](https://zenodo.org/badge/DOI/10.5281/zenodo.23105908.svg)](https://doi.org/10.5281/zenodo.23105908)
[![Article DOI](https://img.shields.io/badge/article%20DOI-10.1038%2Fs41598--025--06693--4-blue.svg)](https://doi.org/10.1038/s41598-025-06693-4)
[![Python](https://img.shields.io/badge/python-3.10%2B-blue.svg)](pyproject.toml)
[![GitHub stars](https://img.shields.io/github/stars/elhajjar1/wrinkleFE?style=social)](https://github.com/elhajjar1/wrinkleFE)

An open-source Python finite element package for predicting strength and stiffness knockdown in composite laminates containing fibre waviness defects.

> ⭐ **Found WrinkleFE useful?** Please [star the repository](https://github.com/elhajjar1/wrinkleFE) and [cite it](#citation) — it's a free academic project, and stars and citations are what keep it supported.

## Try it in your browser

The fastest way to use WrinkleFE is the hosted Streamlit app — no install required:

### [Launch the WrinkleFE Streamlit app](https://wrinklefe.streamlit.app/)

Pick a material, set the wrinkle amplitude / wavelength / morphology, and the
app returns the analytical knockdown, plots, and (optionally) a full FE solve.
Public link, no account needed.

## Features

- **Compression model:** CLT-weighted Budiansky-Fleck kink-band with layup-dependent confinement
- **Tension model:** Three-mechanism (fibre cos^2 theta, Hashin matrix, curved-beam sigma_33 delamination) with thick-ply in-situ correction
- **3D finite element:** Structured hexahedral mesh with LaRC04/05 failure criteria
- **Six morphologies:** Stack, convex, concave, uniform, graded (with configurable decay floor) and tool_flat (one surface flat on the tool, for one-sided surface troughs) — compared side by side in one click from the app (**Compare morphologies**, any of the six) or the CLI (`wrinklefe compare`: stack, convex and concave)
- **Graded averaging:** Through-thickness ply-averaged knockdown for graded wrinkles
- **Movable wrinkle position:** Configurable through-thickness placement (`wrinkle_z_position`, 0.5 = mid-plane)
- **Multi-wrinkle configurations:** Arbitrary N-wrinkle layouts via `AnalysisConfig.wrinkles` — a list of `WrinkleSpec(amplitude, wavelength, width, ply_interface, phase_offset)` — through the analytical, FE, penetration-gate (per-spec, weakest-link) and CZM paths
- **Cohesive-zone delamination (CZM):** Bilinear traction–separation interface elements (`enable_czm=True`) with per-interface damage, energy and crack-length reporting — including continuous cohesive surfaces across adjacent wrinkles for crest-to-crest delamination link-up (see `examples/08_multi_wrinkle_czm_linkup.py`)
- **Resin-pocket material zone:** Graded neat-epoxy lens at the wrinkle crest (modulus + fibre-angle blend, counted once) via `wrinklefe.core.resin_pocket`
- **Tool-flat surfaces & surface resin pockets:** Parts cured against rigid tooling / a caul sheet keep perfectly flat outer surfaces while the fibres undulate internally; the wrinkle troughs fill with neat resin under the flat surface (`enable_surface_resin_pockets`, `surface_pocket_side` = `top`/`bottom`/`both`). FE-only, trough-following, volume-conserving; composes with the crest lens
- **Compaction Vf / ply-thickness gradient:** Under rigid tooling the resin is squeezed out of the compacted regions and pools where the geometry opens up, so ply thickness and local fibre volume fraction vary through the wrinkle (`enable_vf_gradient`, `tool_flat` only — `surface_pocket_side="both"` is the two-caul-plate case). `Vf_local = vf_nominal · h0/h` per element (fibre content conserved), and the local ply is the preset **scaled by the micromechanics Vf ratio** (`wrinklefe.core.compaction`); FE-only, opt-in, and the continuous generalization of the binary surface pockets
- **Progressive-damage FE:** Load-stepping `ProgressiveDamageSolver` to ultimate load with optional crack-band (Bažant–Oh) regularization — the only path that carries UD compression past first-ply failure
- **Penetration gate (θ, D/T, z):** Closed-form two-parameter UD predictor `KD = 1 − (1 − KD_angle(θ))·S(D/T)·P(z)` with calibrated presets; zero FE cost (`wrinklefe.core.penetration_gate`)
- **Linear buckling:** Geometric-stiffness eigenvalue solve (`LinearBucklingSolver`) with a microbuckling knockdown. It is verified structural-buckling infrastructure, **not** the UD wrinkle predictor (the penetration gate is): a homogenised-continuum eigenvalue misses the *fibre-scale* knockdown, and its bifurcation load even *rises* with the wrinkle; see the [modelling findings](docs/wrinkle_modeling_findings.md)
- **16 built-in laminate materials:**
  - carbon/epoxy: AS4/3501-6, IM7/8552, T300/914, T700/2510, T700GC/2510 (AGATE allowables), T800S/M21, M21E/IMA, IM10/8552 and IM6G/3501-6 (the Hsiao & Daniel 1996 wavy-UD study);
  - three dataset-specific carbon cards: `T700_SHI2025`, `UD_CFRP_PILATO2022` and `YC300_LY5052`;
  - glass and aramid: AC318/S6C10 S-glass/epoxy (moulded, plus `AC318_S6C10_vacbag`, the Li 2025 vacuum-bag realization with measured Xc = 335.5 MPa and E1 = 50.8 GPa), S-2 glass/epoxy and Kevlar-49/epoxy;
  - plus an isotropic neat-epoxy card (`EPOXY_S6C10`) for the resin-pocket zone.
- **Applicability warnings:** every analysis reports how far the wrinkle moves the laminate's load path (e/t) and flags a bending-dominated wrinkle. On a bending-dominated whole-thickness wave in compression the analytical knockdown is also capped by an eccentric-column bending bound, and the result says which bound governed (issue #452). `manufacturing_process` flags a laminate made outside the validated processes (wet lay-up, infusion, RTM), and a penetration-gate preset used outside its fitted material or process is flagged too. The warnings reach the summary, the exports and the app; apart from the bending bound, they change no number.
- **Validation ledger:** 13 published datasets (A–M) from 10 research groups; every one but Dataset E is in the ledger with each prediction pinned, so a model change shows up as test drift; see [Validation](docs/validation.md)
- **Inverse goal-seek (acceptance limits):** `wrinklefe.goalseek.find_critical_value` / `wrinklefe critical` answer the disposition question directly — the largest wrinkle amplitude (or any float config field) that still meets a target knockdown or allowable strength. The answer is backed off to the safe side and verified by a real forward run, and the search refuses with an actionable message when the curve admits no unique root
- **Process-parallel sweeps:** `n_workers=N` on both sweep APIs and `wrinklefe sweep --parallel N` fan the independent per-point solves across CPU cores (results identical to and ordered like the sequential run)
- **Comprehensive test suite** covering all modules (run `pytest` to see the current count)

## Installation

WrinkleFE needs Python 3.10 or newer. It is not on PyPI yet, so install
it from GitHub:

```bash
pip install "git+https://github.com/elhajjar1/wrinkleFE.git"
```

Once a release is published on PyPI, `pip install wrinklefe` will work
too.

To contribute, or to run the Streamlit app locally, clone the source:

```bash
git clone https://github.com/elhajjar1/wrinkleFE.git
cd wrinkleFE
pip install -e ".[all]"
```

The 3D cohesive-zone renders (`plot_interface_damage_3d` /
`plot_crack_front_3d`) use PyVista, which pulls in VTK (~150 MB). It is an
optional dependency, so the plain install stays lean and headless-safe;
add the `vtk` extra when you need those plots
(`pip install "wrinklefe[vtk] @ git+https://github.com/elhajjar1/wrinkleFE.git"`; it is already included in
`[all]`).

Verify the install:

```bash
python -c "import wrinklefe; print('WrinkleFE installed successfully')"
```

Run the test suite:

```bash
pytest
```

## Usage

### Streamlit web app

In the [hosted app](https://wrinklefe.streamlit.app/), pick a material
from the sidebar (or, in **Expert mode**, **Custom…** to enter your own
elastic constants and strength allowables), enter a layup in contracted
notation (e.g. `[0/45/-45/90]_3s`), set the wrinkle geometry, and click
**Run analysis**. The app ships with
per-morphology schematic cartoons (including the `tool_flat`
surface-resin-pocket morphology, whose pinned-side and transition-ply
controls appear inline with the Morphology selector), a live wrinkle
preview, and the same analytical + FE pipeline as the Python API. Expert
mode also exposes the through-width `transverse_mode` envelope, and the
sidebar **Config file** section can **download** the current settings as a
portable `AnalysisConfig` JSON — round-tripping with the CLI `--config` /
`--save-config` flags — or **load** a saved case (JSON/YAML) back into the
sidebar. **Share this case (link)** turns the same settings into a
[permalink](#sharing-a-case-as-a-link-permalinks) you can paste into a
chat, a bookmark or a bug report.

A run reports live progress while it solves: the status box tracks the
pipeline phase by phase (building the laminate → analytical predictions →
assembling the FE mesh → solving → evaluating failure → retention
factors), so a long FE or CZM solve is distinguishable from a hung one.
Re-running an identical configuration is served from a small in-session
result cache instead of re-solving; *Reset to defaults* clears it.

A second sidebar button, **Find acceptable limit**, runs the inverse
search on those same inputs — the browser form of
[`wrinklefe critical`](#finding-the-maximum-acceptable-defect). Set a
target knockdown (or an absolute MPa allowable) in the
*Acceptable-limit settings* expander and the app reports the largest
amplitude the laminate can carry and still meet it, together with the
scanned knockdown curve and the per-evaluation ledger. The reported
number is the conservative one — backed off from the root and verified
by a real forward run — and one click applies it back to the sidebar so
you can run the full analysis at the limit. When the search converged
for the *same* configuration a run was made on, the limit is carried
onto the NCR validation summary on the **Export** tab; a limit searched
against different inputs is deliberately left off. If the curve admits
no unique answer the app shows the engine's diagnosis and the measured
curve instead of a number. See
[`DEPLOYMENT_STREAMLIT.md`](docs/internal/DEPLOYMENT_STREAMLIT.md) for the full feature
tour and instructions for self-hosting.

To run the app locally, from a clone:

```bash
pip install -r requirements.txt
pip install -e .
streamlit run app.py
```

The app's interactive 3D views are ordinary library functions, so you can
draw the same figures in a notebook — the mesh surface, the σ contour, the
deformed shape, the failure-index surface and the y-slice scatters:

```bash
pip install "wrinklefe[plotly] @ git+https://github.com/elhajjar1/wrinkleFE.git"
```

```python
from wrinklefe.viz import mesh3d_figure, stress_contour_figure

fig = stress_contour_figure(nodes, elements, stress_per_elem)
fig.show()
```

Plotly is optional: `import wrinklefe.viz` works without it, and only
asking for one of these figures raises an error telling you what to
install. See the [`wrinklefe.viz` API reference](docs/api/viz.rst).

### Python API

```python
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis

config = AnalysisConfig(
    amplitude=0.366, wavelength=16.0, width=12.0,
    morphology="stack", loading="compression",
)
result = WrinkleAnalysis(config).run()
print(result.summary())
```

Runnable scripts for the common workflows — parametric sweeps,
morphology comparison, CZM delamination, the penetration gate,
progressive damage, the crest resin lens, stochastic propagation, export
round-trips, custom materials — live in [`examples/`](examples/); each
states its expected runtime and output, and CI executes them all so they
stay current.

Three pages carry most of what a new user needs:

- [Tutorial](docs/tutorial.md) — a worked walkthrough from a micrograph
  measurement to the NCR validation attachment.
- [Interpreting results](docs/interpreting_results.md) — what each
  headline number means, what it does not mean, which one governs, and
  the severity bands (`src/wrinklefe/io/export.py` remains authoritative
  for those).
- [Units & conventions](docs/units_conventions.md) — mm / MPa / N/mm,
  strain and `Vf` as fractions, angle and coordinate conventions, signs.

The full API reference and user guide are built from [`docs/`](docs/)
with Sphinx (`pip install -e ".[docs]" && sphinx-build -W docs
docs/_build`) and published at
<https://wrinklefe.readthedocs.io>.

`amplitude` (`A`) is the **half-amplitude** [mm]: the peak displacement
of the wrinkled mid-surface from the flat (unwrinkled) reference plane,
so `z(x) = A·cos(2πx/λ)` (modulated by the envelope) and the
peak-to-trough height is `2A`. For a measured wrinkle (e.g. from a
cross-section micrograph or CT slice), `A = (z_max − z_min) / 2`. The
peak fibre misalignment angle scales as `θ_max ≈ arctan(2πA/λ)`, which
drives the Budiansky-Fleck compressive knockdown.

### Wrinkle geometry parameters

All wrinkle length parameters use a single, consistent unit:
**millimetres (mm)** — the same unit as `ply_thickness` and
`domain_length` (the default `amplitude=0.366` mm is exactly two ply
thicknesses of 0.183 mm). Lengths are **not** normalized by thickness.
The longitudinal coordinate `x` runs along the laminate in the fibre
direction; out-of-plane displacement `z(x)` is measured from the flat
(undeformed) mid-surface. Angles are in **radians**.

This table is the canonical reference for every wrinkle-geometry
parameter exposed by `AnalysisConfig`, the CLI, and the Streamlit UI.
The `AnalysisConfig` docstrings, CLI `--help`, and Streamlit `help=`
tooltips all mirror these definitions; the
`tests/test_param_docs_match.py` regression test pins the defaults so
the docs and the dataclass cannot drift.

| Parameter | Units | Default | Definition | Constraint |
|-----------|-------|---------|------------|------------|
| `amplitude` (A) | mm | `0.366` | Half-amplitude: peak displacement of the wrinkled mid-surface from the flat reference, with `z(x) = A·cos(2π(x − x₀)/λ)` (modulated by the envelope) and peak-to-trough height `2A`. For a measured wrinkle, `A = (z_max − z_min)/2`. | ≥ 0 (`0` = flat / no wrinkle) |
| `wavelength` (λ) | mm | `16.0` | Spatial period of the `cos(2π(x − x₀)/λ)` carrier along the longitudinal x-direction (crest-to-crest distance). Wavenumber `k = 2π/λ`. | > 0 |
| `width` (w) | mm | `12.0` | Longitudinal envelope decay length about the centre `x₀`. Exact meaning is profile-dependent: Gaussian 1/e length scale in `exp(−(x−x₀)²/w²)`, tapered flat-top extent (`\|x−x₀\| < w/2`), or triangular half-base (`\|x−x₀\| < w`). Also used as the transverse (y-direction) extent of the wrinkle in 3-D dual-wrinkle / graded mesh deformation. | > 0 |
| `phase` (φ) | rad | `None` | Explicit dual-wrinkle phase offset between the two wrinkle centrelines. `None` derives φ from `morphology` via `MORPHOLOGY_PHASES` (stack φ=0, convex φ=+π/2, concave φ=−π/2). A float overrides the named-morphology phase so arbitrary offsets can be swept (e.g. `0` to `π`). Ignored for single-wrinkle morphologies (`uniform`, `graded`). | finite when set |
| `decay_floor` | dimensionless | `0.0` | Graded morphology only: minimum fraction of the wrinkle amplitude retained at the laminate outer surfaces. `0.0` = full decay to zero amplitude at the surfaces (pure graded); `1.0` = no decay (equivalent to `uniform`). | in `[0, 1]` |
| `amplitude_profile` | name | `"constant"` | Spatially varying in-plane modulation of the wrinkle amplitude `A`, applied on top of the wrinkle's own longitudinal envelope. `"constant"` (default) preserves the legacy uniform `A`; `"gaussian"` multiplies `A` by `exp(−(s/d)²)`; `"linear"` multiplies `A` by `max(0, 1 − \|s\|/d)` (clipped). `s` is the coordinate from the wrinkle centre along `amplitude_profile_axis` and `d` is `amplitude_profile_decay_length`. | one of `constant`, `gaussian`, `linear` |
| `amplitude_profile_decay_length` | mm | `None` | Decay length `d` (mm) for the Gaussian sigma or linear-decay extent. `None` falls back to the wrinkle profile's own `width`. Ignored when `amplitude_profile == "constant"`. | finite and > 0 when set |
| `amplitude_profile_axis` | axis | `"x"` | In-plane axis along which the amplitude modulation runs. Pick `"y"` for an independent transverse tapering of `A` that does not stack with the existing longitudinal envelope on `x`. | one of `x`, `y` |
| `transverse_mode` | name | `"uniform"` | Through-width (transverse, y-direction) wrinkle-surface envelope `f(y)`. `"uniform"` (default) builds the bare x-only wrinkle. The non-uniform modes wrap the profile in a `WrinkleSurface3D` so the crest amplitude varies across the specimen width: `"gaussian_decay"` decays it toward the edges (localized mid-width defect), `"sinusoidal_y"` ripples it across the width, and `"elliptical"` confines it to a mid-width patch. **FE-only** — a non-uniform mode requires `analytical_only=False` and is not yet combinable with multi-wrinkle (`wrinkles`) or `enable_czm` (both rejected at construction). | one of `uniform`, `gaussian_decay`, `sinusoidal_y`, `elliptical` |
| `transverse_span` | mm | `None` | Specimen width `span_y` seen by the transverse envelope. `None` tracks `domain_width` so the envelope always spans the meshed y-extent. Ignored when `transverse_mode == "uniform"`. | finite and > 0 when set |
| `transverse_width` | mm | `None` | Transverse localization half-width `width_y`: the Gaussian 1/e length for `"gaussian_decay"` and the ellipse half-width for `"elliptical"` (ignored by `"uniform"`/`"sinusoidal_y"`). `None` resolves to `span_y / 4` — a localized mid-width patch whose amplitude has fallen to `exp(−4) ≈ 0.018` of the crest at the edges (`gaussian_decay`) or that occupies the central half of the width (`elliptical`). | finite and > 0 when set |
| `surface_transition_plies` | count | `2` | `tool_flat` morphology only: number of plies over which the amplitude ramps linearly from the full-amplitude core to zero at the pinned (tool-flat) surface. A short transition concentrates the full-amplitude trough mismatch into a thin surface band (significant resin pockets), but the crest-side transition elements compress by `amplitude / surface_transition_plies`, so the amplitude is bounded: `amplitude ≤ 0.8 · surface_transition_plies · ply_thickness / nz_per_ply` (enforced at construction; exceed it and the message names both remedies). | integer ≥ 1 |

Peak fibre misalignment: `θ_max ≈ arctan(2πA/λ)` (exact for a pure
cosine; dimensionless because A and λ share the mm length unit). See the
`WrinkleProfile` class docstring in `src/wrinklefe/core/wrinkle.py` for
the full per-profile geometric definitions.

#### Through-width (transverse) variation

By default a wrinkle is treated as uniform across the full specimen width
(`transverse_mode="uniform"`). Real manufacturing wrinkles are usually
*localized* — high amplitude mid-width, fading toward the edges — so a
uniform-width assumption overstates the defect volume and biases the
knockdown conservative. Set `transverse_mode` to `"gaussian_decay"`,
`"sinusoidal_y"`, or `"elliptical"` to run the FE path with a full
`z(x, y)` wrinkle surface (`transverse_span` and `transverse_width` tune
the width envelope). At the same crest amplitude a localized wrinkle
predicts a *milder* knockdown than the uniform baseline. This is FE-only
and single-wrinkle (analytical-only, multi-wrinkle, and CZM runs are
rejected at construction). The CLI exposes it on `analyze
--transverse-mode`, the app in Expert mode. See
`examples/transverse_wrinkle_knockdown.py` for a runnable comparison.

### Thermal / cure-residual loading (`delta_T`)

A laminate leaves the autoclave stress-free at cure temperature and is
used at room temperature. Cooling down locks in residual stress, because
each ply's transverse (matrix) contraction is restrained by its
neighbours' stiff fibres. That residual stress superposes directly onto
the matrix tension/compression states the failure criteria evaluate, so
it is not a second-order detail: on the IM7/8552 quasi-isotropic coupon
below a `delta_T = -155` cool-down adds **+34.3 MPa** of transverse
tension to every ply — 55 % of `Yt = 62.3 MPa`.

> **Sign convention.** `delta_T` is the temperature change **from the
> stress-free (cure) state**: `delta_T = T_service − T_stress_free`.
> **A cure cool-down is NEGATIVE.** A 177 °C cure taken to 22 °C service
> is `delta_T = -155`, *not* `+155` and *not* `22`. A positive value
> means the laminate is hotter than its stress-free state. Get this
> backwards and the residual matrix stress flips from tension to
> compression — a laminate prone to cure microcracking will look safe.

| Parameter | Units | Default | Definition | Constraint |
|-----------|-------|---------|------------|------------|
| `delta_T` | °C | `0.0` | Uniform temperature change from the stress-free (cure) state, `T_service − T_stress_free`. Negative for a cure cool-down. Adds the CLT thermal resultants `N^T`/`M^T` to the ABD solve on the analytical path and the element thermal initial-strain load vector `∫ Bᵀ C ε_th dV` on the FE path; stresses are recovered from the mechanical (not total) strain on both. `0.0` (default) leaves every result bit-identical. | finite, `\|delta_T\| ≤ 1000` |

```python
config = AnalysisConfig(
    amplitude=0.5, wavelength=16.0, width=12.0,
    morphology="graded", loading="compression",
    angles=[0, 45, -45, 90, 90, -45, 45, 0],
    ply_thickness=0.183,
    delta_T=-155.0,       # 177 °C cure -> 22 °C service
)
result = WrinkleAnalysis(config).run()
print(result.summary())          # states the ΔT and its sense
print(result.failure_report.critical_mode)   # -> matrix_tension
```

```bash
wrinklefe analyze --delta-T -155 --angles "[0/45/-45/90]s"
```

**Both paths honour it (issue #273).** The analytical path adds the CLT
thermal resultants to the ABD solve; the FE path assembles the element
thermal initial-strain load vector `∫ Bᵀ C ε_th dV` and subtracts the
thermal strain during stress recovery, so `σ = C̄ (B u − ε_th)`. Because
the wrinkle rotates the fibre frame, the CTE mismatch concentrates in
exactly the elements where the failure criteria are evaluated. A flat
laminate solved this way reproduces the closed-form CLT ply stresses to
better than 0.5 %.

Two deliberate asymmetries: the pristine **retention baseline** is
solved at the same `delta_T` (so the retention factor compares like with
like), while the **measured modulus** is solved at `delta_T = 0` — a
strain-independent reaction offset from a residual load is not a
stiffness change, and folding it in would report a spurious modulus
shift.

### Combined load states — biaxial and in-plane shear

Real wrinkle dispositions are rarely pure uniaxial: skin panels see
compression **plus shear**, pressure shells see biaxial membrane states.
`AnalysisConfig.load_state` takes a `LoadState` and applies its membrane
resultants directly (issue #275):

```python
from wrinklefe.core.laminate import LoadState

config = AnalysisConfig(
    amplitude=0.15, wavelength=12.0, width=8.0, morphology="graded",
    angles=[0, 45, -45, 90, 90, -45, 45, 0], ply_thickness=0.125,
    load_state=LoadState(Nx=-800.0, Nxy=250.0),   # compression + shear
)
result = WrinkleAnalysis(config).run()
print(result.load_state_factor)             # scale to first failure
print(result.load_state_factor_knockdown)   # vs the flat baseline
```

This matters because the transverse and shear components a combined state
produces are what drive the **matrix** failure modes — the cases where the
FE path earns its cost over the analytical one, and exactly what a
uniaxial-only surface could not express.

**Strength under a combined load** is a *proportional load factor*: the
scalar the whole state is multiplied by to reach first failure. Scaling
the state as a whole is what "how much of this load can it take" means
when the load is not a single number, and it reduces to the usual
definition for a uniaxial state. `load_state_factor_knockdown` is the
ratio to the flat baseline solved under the identical state, so it is a
like-for-like knockdown and is invariant to how hard you push.

`load_state=None` (the default) leaves the `applied_strain` path
bit-identical.

**Membrane components only.** `Mx`/`My`/`Mxy` are rejected: the boundary
mapping applies curvature as a prescribed displacement on the same faces a
membrane state loads with traction, so the two cannot be superposed.
`Qx`/`Qy` have no mapping. A load state is also refused with
`analytical_only=True`, `enable_czm`, or `enable_progressive_damage` —
each builds its own uniaxial boundary conditions, so the state would be
silently ignored. Set the temperature on `AnalysisConfig.delta_T`, not on
the load state: one quantity, one owner.

Moisture (`LoadState.delta_C`, and the `beta1/2/3` swelling coefficients
on every material preset) is deliberately **not** exposed on
`AnalysisConfig` yet: nothing in the CLT solve consumes `delta_C`, so a
config field would be a silent no-op.

### Tension analysis

```python
config = AnalysisConfig(
    amplitude=0.366, wavelength=16.0, width=12.0,
    morphology="stack", loading="tension",
    angles=[0, 45, 90, -45, 0, 45, -45, 0, 0, -45, 45, 0, -45, 90, 45, 0],
    ply_thickness=0.152,
)
result = WrinkleAnalysis(config).run()
print(result.summary())
```

### Graded morphology (embedded wrinkle)

```python
config = AnalysisConfig(
    amplitude=0.5, wavelength=15.0, width=11.0,
    morphology="graded", decay_floor=0.0,
    loading="compression",
)
result = WrinkleAnalysis(config).run()
print(result.analytical_knockdown)
```

### Resin pocket + progressive-damage FE (UD compression)

For a unidirectional wrinkle, tag the soft neat-epoxy lens at the crest
and load-step to ultimate with the progressive-damage solver:

```python
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis
from wrinklefe.core.material import MaterialLibrary

lib = MaterialLibrary()
config = AnalysisConfig(
    amplitude=0.366, wavelength=16.0, width=12.0,
    morphology="graded", loading="compression",
    material=lib.get("AC318_S6C10"),
    enable_resin_pocket=True,                       # graded epoxy lens at the crest
    resin_pocket_material=lib.get("EPOXY_S6C10"),   # default if left None
    enable_progressive_damage=True,                 # load-step to ultimate
    progressive_n_increments=15,
)
result = WrinkleAnalysis(config).run()
print(result.progressive_knockdown, result.progressive_strength_MPa)
```

### Tool-flat surfaces & surface resin pockets (FE)

Parts cured against rigid tooling (or under a caul sheet) keep
**perfectly flat outer surfaces**: the fibre undulation is confined to
the interior, and where the outermost undulating ply dips away from the
flat surface the gap fills with **neat resin** — surface-visible pockets
over the wrinkle troughs, thinning to nothing over the crests. WrinkleFE's
default through-thickness decay already leaves the outer surfaces exactly
flat (for `stack`/`convex`/`concave`, or `graded` with `decay_floor=0`);
`enable_surface_resin_pockets` supplies the missing *material* — it tags
the stretched transition elements as fibre-free isotropic resin (a
stiffness hole and a matrix-cracking site where fibre-misalignment
criteria are meaningless).

```python
config = AnalysisConfig(
    amplitude=0.354, wavelength=7.4, width=3.7,
    morphology="graded", loading="compression",   # graded, decay_floor=0 ⇒ flat surface
    material=lib.get("AC318_S6C10"), angles=[0.0] * 15, ply_thickness=0.42,
    enable_surface_resin_pockets=True,             # trough pockets under the flat surface
    surface_pocket_side="both",                    # "top" | "bottom" | "both"
    resin_pocket_material=lib.get("EPOXY_S6C10"),  # default if left None
)
result = WrinkleAnalysis(config).run()
print(result.modulus_retention_global)
```

The pocket geometry is volume-conserving (the tagged resin equals the
integrated kinematic gap between the flat surface and the outermost
undulating ply) and reuses the crest-lens material plumbing, so the two
zones **compose** (per-element maximum) when both are enabled. This is an
**FE-only** effect: the closed-form analytical path keeps using fibre
angles only. A `uniform` morphology (never flat) or `graded` with
`decay_floor > 0` (wavy surface) is rejected with a message naming the
fix.

#### The `tool_flat` morphology (significant pockets)

Under the smooth-decay morphologies above (`stack`/`convex`/`concave`,
or `graded` with `decay_floor=0`) the amplitude tapers **linearly** across
the whole thickness, so the outermost undulating ply barely moves — at the
24-ply defaults the trough gap is only ~0.25 of a *single* ply thickness,
invisible in the preview and mechanically negligible. Those morphologies
model a co-cured *gradual* wrinkle, not a tooling-dominated one.

`morphology="tool_flat"` models the tooling-dominated case directly: a
uniform-amplitude core (like `uniform`), a short linear transition over
`surface_transition_plies` (default 2), and an **exactly-flat pinned
surface** on `surface_pocket_side` (`"top"`/`"bottom"`/`"both"`). The
kinematic mismatch across the transition equals the *full* amplitude, so
the trough pocket is ≈ `amplitude` deep (≈2.7 ply thicknesses at
defaults) — significant and visible. Surface pockets **auto-enable** for
`tool_flat` (they are its defining physics), and the analytical path
takes `uniform`'s kink-band knockdown (M_f = 1.0; the pocket effect is
FE-only) but not its eccentric-column bending bound, which applies only
to a wave through the whole thickness (issue #452).

```python
config = AnalysisConfig(
    morphology="tool_flat", amplitude=0.35, wavelength=16.0, width=12.0,
    surface_pocket_side="both",       # tool-flat face(s), also the pocket side
    surface_transition_plies=3,       # amplitude ramp width (>= 1)
    material=lib.get("IM7_8552"), angles=[0.0] * 24, ply_thickness=0.183,
)
result = WrinkleAnalysis(config).run()   # pockets auto-enabled
print(result.modulus_retention_global)
```

On the crest side the `surface_transition_plies` transition elements
*compress* by `amplitude / surface_transition_plies`, so the amplitude is
bounded — `amplitude ≤ 0.8 · surface_transition_plies · ply_thickness /
nz_per_ply`; beyond it the elements would invert and construction fails
with a message naming both remedies (more transition plies, or a smaller
amplitude).

Those compressed crest-side elements are exactly what the **compaction
Vf gradient** models (issue #379): set `enable_vf_gradient=True` and the
local fibre volume fraction follows the thickness change,
`Vf_local = vf_nominal · h0/h`, so the crest band compacts and stiffens
while the trough turns resin-rich and softens — one continuous field
instead of a binary neat-resin tag:

```python
config = AnalysisConfig(
    morphology="tool_flat", amplitude=0.25, surface_pocket_side="both",
    enable_vf_gradient=True,          # opt-in; tool_flat only
    vf_nominal=None,                  # None -> the card's documented Vf
    vf_max=0.75,                      # compaction cap (square packing)
    material=lib.get("IM7_8552"), angles=[0.0] * 24, ply_thickness=0.183,
    analytical_only=False,
)
```

The local card is the **preset scaled by the micromechanics ratio**
`P_micro(Vf_local) / P_micro(vf_nominal)` for the stiffnesses and CTEs, so
the trend is modelled without importing the micromechanics model's
absolute error. **Poisson ratios and all strengths stay at the preset
values** — no mixing rule predicts a strength from `Vf` (local failure
indices still move, because the local stiffness redistributes stress).
Elements compacted past `vf_max` saturate, counted in a single warning;
the rule carries no lateral resin flow along the ply. When the gradient is
on it supersedes the binary surface-pocket tag (it is its continuous
generalization); the machined crest resin lens composes unchanged.

### Penetration gate (θ, D/T, z) — UD, zero FE cost

The closed-form two-parameter gate predicts a UD knockdown directly from
geometry. Call it on its own with a calibrated preset:

```python
from wrinklefe.core.penetration_gate import penetration_gate_kd, GATE_LI2024_MOULDED

kd = penetration_gate_kd(theta_deg=8.0, dt=0.10, params=GATE_LI2024_MOULDED)
print(kd)
```

Or drive it through `AnalysisConfig.penetration_gate` so
`analytical_knockdown` (and `analytical_strength_MPa`) come from the gate
instead of Budiansky–Fleck (use `GATE_LI2025_VACBAG` with
`AC318_S6C10_vacbag` for the vacuum-bag realization):

```python
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis
from wrinklefe.core.material import MaterialLibrary
from wrinklefe.core.penetration_gate import GATE_LI2024_MOULDED

lib = MaterialLibrary()
config = AnalysisConfig(
    amplitude=0.366, wavelength=16.0, width=12.0,
    morphology="uniform", loading="compression",
    material=lib.get("AC318_S6C10"),
    penetration_gate=GATE_LI2024_MOULDED,
)
result = WrinkleAnalysis(config).run(analytical_only=True)
print(result.analytical_knockdown)
```

When `penetration_gate` is left unset (the default `None`), the
analytical knockdown is unchanged. With a multi-wrinkle configuration
(`AnalysisConfig.wrinkles`, below) the gate evaluates each wrinkle on
its own geometry — `theta_i = arctan(2πA_i/λ_i)`, penetration
`D_i/T = A_i/T`, and the through-thickness position factor `P(z_i)`
from the spec's ply interface — and returns the weakest-link (minimum)
knockdown.

### Multi-wrinkle configurations

Real laminates often carry several wrinkles. Passing a list of
`WrinkleSpec` entries overrides the single/dual-wrinkle dispatch and
places each wrinkle at its own ply interface with its own geometry and
longitudinal position (`phase_offset` shifts a crest by
`phase·λ/2π`):

```python
import numpy as np
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis, WrinkleSpec

config = AnalysisConfig(
    morphology="graded", loading="compression",
    angles=[0.0] * 14, ply_thickness=0.44,
    wrinkles=[
        WrinkleSpec(amplitude=0.75, wavelength=12.9, width=6.45,
                    ply_interface=6, phase_offset=-2.0 * np.pi),
        WrinkleSpec(amplitude=0.75, wavelength=12.9, width=6.45,
                    ply_interface=6, phase_offset=+2.0 * np.pi),
    ],
)
result = WrinkleAnalysis(config).run()
```

The composed displacement and fibre-angle fields feed the FE solve
("compose then differentiate"), the penetration gate scores each
wrinkle and takes the weakest link, and `enable_czm=True` inserts
cohesive surfaces along the full length of every wrinkle-nominated
interface — wrinkles sharing an interface get one continuous surface,
so a delamination can propagate crest-to-crest between neighbours
(`examples/08_multi_wrinkle_czm_linkup.py` demonstrates the link-up
vs far-separated contrast).

### When a CZM solve fails to converge

CZM solves can fail for physical reasons — the tangent goes indefinite
near peak load, the load steps are too coarse for the damage evolution,
or the residual stagnates. When that happens the run still returns
(partial) results with `czm_converged = False`, and the solver
classifies *why* and hands back one actionable hint. `AnalysisResults`
carries `czm_failure_diagnostics` (increment index + load fraction,
iteration count, final residual, residual history, line-search status)
and `czm_failure_hint`; the CLI prints the hint to stderr and the
Streamlit CZM section shows it as an error. The hint names the knob to
turn:

- **still decreasing at the iteration cap** → raise `max_newton_iter` or
  loosen `czm_newton_tol`;
- **stagnated / diverged** → increase `czm_n_load_increments` or reduce
  the applied strain;
- **tangent singular (post-peak snap-back)** → more load increments or a
  smaller strain; displacement/arc-length control is the robust fix on
  the roadmap.

### Uncertainty propagation (probabilistic analysis)

Measured wrinkle geometry is uncertain — amplitude and wavelength come
from a micrograph or C-scan with error. `probabilistic_analysis`
samples `AnalysisConfig` fields from user-supplied distributions
(Latin-hypercube by default, plain Monte-Carlo optional) and runs the
analytical path per sample, turning "the model says 0.64" into "P5–P95
= 0.59–0.86 given my measurement uncertainty" — the form an NCR
disposition rationale actually needs:

```python
from wrinklefe.analysis import AnalysisConfig
from wrinklefe.core.penetration_gate import GATE_LI2025_VACBAG
from wrinklefe.stochastic import probabilistic_analysis

base = AnalysisConfig(
    amplitude=0.75, wavelength=12.9, width=6.45,
    angles=[0.0] * 14, ply_thickness=0.44, morphology="graded",
    penetration_gate=GATE_LI2025_VACBAG,
)
prob = probabilistic_analysis(
    base,
    {"amplitude": ("normal", 0.75, 0.08),
     "wavelength": ("normal", 12.9, 1.0)},
    n_samples=1000, seed=42,
)
print(prob.summary())                     # P5/P50/P95, mean ± std
print(prob.knockdown_percentile(5.0))     # 5th-percentile knockdown
prob.plot()                               # histogram + sensitivity scatter
```

Distributions accept `("normal", mean, std)`, `("uniform", lo, hi)`,
`("lognormal", mu, sigma)` or any frozen `scipy.stats` distribution; a
fixed `seed` makes the whole analysis reproducible, and `n_workers`
reuses the sweep process pool for FE-path sampling. 1000 analytical
samples run in under a second for UD/gate configs (~20 s for a 24-ply
multidirectional layup, or seconds with `n_workers`).

> **Not A-/B-basis values.** The reported percentiles are
> *model-input-propagation* statistics — the deterministic model driven
> by sampled geometry. They are **not** CMH-17 A-/B-basis allowables
> (one-sided tolerance bounds on physical test data with prescribed
> confidence) and must not be presented as basis values in
> certification paperwork.

### Batch parametric sweeps

For exploring how the knockdown varies across a parameter range, use
`WrinkleAnalysis.parametric_sweep` to sweep a single
`AnalysisConfig` field (any numeric field — `amplitude`, `wavelength`,
`width`, `phase`, `applied_strain`, ...):

```python
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis

base = AnalysisConfig(
    amplitude=0.366, wavelength=16.0, width=12.0,
    morphology="stack", loading="compression",
)
results = WrinkleAnalysis.parametric_sweep(
    base, parameter="amplitude", values=[0.1, 0.2, 0.3, 0.4],
    analytical_only=True,
)
for r in results:
    print(f"A={r.config.amplitude:.3f}  KD={r.analytical_knockdown:.4f}")
```

Every sweep point is an independent analysis, so full-FE sweeps
parallelize across processes: pass `n_workers=N` (`0` = all CPU cores)
and the solves fan out over a process pool with results returned in
the same order as `values` — measured 3.6× at 4 workers on an 8-value
FE sweep. Peak memory scales with `n_workers` × the per-solve
footprint, so size the worker count by available RAM for fine meshes.

For multi-parameter cross-product sweeps with JSON output and plots,
use `wrinklefe.sweep.run_sweep`:

```python
import numpy as np
from wrinklefe.sweep import run_sweep, save_sweep_results, plot_sweep_results

sweep = run_sweep({
    "amplitude":  np.linspace(0.183, 0.549, 3),
    "wavelength": np.linspace(8.0, 24.0, 3),
})
save_sweep_results(sweep, "./sweep_output/")
plot_sweep_results(sweep, "./sweep_output/")
```

### Command line

```bash
wrinklefe --help

# Single-parameter sweep (analytical-only is the default, fast)
wrinklefe sweep --parameter amplitude --min 0.1 --max 0.5 --steps 5

# Full-FE sweep across 4 worker processes (--parallel 0 = all cores)
wrinklefe sweep --parameter amplitude --min 0.1 --max 0.5 --steps 8 \
    --no-analytical-only --parallel 4

# Sweep over a saved config's laminate/gate (the headline UD capability):
# --config supplies the base (material, layup, gate); --parameter/--min/
# --max/--steps drive the variation over it. `compare` takes --config too.
wrinklefe sweep --config ud_gate.json --parameter amplitude \
    --min 0.2 --max 0.9 --steps 6

# Inverse: the largest amplitude that still meets a 0.85 knockdown
wrinklefe critical --parameter amplitude --target-knockdown 0.85
```

Through-width (transverse) wrinkle surfaces and the mesh/thickness knobs
are exposed on `analyze` (a non-uniform `--transverse-mode` forces the FE
path; `--ply-thickness` sets the gate's D/T):

```bash
wrinklefe analyze --transverse-mode gaussian_decay --transverse-width 5 \
    --nz-per-ply 2 --ply-thickness 0.25
```

#### Wrinkle-defect capabilities

`analyze` also exposes the wrinkle-defect models. The FE-only features
(`--resin-pocket`, `--surface-resin-pockets`, `--vf-gradient`,
`--progressive`) force the
full FE solve, so they take precedence over `--analytical-only` / `--no-fe`:

```bash
# Two-parameter (theta, D/T) penetration gate — the best UD predictor —
# selecting a calibrated preset (UD-scoped; not for multidirectional laminates)
wrinklefe analyze --gate li2025-vacbag --morphology uniform \
    --angles 0,0,0,0,0,0,0,0 --interface-1 3 --interface-2 4 --amplitude 0.5

# Off-mid-plane wrinkle position (fraction of thickness, in [0, 1])
wrinklefe analyze --morphology graded --wrinkle-z-position 0.7

# Crest resin pocket + progressive-damage ultimate strength (FE)
wrinklefe analyze --resin-pocket --progressive --increments 15

# Surface resin pockets under a tool-flat surface (FE)
wrinklefe analyze --surface-resin-pockets --surface-pocket-side both

# tool_flat morphology — significant pockets (auto-enabled). The
# hyphenated --morphology tool-flat is an alias for tool_flat; amplitude
# is bounded by 0.8 * surface-transition-plies * ply_thickness / nz_per_ply.
wrinklefe analyze --morphology tool-flat --amplitude 0.35 \
    --surface-pocket-side both --surface-transition-plies 3

# Compaction Vf / ply-thickness gradient — two caul plates (FE, tool_flat)
wrinklefe analyze --morphology tool-flat --amplitude 0.25 \
    --surface-pocket-side both --vf-gradient
```

| Flag | Config field | Notes |
| --- | --- | --- |
| `--wrinkle-z-position Z` | `wrinkle_z_position` | Fraction of thickness in `[0, 1]` (0.5 = midplane); graded morphology **and the penetration gate on every morphology** — see note below |
| `--gate {li2024-moulded,li2025-vacbag}` | `penetration_gate` | Calibrated `GateParameters` preset; UD-scoped |
| `--process {prepreg-autoclave,…,wet-layup,unspecified}` | `manufacturing_process` | How the laminate was made. A process outside the validated scope (anything but the three prepreg routes) adds a warning, as does a `--gate` preset used outside its fitted material or process; no number changes (issue #445) |
| `--resin-pocket` | `enable_resin_pocket` | Crest resin lens (FE-only) |
| `--surface-resin-pockets` / `--surface-pocket-side {top,bottom,both}` | `enable_surface_resin_pockets` / `surface_pocket_side` | Tool-flat surface pockets (FE-only). Auto-enabled for `--morphology tool-flat` |
| `--surface-transition-plies N` | `surface_transition_plies` | `tool_flat` only: amplitude-ramp width (≥ 1, default 2); bounds the amplitude |
| `--vf-gradient` | `enable_vf_gradient` | Compaction Vf / ply-thickness gradient (FE-only, `tool_flat` only); supersedes the binary surface pockets |
| `--vf-nominal V` / `--vf-fiber F` / `--vf-matrix M` | `vf_nominal` / `vf_fiber` / `vf_matrix` | Ratio anchor and constituents; omit to use the material card's documented values |
| `--progressive` / `--increments N` | `enable_progressive_damage` / `progressive_n_increments` | Load-stepping ultimate strength (FE-only) |

> **`--wrinkle-z-position` is not gate-neutral.** Outside the gate it is
> genuinely inert for `stack` / `convex` / `concave` / `uniform` — those
> morphologies fix their own through-thickness behaviour. But the gate's
> position factor `P(z) = (2·min(z, 1−z))^position_q` is applied on
> **every** morphology, and with `li2025-vacbag` (`position_q = 5.26`) it
> is steep: on the 14-ply UD case above, moving z from 0.5 to 0.9 takes the
> knockdown from 0.6427 to 0.9999 — a **120 MPa less conservative**
> allowable. Leave it at the midplane default unless you have a micrograph
> that says otherwise.

The long tail of finer knobs (custom `GateParameters`, resin-pocket
geometry scales, progressive load-ramp targets) stays reachable through
`--config` (see below). `sweep` sweeps any numeric `AnalysisConfig` field,
including `wrinkle_z_position` (add `--config` to sweep over a saved setup):

```bash
wrinklefe sweep --parameter wrinkle_z_position --min 0.2 --max 0.8 \
    --steps 4 --morphology graded
```

### Iterative solver controls

`solver="iterative"` runs conjugate gradient with an incomplete-LU (ILU)
preconditioner — memory-efficient for large meshes (>100 k DOFs). Its
knobs are `AnalysisConfig` fields, reachable through `--config`:

| Config field | Default | Meaning |
| --- | --- | --- |
| `iterative_rtol` | `1e-10` | CG relative-residual convergence tolerance (> 0) |
| `iterative_maxiter` | `10000` | CG iteration cap (≥ 1) |
| `ilu_drop_tol` | `1e-4` | ILU drop tolerance — the main quality/memory knob; larger = sparser, cheaper, weaker (≥ 0) |
| `ilu_fill_factor` | `None` | Upper bound on ILU fill; `None` keeps SciPy's default (≥ 1 when set) |
| `preconditioner` | `"ilu"` | `"ilu"`, `"jacobi"` (diagonal — much lower memory for huge meshes), or `"none"` |

If ILU construction fails (out of memory, or a structurally/numerically
singular factor), the solver emits a **`WARNING`** and falls back to
the diagonal (Jacobi) preconditioner — an
ill-conditioned matrix can make that fallback orders of magnitude slower,
so the warning names the original error and points at
`preconditioner="jacobi"` (to silence it) or the direct solver. Only the
narrow set of exceptions `spilu` uses to signal a real factorisation
failure is caught; any other error propagates. On non-convergence the
`RuntimeError` names the active preconditioner, the iterations used, the
cap, and the final relative residual.

### Saving and reusing a configuration

`analyze` can persist and reload a full `AnalysisConfig`. `--save-config`
writes the *effective* configuration (after any `--config` file and CLI
overrides are applied); `--config` reloads it. Any flag given on the same
command line as `--config` overrides the file value, while flags left off
keep the file's value:

```bash
# Save the effective config to a file, then reuse it verbatim
wrinklefe analyze --amplitude 0.4 --morphology concave --save-config case.json
wrinklefe analyze --config case.json

# Reuse the file but override one parameter (0.9 wins over the file's value)
wrinklefe analyze --config case.json --amplitude 0.9
```

The same round-trip is available programmatically via
`AnalysisConfig.to_dict()` / `from_dict()` and the
`save_json` / `load_json` (and extension-dispatching `save` / `load`)
helpers. The JSON pins a `config_version` field; loading a file with an
unknown key or a mismatched version fails loudly. YAML is supported when
PyYAML is installed (it is not a required dependency). Library materials
serialise by name, custom materials inline, and penetration-gate presets
serialise by their registry name.

### Sharing a case as a link (permalinks)

The app's **Share this case (link)** expander, at the bottom of the
sidebar, gives the current inputs as one URL. Opening it — in a fresh
browser, on someone else's machine, from a bookmark — loads that exact case
into the sidebar: layup, material (including a custom one), geometry,
morphology, loading, mesh and CZM settings. It is the browser counterpart to
a config file, and the two are built from the same inputs, so a link and a
downloaded `.json` of the same sidebar describe the same case.

```
https://wrinklefe.streamlit.app/?cfg=eNo1jEsOgyAQhu_yrwlB21o7lzFERyQ...
```

A few properties worth knowing:

- **Inputs, not results.** A link reopens the case; click **Run analysis**
  to get the numbers. To hand over a finished run, use the full-result
  archive on the Export tab
  ([`.wfr`](#archiving-a-full-result-wfr--reload-without-re-solving)).
- **Short.** Only inputs that differ from the defaults are carried, then
  compressed: a two-parameter case is about 170 characters, and a config
  with every sidebar section changed and a custom material about 740 —
  well inside the limits of chat clients and mail gateways.
- **Bookmarkable.** The `?cfg=` parameter stays in the address bar, so a
  reload or a bookmark reopens the case. It is read **once** per visit: once
  you start editing, your changes are not overwritten. Re-copy the link from
  the expander after changing inputs.
- **Safe to open.** The payload is treated as untrusted input — size-capped,
  bounded on decompression, and validated through `AnalysisConfig.from_dict`
  like any other config. A truncated or tampered link shows a plain-English
  message at the top of the sidebar and leaves the defaults in place; it
  never produces an error page.
- **Says when it may have drifted.** Because only non-default inputs travel,
  a link means "these values, plus the defaults". If a later release changes
  a default, opening an old link still works but warns that any input it
  did not name now uses the current default.

The codec is an ordinary library module, so the same links can be built or
read from Python:

```python
from wrinklefe.analysis import AnalysisConfig
from wrinklefe.io import decode_config, permalink_url

url = permalink_url(AnalysisConfig(amplitude=0.5), "https://wrinklefe.streamlit.app")
config = decode_config(url.split("cfg=", 1)[1]).config
```

### Probabilistic analysis from the command line

`wrinklefe stochastic` propagates input distributions through the analysis
to report percentile knockdowns — the "P5 knockdown under measurement
uncertainty" number — from a `--config` base plus one or more repeatable
`--distribution FIELD:DIST:P1:P2` specs (`DIST` is `normal`, `uniform`, or
`lognormal`). A fixed `--seed` makes the whole run reproducible; results
write to JSON (percentiles + per-sample arrays) and/or a per-sample CSV.

```bash
wrinklefe stochastic --config case.json \
    --distribution amplitude:normal:0.5:0.05 \
    --distribution wavelength:uniform:12:20 \
    --n-samples 1000 --seed 42 --method lhs \
    --output-json prob.json --output-csv prob.csv
```

These are model-**input-propagation** statistics, **not** CMH-17 A-/B-basis
allowables (the printed summary carries the disclaimer).

### Finding the maximum acceptable defect

`wrinklefe critical` inverts the analysis: instead of "given this wrinkle,
what is the knockdown?", it answers "given this allowable, what is the
largest wrinkle we can accept?" — the number that goes into an inspection
criterion or an NCR disposition. Pass `--target-knockdown` (a knockdown
factor) or `--target-strength` (an absolute MPa allowable), optionally
`--parameter` (any float `AnalysisConfig` field; default `amplitude`),
`--objective`, `--bracket LO HI`, `--max-value`, `--scan-points` and
`--rtol`. `--save-config` writes the config at the critical value — the
run to attach to the disposition — and `--output-json` / `--output-csv`
carry the scan curve and the full evaluation ledger.

```bash
# Compression: largest amplitude meeting a 0.85 knockdown
wrinklefe critical --parameter amplitude --target-knockdown 0.85 \
    --save-config limit.json

# Tension, expressed as an absolute allowable instead
wrinklefe critical --loading tension --target-strength 1020
```

**The returned value satisfies the criterion under a real forward
evaluation, not merely to within the root tolerance.** `brentq` converges
to a root, not to the side of it that satisfies the inequality, so the
search backs the answer off toward safety and verifies it with an extra
solve; the printed `Criterion satisfied:` line is that check.

The search runs on the analytical path by default (about 25 forward
evaluations, well under a second; `--bracket LO HI --scan-points 3` cuts
it to about 15). `--no-analytical-only` forces the full FE pipeline, where
each evaluation is 4 linear solves — minutes, not seconds. Root-finding is
inherently sequential, so there is no `--parallel`.

Neither monotonicity nor uniqueness is assumed: the direction is measured
from a log-spaced scan of the search range, and when the target is never
crossed (`no_crossing`), never met (`target_unreachable`), reached by more
than one root (`non_monotonic` — e.g. `morphology="graded"` against
`wavelength`, whose curve is U-shaped), or unresolvable because the
parameter is inert for the configuration (`flat`), the search exits 1 with
a message quoting the measurements behind the refusal — never a scipy
traceback. "No crossing in range" is not an error: it means no defect size
in range fails your criterion.

The same search is the **Find acceptable limit** button in the
Streamlit app (see [Streamlit web app](#streamlit-web-app)).

### Exporting results to CSV / JSON

Numeric outputs (load factor, per-ply failure index, knockdown factors,
stress-field summary) can be written to a schema-versioned JSON or a
Pandas-friendly per-ply CSV for downstream comparison and plotting in
Excel, Pandas, or shared Jupyter notebooks:

```python
from wrinklefe.analysis import AnalysisConfig, WrinkleAnalysis
from wrinklefe.io.results import export_results_csv, export_results_json

result = WrinkleAnalysis(AnalysisConfig()).run()

export_results_json(result, "results.json")  # schema-versioned JSON
export_results_csv(result, "per_ply.csv")    # per-ply tabular CSV
```

The JSON output is deterministic (`sort_keys=True`), pins a top-level
`schema_version` field, and reduces large numpy arrays (e.g.
per-Gauss-point stress fields) to summary statistics so the file stays
compact. The CSV is one row per ply with columns `ply_index,
angle_deg, max_FI, min_RF, critical_mode, critical_criterion`, suitable
for `pandas.read_csv` or `csv.DictReader`.

Every JSON export carries a `provenance` block recording the installed
WrinkleFE version (never a hardcoded literal), the Python/numpy/scipy
versions and the platform — so a result file can be audited and
reproduced against the validation ledger. The schema-versioned export
adds a solver snapshot and leaves out the write time, so identical
input gives a byte-identical file; the NCR validation summary
(`build_analysis_summary`) embeds the block with a UTC timestamp.

Both also carry a `validation_margin` block (issue #455). It gives the
factor that would have kept 95% of the validation cases in the
analytical knockdown's regime conservative: the same model path, load,
laminate type and bending regime. With it come the regime and the
number of cases behind it, or "insufficient evidence" below 10 cases.
It is reporting only; see `docs/interpreting_results.md`.

The Streamlit web app exposes the same exports as **Download results as
JSON** and **Download per-ply results as CSV** buttons on the Export
tab.

On the command line, `analyze`, `sweep`, `compare` and `critical` write
JSON with `--output-json`. Pass `--json-schema structured` to get the
schema-versioned document above (one per run). The default is still the
older layout (top-level `wrinklefe_version`, written by
`wrinklefe.io.export.export_results_json`), which is **deprecated since
1.3 and will be removed in 2.0**; the CLI prints a note each time it
writes it, and the Python functions raise a `DeprecationWarning`. The two
layouts share no result paths, so a consumer written against one reads
nothing from the other.

### Archiving a full result (`.wfr`) — reload without re-solving

The JSON and CSV exports are report tiers: they reduce per-Gauss-point
fields to `{min, max, mean, p95, n}` so the file stays a few KB. That is
the right trade for a report and the wrong one for post-processing —
once the Python session ends, the expensive part of an FE or CZM run
(stress and strain per Gauss point, failure-index and mode fields,
cohesive damage, separation and traction) is gone, and the only way to
plot a different slice tomorrow is to solve again.

`save_results` / `load_results` are the archive tier underneath them:

```python
from wrinklefe.io import load_results, save_results
from wrinklefe.viz import plot_displacement_3d

save_results(result, "run.wfr")

# ... a week later, in a different process ...
restored = load_results("run.wfr")
plot_displacement_3d(restored.field_results)   # viz/ takes it unchanged
```

A `.wfr` is a single compressed `.npz`: each array is stored losslessly
under its dotted path in the result graph
(`field_results.stress_local`,
`failure_report.ply_failure_indices.larc05`), and everything else —
scalars, the config, the provenance block, and the structure needed to
rebuild the objects — goes in a JSON manifest. Round-tripping reproduces
every array bit-for-bit and every scalar exactly, and it preserves the
object graph rather than duplicating it, so `field_results.mesh is
results.mesh` still holds on the reloaded result.

**No pickle.** Arrays are written as plain numeric or unicode dtypes and
read back with `allow_pickle=False`, so opening an archive cannot execute
anything and does not depend on this package's class layout — a `.wfr` is
inspectable with nothing but `numpy` and `json`. The manifest pins an
`ARCHIVE_FORMAT_VERSION`, and a newer major version raises
`ArchiveFormatError` on load rather than guessing: a half-read result is
worse than a refusal, because the numbers would still look plausible.

From the command line, `--save-results` writes the archive alongside any
other output:

```bash
wrinklefe analyze --amplitude 0.6 --wavelength 20 --fe \
    --save-results run.wfr
```

The Streamlit app offers the same file as **Download full result archive
(.wfr)** on the Export tab. `examples/17_archive_reload_results.py`
walks the round-trip end to end, including the bit-identity and
plot-identity checks.

### Building a ply from its constituents (fibre volume fraction)

The 16 built-in laminate cards above are fixed. When the question is
*"what if my fibre volume fraction is 0.55 rather than 0.60?"*,
`wrinklefe.core.micromechanics` builds the ply from its constituents —
a fibre, a neat resin, and a Vf — instead:

```python
from wrinklefe.core.material import MaterialLibrary, OrthotropicMaterial
from wrinklefe.core.micromechanics import FIBER_PRESETS, MATRIX_PRESETS

lib = MaterialLibrary()

ply = OrthotropicMaterial.from_constituents(
    FIBER_PRESETS["IM10"],            # Hexcel HexTow IM10 carbon
    MATRIX_PRESETS["EPOXY_8552"],     # HexPly 8552 neat resin
    0.55,                             # fibre volume fraction
    strengths_from=lib.get("IM10_8552"),   # allowables carried over as-is
)
ply.E1, ply.E2, ply.G12, ply.nu12
```

The rules are the standard ones, each cited in the module docstring:
Voigt rule of mixtures for `E1`, `nu12` and `nu23`, Halpin–Tsai for `E2`
(ξ = 2) and `G12` (ξ = 1), transverse isotropy for `G23`, and Schapery
for the thermal expansion coefficients. `FIBER_PRESETS` covers the nine
fibres behind the library systems (AS4, T300, T700S, IM7, IM10, IM6G,
T800S, S-2 glass, Kevlar 49) and `MATRIX_PRESETS` the two resins with
published neat-resin data (3501-6, 8552); any isotropic card already in
the library — `EPOXY_S6C10`, say — can be used as the matrix through
`MatrixProperties.from_material`. Every constituent constant is sourced
in a comment: `E1f`/`alpha1f` from the manufacturer data sheet, the
transverse set from Daniel & Ishai (2006) Table A.2/A.3.

**Two limitations, both deliberate.**

1. **Strengths are not predicted.** No mixing rule here maps Vf to an
   allowable, so `Xt`/`Xc`/`Yt`/`Yc`/`S12`/… are carried over unchanged
   from `strengths_from` (or left at the defaults) and do *not* track
   Vf. Longitudinal strength is set by fibre-strength statistics and
   misalignment, transverse and shear strengths by the matrix and the
   fibre–matrix interface; a Vf-scaled strength model would be quietly
   wrong.
2. **The elastic predictions are approximate.** Rebuilt from
   constituents at the Vf each preset documents in its own source
   comment, the model lands within 12 % on `E1`, 26 % on `nu12`, 33 %
   on `E2` and 32 % on `G12` — except Kevlar-49/epoxy `G12`, which
   Halpin–Tsai over-predicts by 84 % (published aramid `G12f` values
   scatter by an order of magnitude; that case is a recorded `xfail`
   rather than a tuned fibre constant). `nu23` is under-predicted
   throughout. `tests/test_micromechanics.py` pins every one of those
   deviations. Use the model for the *trend* — how properties move with
   Vf, anchored on a measured ply — not as a source of allowables.

In the analysis pipeline these rules are used only as a ratio, by the
compaction Vf gradient (`enable_vf_gradient`, issue #379, above), where
the local Vf varies from element to element.

## Validation

### What the in-repo tests check

The integration tests under `tests/test_integration/` exercise the full
`WrinkleAnalysis` pipeline and assert physical-sanity properties of the
analytical knockdown rather than reproducing absolute experimental
strengths:

- `test_elhajjar_validation.py` (10 tests): zero-amplitude returns
  knockdown ≈ 1, knockdown decreases monotonically with amplitude,
  morphology ordering convex > stack > concave in compression, knockdown
  stays in `(0, 1]`, and strength equals `Xc * knockdown`.
- `test_tension_validation.py` (13 tests): tension pipeline completes,
  uses `Xt` (not `Xc`), three-mechanism (`kd_fiber`, `kd_matrix`,
  `kd_oop`) decomposition is populated with the controlling mode, and
  tension knockdown is no more severe than compression for the same
  defect.

These act as regression guards on the analytical model. The measured
knockdowns live in the committed validation ledger
([`tests/test_validation/ledger.json`](tests/test_validation/ledger.json)):
85 cases in 16 entries covering Datasets A–D and F–M, each with its
analysis recipe and a pinned prediction (Dataset E, Li et al. 2024, is
scored by the validation scripts only). One command regenerates the
full per-case predicted-vs-measured table, with drift detection against
the pinned baselines:

```bash
python scripts/validate.py
```

For **Li et al. (2025)** (Dataset F, 6 single-wrinkle UD S-glass
compression cases) the ledger scores three predictors per case: the plain
Budiansky–Fleck angle floor, the closed-form **modulus** knockdown, and
the calibrated **penetration gate** — the UD strength path (issue #161)
that is sensitive to amplitude and through-thickness position
independently of the peak angle. On the S-M-2/4/5 trio (identical 20°
peak angle, amplitude 1.5/1.0/0.5 mm, measured KD 0.63/0.94/1.00 — a
~60 % strength spread invisible to any angle-only model) the gate lands
within **2.2 % / 0.6 % / 0.3 %**; over all six cases (including the
near-surface S-A-2 via the position factor) the mean absolute KD error
is **0.035**, with every case inside the ±20 % parity band. These
acceptance criteria are pinned as permanent regression tests in
[`tests/test_validation/test_ledger.py`](tests/test_validation/test_ledger.py).

### Quantitative validation against experiment

Comparison of the analytical predictions against published experimental
data is documented in the accompanying paper:

- Elhajjar, R. (2025). *Fat-tailed failure strength distributions and
  manufacturing defects in advanced composites.* Scientific Reports,
  15:25977. https://doi.org/10.1038/s41598-025-06693-4

Additional datasets referenced by the model calibration (Mukhopadhyay
et al., 2015; Li et al., 2026) are cited in [References](#references)
below. Mukhopadhyay's case data are in the validation ledger (Dataset C,
issue #462); Li et al. (2026)'s raw data are not included (tracking
issue: #22), so its case-level pass/fail table cannot be reproduced from
this repository.

For a consolidated predicted-vs-experimental view, the script
[`validation/plot_all_validation.py`](validation/plot_all_validation.py)
regenerates two parity figures against a ±20% band. The headline,
`validation/fig_all_validation_parity.png`, plots the analytical
prediction the package reports by default against each case's
ultimate strength, one point per case (Datasets A–F and H–M), with
each dataset predicted by the model that physically applies to it
(Budiansky–Fleck / three-mechanism for the multidirectional cases,
the penetration gate for the UD glass cases E/F). 61 of its 74 cases
fall inside the band and none is more than 20% unsafe. The 13 outside
it are all conservative: five of Dataset H's one-sided troughs (−22%
to −71%), four kink-band predictions on UD carbon (Datasets J and K,
−22% to −23%, with the UD strain budget read from the material card,
issue #451), three Dataset M tension coupons (−22% to −23%) and one of
Wang's wrinkles (−27%).

The second figure, `validation/fig_fe_onset_parity.png`, holds the
estimates that are not the default prediction or not of the ultimate
strength: first-ply FE retention on Datasets H to L, Dataset M's
first-ply FE onset against its measured damage onset, and
Mukhopadhyay's analytical onset estimate. 16 of its 31 points fall
inside the band. The 8 more than 20% unsafe are all first-ply FE onset
estimates on Dataset M, +23% to +154% against the measured damage
onset, as the FE warning and the guidance state.

No analytical point falls outside the band on the unsafe side. On
Dataset H's one-sided troughs, scored on the specimens' `tool_flat`
geometry (issue #461), the analytical path is conservative except +5%
on the deepest multidirectional trough, which it flags as
bending-dominated (run FE: on that half the FE retention, on the
specimen's own geometry, is within −8% to −1%). The most severe
whole-thickness waves (in Datasets A and J) carry an eccentric-column
bending bound and sit at +10% and +15% (issue #452). The earlier unsafe
points and the changes that removed them (issues #452, #454 and #461)
are recorded in [CHANGELOG.md](CHANGELOG.md) and
[VALIDATION.md](docs/internal/VALIDATION.md).

Dataset L, a low-quality wet lay-up laminate, is **quarantined**
(issue #445): it is plotted with hollow markers and left out of that
count. On its own, 8 of its 36 points fall inside the band; the
kink-band is conservative on all 18 cases, and the FE over-predicts its
embedded (mid-plane) wrinkles by up to +100%.

> **FE strength is indicative, not an allowable.** The analytical
> knockdown (the penetration gate for UD) is the validated strength path.
> The FE strength outputs — `retention_factors` and
> `progressive_knockdown` — were first checked on six UD glass/epoxy
> specimens, and there they missed on the unsafe side exactly where it
> matters most: FE LaRC05 over-predicted retained strength
> on the two most severe wrinkles, by up to +39% (while under-predicting
> the four milder ones), and the progressive-damage model over-predicted
> the most severe wrinkle by +31% (+40% on a refined mesh). On
> multidirectional laminates in compression FE LaRC05 has erred both
> ways: conservative (−9% to −1%) on Shi et al. 2025 and Thor et al.
> 2021 (Datasets H and J), but unsafe on Datasets A, C and D, by up to
> +92% (the head-to-head in
> `validation/strength_error_summary.py`, section 4). In tension the
> first-ply FE is an onset estimate: against Dataset M's measured damage
> onset it was unsafe on 9 of 12 coupons, by up to +154%. For an
> embedded wrinkle under flat surfaces no FE strength path is offered at
> all: the two size-aware candidates (progressive peak load,
> cohesive-zone delamination onset) stayed unsafe on every such dataset,
> by up to +109% and +139% (issue #453). The app and the NCR summary show this
> warning beside the FE numbers; the per-case figures are in
> [Interpreting results](docs/interpreting_results.md#how-far-to-trust-each-number)
> and `python validation/strength_error_summary.py` regenerates them.

### Stiffness (modulus) knockdown

Besides strength, WrinkleFE reports a **stiffness** knockdown of the axial
Young's modulus two ways: the FE `modulus_retention` (wrinkled vs pristine
from the linear static solve, any layup) and the closed-form
`analytical_modulus_knockdown` (a CLT series-average of the off-axis
lamina modulus over the wrinkle profile, any layup, no FE solve). The
script [`validation/validate_modulus.py`](validation/validate_modulus.py)
scores the FE value and the script's own CLT estimate against the UD
datasets that report a *measured modulus* — **F** (Li 2025, S-glass),
**G** (Hsiao & Daniel 1996, carbon — the `IM6G_3501_6` card), and the
indicative **E** (Li 2024): the script's estimate lands at 3.9 % MAE (F)
/ 1.2 % (G) and the FE at 7.2 % (F) / 5.7 % (G). The shipped
`analytical_modulus_knockdown`, pinned for every ledger case, is at
7.5 % MAE on F (see `docs/internal/VALIDATION.md`). The data
and both models agree that stiffness is far more wrinkle-tolerant than
strength: the modulus knockdown stays ≈0.81–0.98 for the S-glass cases
and only reaches ≈0.52–0.57 for a carbon uniform wrinkle at θ = 15°. The
script
[`validation/plot_modulus_validation.py`](validation/plot_modulus_validation.py)
renders the comparison as `validation/fig_modulus_validation.png` —
knockdown-vs-angle and a predicted-vs-experimental parity plot across all
three datasets.

## Supported morphologies

WrinkleFE ships six wrinkle morphologies (defined in
`src/wrinklefe/core/morphology.py`). They differ along *two* independent
axes: **how many wrinkles** are placed in the laminate and **how the
amplitude varies through the thickness**. The first three names below
are *dual-wrinkle* modes distinguished by the phase offset φ between
two adjacent wrinkle centrelines (the through-thickness amplitude
follows a linear taper from the wrinkle interface plies down to zero
at the laminate outer surfaces). The last three are *single-wrinkle*
modes that swap that taper for a different through-thickness profile.

| Morphology   | # wrinkles | Phase φ | Through-thickness amplitude            | M_f (compression) | When to use                                                                                                                                          |
|--------------|------------|---------|-----------------------------------------|-------------------|------------------------------------------------------------------------------------------------------------------------------------------------------|
| `stack`      | 2          | 0       | Linear decay, 1 at interface → 0 at surfaces | 1.0 (baseline)    | Two aligned wrinkles, peaks-over-peaks. The dual-wrinkle reference case used to scale `convex` / `concave`.                                          |
| `convex`     | 2          | +π/2    | Linear decay, 1 at interface → 0 at surfaces | < 1               | Two phase-shifted wrinkles whose interface bulges outward. *Least* damaging dual-wrinkle case in compression.                                       |
| `concave`    | 2          | −π/2    | Linear decay, 1 at interface → 0 at surfaces | > 1               | Two phase-shifted wrinkles whose interface pinches inward. *Most* damaging dual-wrinkle case in compression — design-driving.                       |
| `uniform`    | 1          | n/a     | Full amplitude on **every** ply (no decay)   | 1.0 (no pairing)  | A single through-thickness-wide wrinkle — every ply wavy with the same A. Conservative bound and sanity-check baseline.                              |
| `graded`     | 1          | n/a     | Linear decay from mid-ply to surfaces, with floor `decay_floor` ∈ [0, 1] | 1.0 (no pairing) | An embedded wrinkle that fades toward the surface plies. `decay_floor=0` is pure graded; `decay_floor=1` collapses to `uniform`.                    |
| `tool_flat`  | 1          | n/a     | Full-amplitude core, linear ramp over `surface_transition_plies`, **flat** at the pinned surface | 1.0 (no pairing) | A wrinkle cured against rigid tooling / a caul sheet: the pinned-flat surface fills the wrinkle troughs with **significant** resin pockets (auto-enabled). Also the morphology for a one-sided surface trough over a flat tool face (issue #461). Shares `uniform`'s kink-band knockdown but not its bending bound (it is not a whole-thickness wave); amplitude is bounded (see [The `tool_flat` morphology](#the-tool_flat-morphology-significant-pockets)). |

### `stack` vs `uniform` — what's the difference?

These two get conflated because both have `M_f = 1.0`, but they model
very different defects:

- **`stack`** places **two** wrinkles at adjacent interfaces with
  φ = 0 (aligned crests). Through the thickness the wrinkle decays
  linearly from the interface plies to zero at the outer surfaces —
  surface plies are flat.
- **`uniform`** places a **single** wrinkle and disables the
  through-thickness decay — every ply, including the outer surfaces,
  is displaced by the full profile.

For the same `amplitude` / `wavelength`, `apply_to_nodes` therefore
produces *different* deformed meshes: `stack` has a wrinkle
concentrated near the interface plies (and flat top/bottom plies),
while `uniform` has a wrinkle of the same amplitude at every single
ply. The `M_f = 1.0` coincidence is purely the analytical knockdown
parameter — the FE geometry, the per-ply fibre-angle field, and the
predicted ply-by-ply failure are not the same.

## Supported failure criteria

The criteria below live in `src/wrinklefe/failure/` and can be selected
through `FailureEvaluator` or used independently:

- **LaRC04/05** (`larc05.py`) — Pinho/Camanho 3-D criterion with
  fibre-kinking under compression (Xc-calibrated intrinsic misalignment and
  a kink-plane search, so a pristine ply kinks at exactly `-Xc`), in-situ
  matrix strengths, and a fracture-plane search. Default for the FE solve.
- **Tsai-Wu** (`tsai_wu.py`) — 3-D tensor-polynomial criterion with a
  configurable interaction coefficient.
- **Tsai-Hill** (`tsai_hill.py`) — 3-D extension of the classical
  quadratic Tsai-Hill index.
- **Hashin** (`hashin.py`) — 3-D Hashin criterion with separate
  fibre-tension/-compression and matrix-tension/-compression modes.
- **Puck** (`puck.py`) — action-plane (Mode A/B/C) inter-fibre-failure
  criterion with simplified fibre failure.
- **Maximum Stress** (`max_stress.py`) and **Maximum Strain**
  (`max_strain.py`) — non-interactive checks against the principal
  material-frame allowables.
- **Budiansky-Fleck kink-band** (`kinkband.py`) — analytical compression
  knockdown with an optional interlaminar damage coupling
  (`InterlaminarDamage`); this is the model exposed in the
  `analytical_knockdown` field of `AnalysisResults`.
- **Progressive damage** (`progressive.py`) — `PlyDiscount` and
  `ContinuumDamage` post-failure stiffness reduction models that wrap
  any of the criteria above.

## How It Works

The full mechanism-by-mechanism derivation lives on the
[Theory: physics & mechanics](https://wrinklefe.readthedocs.io/en/latest/theory.html)
page (`docs/theory.md`). The essentials:

### Wrinkle kinematics

A wrinkle is reduced to its peak fibre-misalignment angle

```
theta_max = arctan(2*pi*A / lambda)
```

(half-amplitude `A`, wavelength `λ`) and — for unidirectional laminates —
its through-thickness penetration `D/T = A/T` (`T` = laminate thickness).
`theta_eff = M_f * theta_max` folds in the morphology factor `M_f`
(`stack` = 1, `convex` < 1, `concave` > 1).

### Compression — CLT-weighted Budiansky–Fleck kink-band

```
KD_lam      = f_0 * KD_BF + (1 - f_0)
KD_BF       = 1 / (1 + r + c_AF * r^2),   r = theta_eff / gamma_Y_eff
gamma_Y_eff = max(0.032 + 0.050 * f_confined
                        - 0.010 * max(n_block_max - 1, 0),  0.016)
```

`f_0` is the axial-stiffness fraction carried by the 0° plies (the plies
that kink); the `(1 - f_0)` term is the off-axis plies riding through at
full strength. The matrix shear-yield strain `gamma_Y_eff` **rises** with
the confinement `f_confined` (off-axis neighbours bracing the 0° plies
against kink-band rotation) and **falls** with the longest run of
consecutive 0° plies `n_block_max` (blocked 0° plies kink more easily),
floored at half the UD value. So a dispersed `[0/45/90/-45]s` resists
wrinkle knockdown far better than a blocked `[0_4/90_4]s`. The optional
Argon–Fleck quadratic term `c_AF` (`kink_band_quadratic_coeff`) defaults
to `0` — the pure linear Budiansky–Fleck floor.

UD-dominated laminates (at least 75% 0° plies) are the exception on
every morphology except `graded` (issue #451). The confinement model
does not read the material and left UD carbon 58–78% conservative, so
these laminates take the strain budget from the material card instead:

```
gamma_Y_eff = g = S12 / Xc,   KD_BF = g / (g + theta_eff)   (c_AF = 0)
```

This is Budiansky–Fleck with an intrinsic misalignment, normalised by
the pristine `Xc`. The `graded` path keeps the confinement budget,
which was calibrated together with its profile average.

### Tension — three-mechanism minimum, CLT-weighted

```
KD_lam = f_0 * min(cos^2(theta), KD_matrix, KD_oop) + (1 - f_0)
```

The 0° ply knockdown is the *most severe* of three competing mechanisms:
fibre load-rotation (`cos^2(theta)`), in-situ matrix cracking (a
Hashin/LaRC `σ22`–`τ12` interaction with an in-situ shear strength that
is never below the lamina `S12`, `KD_matrix`), and a curved-beam
out-of-plane delamination check (`KD_oop`: the wrinkle curvature drives
an interlaminar `σ33` at the crest and `τ13` at the flanks, each bounded
by the interface strength, so `KD_oop >= 1/sqrt(2)`: past that strength
the interface delaminates and sheds the traction, and the ply fails later
by fibre or matrix). A Benzeggagh–Kenane mixed-mode delamination-*onset*
knockdown is reported alongside; it is unvalidated, and the output flags
it as such (issue #463). The tension knockdown is not floored by
the compression (kink-band) value (issue #454).

### Graded morphology

For the `graded` morphology the knockdown is averaged over the wrinkle
profile in both the longitudinal (`x`) and through-thickness (`z`)
directions. The **compression** path weights each ply by a Gaussian
through-thickness envelope centred at `wrinkle_z_position` (decay scale
`max(λ/2, A)`); the **tension** path uses an analogous linear taper,
applied to each ply's angle and to its amplitude (so its curvature and
interlaminar stress) alike. In
both, `decay_floor` sets the surface-ply amplitude: `0` is a fully
embedded wrinkle that fades to flat at the surfaces, `1` collapses to
`uniform`.

## References

- Elhajjar, R. (2025). Scientific Reports, 15:25977. https://doi.org/10.1038/s41598-025-06693-4 (Datasets A and B)
- Li, Y. et al. (2026). Composites Part A, 205:109719.
- Mukhopadhyay, S., Jones, M.I. & Hallett, S.R. (2015). Composites Part A, 73:132-142. https://doi.org/10.1016/j.compositesa.2015.03.012 (Dataset C, compression)
- Mukhopadhyay, S., Jones, M.I. & Hallett, S.R. (2015). Composites Part A, 77:219-228. https://doi.org/10.1016/j.compositesa.2015.07.007 (Dataset C, tension)
- Wang, Y., Zhu, S., Li, H., Zhou, L. & Yi, W. (2021). Polymers, 13(19):3204. https://doi.org/10.3390/polym13193204 (Dataset D)
- Li, X., Ge, J., Chen, G., Zhang, B. & Liang, J. (2024). Composites Science and Technology, 256:110762. https://doi.org/10.1016/j.compscitech.2024.110762 (Dataset E)
- Li, Y., Li, X., Ge, J. & Liang, J. (2025). Polymer Composites, 46(16):15176-15187. https://doi.org/10.1002/pc.30121 (Dataset F)
- Hsiao, H.M. & Daniel, I.M. (1996). Composites Science and Technology, 56(5), 581-593. (Dataset G)
- Shi, J., Yang, G., Sun, N., Zheng, J., Qian, J., Wang, W. & Song, K. (2025). Materials, 18:4503. https://doi.org/10.3390/ma18194503 (Dataset H)
- Luong, H.M., Trevarthen, J., Butler, R., Srisuriyachot, J. & Lunt, A.J.G. (2025). Composites Part B, 307:112818. https://doi.org/10.1016/j.compositesb.2025.112818 (review; knockdown context)
- Calvo, J.V., Quiñonero-Moya, A.R., Feito, N., Miguélez, M.H. et al. (2023). Composite Structures, 323:117498. https://doi.org/10.1016/j.compstruct.2023.117498 (Dataset I)
- Thor, M., Mandel, U., Nagler, M., Maier, F., Tauchner, J., Sause, M.G.R. & Hinterhölzl, R.M. (2021). International Journal of Material Forming, 14:19-37. https://doi.org/10.1007/s12289-020-01540-5 (Dataset J)
- Pilato, N.R., Butler, R., Trevarthen, J. & Lunt, A.J.G. (2022). Materials & Design, 222:111055. https://doi.org/10.1016/j.matdes.2022.111055 (Dataset K)
- Krishnappa, S. & Gururaja, S. (2024). Composites Part B, 284:111688. https://doi.org/10.1016/j.compositesb.2024.111688 (Dataset L)
- Elhajjar, R.F. & Shams, S.S. (2016). Composites Science and Technology, 122:140-148. https://doi.org/10.1016/j.compscitech.2015.11.026 (Dataset M)
- El-Hajjar, R.F. & Petersen, D.R. (2011). Composite Structures, 93(9):2400-2408. https://doi.org/10.1016/j.compstruct.2011.03.029 (Dataset M control strengths)
- Budiansky, B. & Fleck, N.A. (1993). J. Mech. Phys. Solids, 41(1), 183-211.
- Pinho, S.T. et al. (2005). NASA-TM-2005-213530.
- Camanho, P.P. et al. (2006). Composites Part A, 37(2), 165-176.
- Jin, L. et al. (2026). Thin-Walled Structures, 219:114237.

## License

MIT - see [LICENSE](LICENSE)

## Changelog

Notable changes between versions — including any that shift predictions,
flagged under a **Numerical results** heading — are recorded in
[CHANGELOG.md](CHANGELOG.md).

## Contributing

See [CONTRIBUTING.md](CONTRIBUTING.md). Measured strengths of wrinkled
laminates from published work are especially welcome: the
[**Contribute a validation dataset**](https://github.com/elhajjar1/wrinkleFE/issues/new?template=validation_dataset.yml)
form asks for what a validation-ledger entry needs, and every
contributed dataset is cited by name in the validation record. Bug
reports and usage questions go through the
[issue forms](https://github.com/elhajjar1/wrinkleFE/issues/new/choose).

## Citation

If you use WrinkleFE in your research, please cite it. The quickest way is the
**"Cite this repository"** button on the
[GitHub page](https://github.com/elhajjar1/wrinkleFE) — it's generated from
[`CITATION.cff`](CITATION.cff) and exports APA or BibTeX. The full software
citation:

> Elhajjar, R. (2026). WrinkleFE: An open-source finite element package for strength prediction of wrinkled composite laminates (Version 1.4.0) [Computer software]. Zenodo. https://doi.org/10.5281/zenodo.23105908

```bibtex
@software{elhajjar2025wrinklefe,
  author = {Elhajjar, Rani},
  title = {{WrinkleFE}: An Open-Source Finite Element Package for Strength
           Prediction of Wrinkled Composite Laminates},
  year = {2026},
  version = {1.4.0},
  publisher = {Zenodo},
  doi = {10.5281/zenodo.23105908},
  url = {https://github.com/elhajjar1/wrinkleFE},
  note = {University of Wisconsin-Milwaukee}
}
```

> **Software DOI:** [`10.5281/zenodo.23105908`](https://doi.org/10.5281/zenodo.23105908) is the
> *concept* DOI: it covers every version and always resolves to the latest
> release, so it is the right default. To pin the exact version you used,
> cite that release's own DOI from the
> [Zenodo record](https://doi.org/10.5281/zenodo.23105908) instead (1.2.0 is
> [`10.5281/zenodo.23105909`](https://doi.org/10.5281/zenodo.23105909)). The DOI badge for
> the *Scientific Reports* article is a separate citation, for the
> experimental data below.

Please also cite the underlying experimental validation data:

> Elhajjar, R. (2025). Fat-tailed failure strength distributions and manufacturing defects in advanced composites. *Scientific Reports*, 15, 25977. https://doi.org/10.1038/s41598-025-06693-4
