LEAPS, Documentation
Theory, engine internals, validation cases, sign conventions, distress models and a usage guide for the Linear Elastic Analysis of Pavement Structures workbench. Everything here runs client-side; the numerical engine (LEAF-JS) is a from-scratch multilayer elastic solver.
1 · Overview
LEAPS computes the stress, strain and displacement fields inside a pavement idealized as a stack of homogeneous, isotropic, linear-elastic layers resting on a semi-infinite elastic halfspace. Each layer is defined by a thickness h, an elastic modulus E and a Poisson's ratio ν; the bottom layer is infinitely deep. Surface loads are applied as a uniform circular contact pressure (a tire, the classical layered elastic analysis idealization used throughout mechanistic-empirical pavement design), as a concentrated point force (Boussinesq), or as a line load spread along a straight segment. All three carry the same total force, so the switch changes only how that force is spread.
The tool is meant for intuition, teaching and rapid mechanistic checks, not as a replacement for a validated design package. Its core assumptions are:
- Each layer is homogeneous, isotropic and linearly elastic (no stress dependence, no visco-elasticity, no temperature gradient).
- Layers extend infinitely in the horizontal plane; the lowest layer is a halfspace.
- Loads are static and vertical, applied over a circular contact area, at a point, or along a line.
- Full continuity, frictionless slip or a linear shear-spring is enforced at each interface.
2 · Multilayer elastic theory
For an axisymmetric vertical load, the elasticity equations admit a solution through a single biharmonic Love stress function φ(r, z) that satisfies ∇⁴φ = 0. Taking the zeroth-order Hankel transform in the radial coordinate maps the biharmonic equation to an ordinary differential equation in the transform variable m, whose general solution in layer i is
The four constants Aᵢ … Dᵢ per layer are fixed by boundary and interface conditions. Stresses and displacements are recovered from φ̄ by differentiation in Hankel space and an inverse Hankel transform (an integral over m weighted by Bessel functions J₀, J₁).
Boundary conditions
- Surface (z = 0): prescribed normal pressure from the load transform, zero shear.
- Each interface: continuity of vertical stress and vertical displacement always; the treatment of shear stress and radial displacement depends on the bond (see §4).
- Depth (z → ∞): the growing exponentials in the halfspace must vanish, so
C_N = D_N = 0.
These assemble into a 4N − 2 linear system per transform parameter m for an N-layer system, solved independently for every quadrature abscissa.
e^(+m h) terms that overflow for thick layers or high m·h. LEAF-JS scales each layer's exponentials by its own e^(−m·z) reference and uses partial pivoting, so no modulus ratio or thickness, however extreme, can overflow the double-precision range.
3 · Numerical method
Each response is an infinite integral of an oscillatory, slowly-decaying integrand. LEAF-JS evaluates it with:
- Panel-wise Gauss-Legendre quadrature between consecutive zeros of the relevant Bessel function, so each panel integrates a smooth half-wave.
- Wynn ε-acceleration of the sequence of panel partial sums, which converts the slowly-converging alternating tail into a rapidly-convergent estimate.
- Decay-scale refinement for deep evaluation points, where the integrand's envelope decays fast and panels are re-spaced accordingly.
- Closed-form surface tails: at
z = 0the integrals reduce to Weber-Schafheitlin forms evaluated with complete elliptic integralsK,E(arithmetic-geometric mean). Surface deflection of a halfspace is therefore returned to machine precision rather than by truncated quadrature.
The solver tolerance (Solver panel) sets the convergence target of the accelerated tail; the contour resolution sets the size of the sampling grid used only for the stress-bulb visualization, not for the reported key responses.
What makes it fast enough to be interactive
A contour grid is 2,600 points and each point is a quadrature, so the whole picture is a few million evaluations of the layered kernels. Three things keep that inside a couple of seconds, and none of them costs accuracy.
- Shared quadrature nodes. The panel breakpoints come from the zeros of the Bessel families and the exponential decay scale, and both the decay length and the radial distance are quantized to powers of two. Neighboring points then land on exactly the same nodes and every 4N-2 coefficient solve is computed once and reused. Rounding the radius up makes the panels finer than the oscillation needs, so it cannot cost accuracy; measured on a 61x43 grid under a dual tandem it took the linear solves from 34,792 to 3,032.
- Half the grid, mirrored. When the gear is its own mirror image in x, which is every layout the gear generator builds, only half the columns are solved and the other half is a reflection: the stresses and strains carrying one x index change sign, everything else does not, and the principal values are untouched because a reflection is orthogonal. Exact, not approximate, and it halves the most expensive job in the tool.
- A coarse pass first. Where the full grid is more than a moment away, a grid at about a third of the resolution is solved first and drawn while the full one runs. The picture arrives in about a second and sharpens rather than appearing all at once two seconds later.
The contour is also not solved at all while the overlay is switched off, and the status bar names what each pass cost, so none of this has to be taken on trust.
4 · Interfaces
Between adjacent layers, vertical stress and vertical displacement are always continuous. The shear behavior is set by a dimensionless slip value, the same input WinJULEA takes, where 0 is fully bonded and 1 is frictionless. Both endpoints are exact. Anything between is Goodman's shear spring:
| Interface | Shear condition | Typical use |
|---|---|---|
| Fully bonded | Shear stress and radial displacement continuous | Well-constructed flexible layers, tack-coated lifts |
| Frictionless | Zero shear; radial displacement free to slip | Unbonded PCC/base, bond-breaker, debonded overlay |
| Shear spring | Shear stress = k · (radial slip); k in MPa/mm | Partial slip, aged tack, transition between the two limits |
The spring stiffness k spans the two limits continuously: k → ∞ recovers the fully-bonded solution, k → 0 recovers the frictionless one. Both limits are checked in the validation suite (§9).
Turning a dimensionless slip into a stiffness necessarily introduces a length. The literature, BISAR's shear-spring compliance, Uzan's friction parameter, scales it on the loaded radius, which is why a dimensionless friction parameter is documented as depending on the diameter of the load circle. LEAPS uses
k(s) = (Glower / aref) · (1 − s) / swith aref the reference contact radius, and prints the resulting stiffness beside the control. An intermediate value is a defensible partial bond, but it is not guaranteed to reproduce another program's intermediate value digit for digit, because each program normalizes its own way. Enter k directly when the number matters.
5 · Loads & gears
The load model selector in the toolbar chooses how each load is applied. Every model carries the same total force F, so switching between them is a direct experiment in Saint-Venant's principle: near the contact the three answers differ by a factor of two or more, and by a couple of contact radii down they agree to a fraction of a per cent.
| Model | What is applied | Where it is singular |
|---|---|---|
| Circular imprint | Uniform pressure p over radius a. The LEA standard, and the only model to use when comparing against another layered-elastic program. | Nowhere, finite everywhere. |
| Point load | A concentrated force F (Boussinesq). | At the load itself. |
| Line load | Force per unit length q = F/L along a segment of length L and azimuth θ. | Anywhere on the segment, at the surface. |
A singular point is reported as such rather than given the finite number a quadrature happens to return there. That the concentrated idealizations have no answer at the load is not a defect in the solver, it is the reason design methods spread a wheel over a contact area in the first place.
Load, pressure and area are one equation
A circular contact is described by three numbers and they are not independent:
F = p · A, A = π a², a = √( F / (π · p) )Any two fix the third, and which two you have depends on where the number came from: a legal axle load and a placarded inflation pressure, or a measured footprint and a scale reading. The Solve for switch names the one LEAPS computes. That quantity becomes a read-only readout rather than an input, so there is no state in which three numbers can sit on screen that do not satisfy the relation, and no question about which of them the solver believed. Contact radius is offered alongside area because a footprint is quoted both ways.
The same switch governs the gear generator, so a gear built from load and pressure and a gear built from load and area are the same gear.
Multiple wheels are superposed: LEAPS solves each wheel independently in its own axisymmetric frame, then rotates every stress/strain tensor into the global frame before summation, so shear coupling between offset wheels is captured exactly. The Gear generator builds standard layouts, all parameterized by wheel load, pressure, dual spacing Sd (transverse gap between paired tires) and tandem spacing St (longitudinal gap between axles):
| Gear | Wheels | Layout |
|---|---|---|
| Single | 1 | One wheel at the origin |
| Dual | 2 | Pair at ±Sd/2 |
| Dual tandem | 4 | Two duals, axles Sd apart at ±St/2 |
| Tridem | 3 | Three single wheels, axles St apart |
| Dual tridem | 6 | Three duals (±Sd/2) on axles at −St, 0, +St |
Individual wheels can be repositioned and re-loaded in the wheel table after a gear is built, so irregular configurations are fully supported.
6 · Conventions & units
- Stresses are tension-positive. Compression under a wheel is therefore negative.
- z is measured downward from the surface; deflection w is positive downward.
- The section is a vertical cut in the x-z plane at the plan y-coordinate set by the draggable section line.
- Shear strains are engineering strains, γ = τ/G.
- Internally everything is
mm · N · MPa, and SI display is that same system, so every SI conversion factor is exactly one. The English system is the column set WinJULEA prints: inches, pounds and psi, with a modulus in psi rather than ksi and a displacement in inches rather than mils, so a result reads straight across from the other program. - LEAPS asks which system you want once, before anything is typed. Every number in the workspace is in it, and an analysis read in the wrong system is not wrong by a factor anyone notices in a plot. The choice is remembered and stays changeable from the toolbar.
| Quantity | SI display | US display |
|---|---|---|
| Length / thickness | mm | in |
| Stress / pressure | MPa | psi |
| Modulus | MPa | psi |
| Wheel load | N | lb |
| Deflection | mm | in |
| Strain | µε | µε |
| Interface stiffness | MPa/mm | psi/in |
7 · Response quantities
At every evaluation point LEAPS reports the full state:
- Stress tensor: σx, σy, σz and the shear components τxy, τxz, τyz.
- Principal stresses σ1 ≥ σ2 ≥ σ3, the maximum shear τmax = (σ1 − σ3)/2, and the von Mises stress.
- Strain tensor: εx, εy, εz (via the 3-D elastic compliance of the local layer).
- Displacement: ux, uy and the vertical deflection w = uz.
The Key responses rail highlights the classic design drivers: maximum surface deflection, horizontal tensile strain at the bottom of the asphalt (fatigue), vertical compressive strain at the top of the subgrade (rutting), tensile stress at the bottom of bound layers, and peak interface shear. Those cards carry a tint of their own so they read as conclusions rather than as more of the table, and the one that governs is outlined.
Reading the results table
The table is transposed the way a layered-elastic program prints it: one column per evaluation point, one row per quantity. Every cell also carries a bar scaled to the largest magnitude in its own row and tinted by sign, and the largest cell of each row is marked. A row of numbers in scientific notation is hard to scan and easy to misread by a decade; the bar is there so the shape of a response along the section can be seen before any single value is read. Because each row is scaled to itself, bars are comparable across a row and never across rows.
A response that is zero by symmetry comes back from the quadrature as roundoff rather than as zero: on a wheel axis the shear is 2e-18 where the answer is exactly nought. Anything twelve orders below the largest value in its own row is printed as 0, in the table and in the exported analysis alike. Twelve, not six: the engine agrees with an independent solver to about five significant figures, so a number that small is not a small response, it is the last bits of the double. The project file keeps the raw value.
Greek letters, subscripts and superscripts are typeset rather than spelled: the table, the field picker, the colorbar and the canvas all draw from one notation table, so τxz is the same glyph everywhere it appears and the CSV twin is a plain-text transliteration of the same entry. Real equations are set with KaTeX, each carrying a plain-text twin that shows until it loads.
8 · Distress models
The Performance tab converts the two critical mechanistic strains into allowable load repetitions using the Asphalt Institute transfer functions. Fatigue cracking is driven by the horizontal tensile strain εt at the bottom of the asphalt and the asphalt modulus E (in psi):
Subgrade rutting is driven by the vertical compressive strain εv at the top of the subgrade:
The governing life is the smaller of the two. N_f and N_r are allowable repetitions of the analyzed load, not calendar life, no seasonal, reliability or shift factors are applied.
9 · Validation suite
LEAF-JS ships with a 112-check regression suite (run under Node during development; a subset runs in your browser at every startup and drives the status-bar badge). The checks cover:
| Category | What is verified |
|---|---|
| Boussinesq closed forms | σz, σr, radial/vertical displacement of a point and circular load on a homogeneous halfspace, at the surface and at depth. |
| Multilayer reduction | An N-layer stack with equal moduli reproduces the single-halfspace solution. |
| Interface limits | Shear-spring stiffness k → ∞ matches the fully-bonded solution; k → 0 matches the frictionless one. |
| Vertical equilibrium | Integrated vertical stress on horizontal planes balances the applied load. |
| Superposition symmetry | Symmetric multi-wheel gears yield symmetric response fields. |
| Burmister two-layer | Surface deflection factor F₂ matches Burmister's tabulated two-layer solution. |
| Extreme contrast | A 1000:1 modulus contrast with a 3 mm interlayer solves stably without overflow. |
| Load idealizations | A point load is the limit of a shrinking circle at constant total force; a line load of vanishing length is a point load; all three carry the applied load through any horizontal plane; a rotated line load gives a rotated field. |
| Slip values | Slip 0 reproduces the bonded solution exactly and slip 1 the frictionless one; deflection rises monotonically between them. |
| Principal values | σ₁ ≥ σ₂ ≥ σ₃ everywhere, the eigenvalues carry the tensor's own invariants, and the principal strains are Hooke's law on the principal stresses. |
Two quadrature defects were found by cross-checking the engine against an independent implementation, and both were converged and wrong, stable at every tolerance, so no amount of tightening would have shown them. The panel mesh near m = 0 now resolves the interface reflections rather than the evaluation depth, and the first panel is graded geometrically toward zero. A surface deflection that came out 0.8 % high on a three-layer section now agrees with the reference to 4e−5. The lesson is recorded here because it generalizes: a single solver cannot check its own quadrature.
The startup self-check compares the engine against the Boussinesq closed forms; a green badge means those passed on your machine, in your browser, this session.
10 · Using LEAPS
- Start: Choose units and a template. Set layer thicknesses, moduli and interfaces. The bottom layer is a half-space.
- Loads: Choose a load model. For circular contacts, enter two of force, pressure and area. Press Build after editing gear settings.
- Points: Add evaluation points in the panel or double-click the section.
- Run: Resolve input warnings, then review the table, contours and key responses.
- Study: Sweep a parameter in Design study. Set a target life to find a crossing within the sampled range.
Tension is positive; depth and vertical displacement are positive downward. Slip 0 is bonded, 1 is frictionless. Point and line loads are singular at the loaded surface.
Views
In Section, drag an interface to change thickness, click its chip to change bonding, or click a layer to select it. In 3D, drag to orbit, scroll to zoom and shift-drag to pan. Click the compass to restore isometric view. Drag the cut plane or the plan inset's dashed line to move the section.
Keyboard & interaction
| Key / action | Effect |
|---|---|
F | Fit the section view |
R | Run the analysis |
Ctrl/⌘ + Z / Y | Undo / redo |
Delete | Remove the selected evaluation point |
| Scroll / drag on section | Zoom / pan the viewport |
| Double-click section | Drop an evaluation point |
| Drag an interface line | Set the thickness of the layer above it |
| Click an interface chip | Cycle bonded, partial, free |
| Click a layer | Select it here and in the structure panel |
| Drag in the 3D view | Orbit the camera (shift-drag pans) |
| Drag the cut plane in 3D | Move the section line through the gear |
Exports
Export full analysis saves inputs and responses in the displayed units. Save a .leaps.json project to reopen the model. Export contours or studies as CSV and the section as PNG. Work also autosaves in this browser.
11 · References
- Burmister, D. M. (1945). The general theory of stresses and displacements in layered soil systems. Journal of Applied Physics, 16.
- Huang, Y. H. (2004). Pavement Analysis and Design (2nd ed.). Pearson.
- The Asphalt Institute (1982). Research and Development of the Asphalt Institute's Thickness Design Manual (MS-1), Research Report 82-2.
- ARA / NCHRP (2004). Guide for Mechanistic-Empirical Design of New and Rehabilitated Pavement Structures (MEPDG), NCHRP 1-37A.
- FAA (2021). Airport Pavement Design and Evaluation, AC 150/5320-6G, and the LEAF layered elastic engine used in FAARFIELD.
- Watson, G. N. (1944). A Treatise on the Theory of Bessel Functions. Cambridge University Press (Weber-Schafheitlin integrals).