Formula reference

Every engineering formula the program computes — with variables, units, and where it lives in the source.

Units throughout: speed km/h, lengths m, rates as fractions (0.06 = 6%), K-values m per 1% of algebraic grade change.

Horizontal geometry

Curvature & bearing

κ = 1 / R · heading = (sin b, cos b) · right perp = (cos b, −sin b)

κ signed curvature (1/m); κ>0 turns right, κ<0 left, κ=0 tangent. b = bearing (deg, 0°=North, clockwise).

src/geom/include/yu3d/geom/HorizontalAlignment.hpp

Tangent: P(s) = start + s·(sin b, cos b); bearing constant.
Circular: constant κ = 1/R; bearing change Δb = s/R (rad).
Clothoid: κ(s) = κ₀ + (κ₁ − κ₀)·s/L (linear curvature) — see below.

Clothoid (Euler spiral) — derivation & integration

A clothoid is the transition whose curvature changes linearly with arc-length, so lateral acceleration ramps smoothly. yu3d stores it by its start point, start bearing, start/end curvatures (κ₀, κ₁, signed) and length L.

κ(s) = κ₀ + (κ₁ − κ₀) · s / L

Signed curvature at arc-length s along the spiral (0 ≤ s ≤ L). A “tangent end” is just κ = 0.

Heading is the integral of curvature (heading turns at the rate of the curvature), plus the start bearing θ₀:

θ(s) = ∫₀ˢ κ(t) dt + θ₀ = θ₀ + κ₀·s + (κ₁ − κ₀)·s² / (2L)

Heading varies quadratically in s (a closed form). This one has an exact expression — the position below does not.

Position integrates the unit heading vector. In the project's bearing convention (0° = North, heading = (sinθ, cosθ)):

x(s) = x₀ + ∫₀ˢ sin θ(t) dt · y(s) = y₀ + ∫₀ˢ cos θ(t) dt

Because θ(t) is quadratic, these are Fresnel-type integrals. They have no elementary closed form for general κ₀, κ₁ (only the special tangent-start case κ₀ = 0 reduces to the standard Fresnel integrals C(s), S(s)).

How the integrals are taken

yu3d evaluates x(s), y(s) by adaptive Gauss–Kronrod quadrature — Boost.Math gauss_kronrod<double, 15> (the G7–K15 rule: a 7-point Gauss estimate embedded in a 15-point Kronrod estimate), recursed to max_depth = 10.

The 15-point Kronrod value is the result; the difference from the embedded 7-point Gauss value is a cheap error estimate.
Where the estimate exceeds tolerance (a long or very tight spiral whose integrand oscillates), the interval is bisected and recursed (up to depth 10) — that's the “adaptive” part.
For typical road spirals (L ≈ 10–200 m, |Δκ| < 0.05) the 15-point rule is already far below drawing tolerance with no recursion.

Curvature at s is returned exactly from the linear law; only the position needs quadrature. This is why a clothoid is more expensive to sample than a tangent or arc (which have closed-form positions).

Sample at s: position via G7–K15 of (sinθ, cosθ) over [0, s]; bearing = θ(s); curvature = κ(s)

src/geom/src/HorizontalAlignment.cpp · HClothoid::theta_rad / HClothoid::evaluate (Boost gauss_kronrod)

Transition family (Bloss, cosine, sine, Helmert)

The clothoid is the linear member of a family: every transition ramps the curvature from κ₀ to κ₁ over the length L, and only the shape of the ramp differs. Write u = s/L (0…1) and let f(u) be the shape function:

κ(u) = κ₀ + (κ₁ − κ₀)·f(u)

f(0) = 0 and f(1) = 1 for every member, so all of them start and end at the required curvature.

src/geom/src/HorizontalAlignment.cpp

Law	f(u)	f′(0), f′(1)	Character
Clothoid	u	1, 1	curvature changes at a constant rate; a corner in the rate at each end
Bloss	3u² − 2u³	0, 0	rate starts and ends at zero — no jolt at the joints
Cosine	(1 − cos πu) / 2	0, 0	smooth like Bloss, slightly gentler in the middle
Sine	u − sin(2πu)/(2π)	0, 0	smoothest ends, steepest middle
Helmert	2u² (u ≤ ½) 1 − 2(1−u)² (u > ½)	0, 0	two parabolas meeting at mid-length; the classic German railway ramp

The invariant worth knowing: whatever the law, the total deflection over the transition is the same:

Δθ = L · (κ₀ + κ₁) / 2

The area under κ only depends on the endpoints, so swapping one law for another moves the curve's SHAPE but not where it ends up pointing. Which is exactly why a mis-mapped law is invisible at the joints and shows up in the middle.

Cubic parabola

The classical railway transition, defined in Cartesian form rather than by arc length:

y = x³ / (6·R·X)

X = the transition's extent along the entry tangent, R = the radius reached at its end. Curvature therefore grows with x, and the arc length is obtained by integrating √(1 + y′²); mapping a station back to x is a Newton inversion.

src/geom/src/HorizontalAlignment.cpp

Because the shape is an approximation by construction, the curvature at the end stays slightly below 1/R — that is the classical definition, and yucad3d keeps it rather than “fixing” a curve the standards define this way.

Minimum curve radius

R_min = V² / (127 · (e_max + f))

V design speed (km/h); e_max max superelevation; f side-friction factor (speed-dependent); 127 = unit/gravity constant.

src/geom/src/DesignStandard.cpp · min_radius()

Superelevation & runoff

Design rate by radius (AASHTO Method 2)

e(R) = e_max · R_min / |R|, clamped to [0, e_max]

e and side-friction both vary in proportion to 1/R between R_min (full e) and the no-superelevation radius. Sign banks into the curve.

src/geom/src/SuperelevationDesign.cpp · design_rate_magnitude()

Runoff length

L_r = w · |e| / Δ

w = width of rotated lanes = lane_width × rotated_lanes (m); |e| design rate; Δ max relative gradient from the standard (e.g. 0.005 = 0.5% at 80 km/h).

src/geom/src/SuperelevationDesign.cpp · compute_superelevation_runoff()

Runoff placement (no-spiral curve)

At a direct tangent↔curve boundary the runoff is split across the PC/PT:

entry: 0 at (B − f·L_r) → e at (B + (1−f)·L_r)
exit: e at (B − (1−f)·L_r) → 0 at (B + f·L_r)

B = boundary station; f = runoff_tangent_fraction ≈ 0.67 (~2/3 on the tangent). Where a clothoid exists instead, the cant ramps linearly along the spiral.

src/geom/src/SuperelevationDesign.cpp

Banking a point

dz = dz_raw + offset · e (axis-of-rotation points only)

offset = lateral distance from CL (m, + right); e = superelevation at the station. Non-AOR (auxiliary) points are not rotated.

src/geom/src/FlowGraph.cpp · emit_point()

Effective rate at a station

e(s) = interp over the diagram control points (clamped)

The diagram is the manual override if set, else the standard-aware auto runoff diagram when auto-super is on, else 0.

src/geom/src/Corridor.cpp · superelevation_at()

A simple legacy helper also exists — proportional to curvature, e = −κ × 6.0 clamped to ±8% (superelevation_rate()) — but the corridor's auto banking uses the Method-2 runoff diagram above.

Curve widening

w = n · (R − √(R² − L²)) + 0.1 · V / √R

n lanes; R radius (m); L wheelbase (m); V speed (km/h). First term = mechanical off-tracking, second = dynamic allowance. Capped at max_widening; ~0 above max_radius.

src/geom/src/CurveWidening.cpp

The widening tapers over a transition length and is interpolated per station, the same way superelevation runs off.

Vertical curves & sight distance

Parabolic vertical curve

z(x) = z_VPC + g_in·x + (g_out − g_in)/(2L) · x²

x = distance from VPC (0..L); z_VPC = e_VPI − g_in·(L/2); g grades (decimal); L curve length. Grade: g(x) = g_in + (g_out−g_in)·x/L.

src/geom/src/VerticalAlignment.cpp

Circular vertical curve (exact)

Some standards specify a circular arc rather than a parabola. yucad3d evaluates it exactly — no small-angle shortcut — which matters on steep grades where the two differ visibly.

φ = arctan g · τ = L / (cos φ_in + cos φ_out) · R = τ / tan(Δφ / 2)

φ = the grade as an ANGLE; τ = the tangent length from the VPI, equal on both sides. The horizontal span is therefore split ASYMMETRICALLY: a = τ·cos φ_in before the VPI, b = τ·cos φ_out after it.

src/geom/src/VerticalAlignment.cpp

sin φ(s) = sin φ_in + (s − s_VPC) / R · grade(s) = tan φ(s)

Walking the arc: the SINE of the slope angle changes linearly with distance along the curve — that is what makes it a circle, and what a parabola only approximates.

Not implemented on purpose: the vertical clothoid that IFC 4.3 also allows. Its exact z(s) definition is not published unambiguously, so importing one is refused with a reason rather than approximated by a parabola.

Algebraic grade change & K

A = (g_out − g_in) × 100 (percent) · K = L / |A|

A > 0 crest, A < 0 sag, |A| < 0.1% no curve needed. K = curve length per 1% grade change.

src/algo/src/DesignAudit.cpp

Stopping sight distance (SSD)

SSD = 0.278 · V · t + V² / (254 · (a/g ± grade))

V speed (km/h); t perception/reaction time (s); a deceleration (m/s²); g = 9.81; grade signed decimal (downgrade lengthens SSD). 0.278 = km/h→m/s; 254 = 3.6²/(2g).

src/geom/src/DesignStandard.cpp · stopping_sight_distance()

Crest & sag K-values

K_crest = SSD² / (200 · (√h_eye + √h_obj)²)

Sight line over a summit; h_eye driver eye height, h_obj object height (m).

src/geom/src/DesignStandard.cpp · crest_k()

K_sag = SSD² / (200 · (h_headlight + SSD · tan β))

Headlight criterion in a valley; h_headlight beam height (m), β beam up-angle (deg).

src/geom/src/DesignStandard.cpp · sag_k()

Railway cant

Equilibrium cant

E = G · V² / (127 · R)

E = cant in mm; G = dynamic gauge (mm, 1500 for standard gauge); V = speed (km/h); R = radius (m). The applied cant is E capped at the line's maximum; the shortfall (E − applied) is the cant DEFICIENCY a faster train feels.

src/modules/rail/src/Cant.cpp

Computed cant follows the curvature, so it ramps naturally through a transition and steps at a bare tangent–arc joint. Sign convention: positive = right-hand curve (the left rail is raised).

Stored cant profiles

A cant profile that came from a file (IFC, LandXML) is data, not a calculation, and yucad3d keeps it verbatim: segments of (start station, length, start/end cant, law). The laws are the same shape family as the horizontal transitions:

cant(u) = c₀ + (c₁ − c₀)·f(u)

f as in the transition table above: constant (f ≡ 0), linear, Bloss, cosine, sine, Helmert. Gaps between segments are unspecified rather than assumed to be zero.

src/geom/src/CantProfile.cpp

Where a stored profile exists it OUTRANKS the formula — the file's numbers are the engineer's decision. The Viennese bend is deliberately not implemented: its published definition is a 7th-order form we have not been able to verify, so importing one is refused rather than guessed.

Clearance envelope

The gauge outline is widened for the curve and tilted by the cant before anything is tested against it.

d_i = L_b² / (8R) d_a = (L_v² − L_b²) / (8R)

d_i = widening on the INSIDE of the curve (the bogie centres cut the corner); d_a = widening on the OUTSIDE (the ends overhang); L_b = bogie centre distance (m); L_v = vehicle length (m); R = |radius| (m). Both zero on a tangent, and both skipped for a gauge that names no reference vehicle.

src/geom/include/yu3d/geom/VehicleSweep.hpp — the one definition; applied numerically by the clearance sweep and as a station-varying expression by the platform cross-section

θ = asin(c / G) o′ = o·cosθ + h·sinθ h′ = −o·sinθ + h·cosθ

c = applied cant (mm, signed); G = rail head distance the cant is measured over (mm, 1500 for standard gauge); o, h = offset right of centre and height above top of rail (m). Positive cant raises the LEFT rail. The pivot is the track centre at rail level.

src/modules/rail/src/ClearanceGauge.cpp

These are the geometric terms every gauge standard builds on. The tolerance and quasi-static terms a specific norm adds belong in the gauge file, and the kinematic case — suspension travel, wear, vehicle roll — is not modelled at all: see Alignments for what a clean sweep does and does not prove.

Surfaces (TIN)

Delaunay triangulation

Incremental Bowyer–Watson insertion: for each new point, remove every triangle whose circumcircle contains it (the in-circle predicate), then re-triangulate the cavity. Breaklines are enforced by Lawson edge flips.

orient2d(a,b,c) = (b.x−a.x)(c.y−a.y) − (b.y−a.y)(c.x−a.x)

> 0 CCW, < 0 CW, = 0 collinear. The in-circle test lifts points to z = x²+y² and tests a 3D determinant.

src/geom/src/Surface.cpp

Elevation by barycentric interpolation

z = u·z_a + v·z_b + w·z_c, u+v+w = 1

u, v, w barycentric weights of (x,y) in the containing triangle (a,b,c). Outside the hull / inside a clip hole → no value.

src/geom/src/Surface.cpp · elevation_at()

Cut / fill volumes

cut = Σ(z_exist − z_design)·cell² (where exist > design)
fill = Σ(z_design − z_exist)·cell² (where design > exist)
net = fill − cut

A regular XY grid over the surfaces' overlap; elevations compared at cell centres (midpoint rule). Accuracy O(cell²). The overlap is rarely a whole number of cells, so the cells along the far edges are clipped to it and each contributes only the area it really covers — the reported volume is for the whole overlap, not for a grid that stops short of it.

src/algo/src/CutFill.cpp

Georeferencing & scale factors

combined = k_grid · k_height · k_height = R / (R + H)

k_grid projection distortion (from PROJ); R Earth radius (≈ 6 371 008 m); H elevation (m). grid = ground × combined.

src/core/include/yu3d/core/Geodesy.hpp

A line's combined factor is the mean of its endpoints' factors.

Design-standard constants

See Design standards for the per-standard table (e_max, cross slope, perception time, deceleration, sight heights). Source: src/geom/src/DesignStandard.cpp.

Side friction f (AASHTO, km/h → f)

30→0.28, 40→0.23, 50→0.19, 60→0.17, 70→0.15, 80→0.14, 90→0.13, 100→0.12, 110→0.11, 120→0.09, 130→0.08 (linearly interpolated).

Max relative gradient Δ (AASHTO, km/h → Δ)

30→0.75%, 40→0.70%, 50→0.65%, 60→0.60%, 70→0.55%, 80→0.50%, 90→0.47%, 100→0.44%, 110→0.41%, 120→0.38%, 130→0.35%.

Worked examples

Minimum radius @ 80 km/h (AASHTO)

R_min = 80² / (127 × (0.08 + 0.14)) = 6400 / 27.94 ≈ 229 m

e_max = 0.08, f(80) = 0.14. A curve sharper than ~229 m is flagged at 80 km/h.

Superelevation runoff length

L_r = w · e / Δ = (3.6 × 1) · 0.06 / 0.005 = 43.2 m

One 3.6 m lane rotated, design rate e = 6%, Δ(80) = 0.5%. ~2/3 (29 m) sits on the tangent before the PC, ~1/3 (14 m) into the curve.

Stopping sight distance @ 80 km/h (AASHTO, level)

SSD = 0.278·80·2.5 + 80²/(254·(3.4/9.81)) ≈ 55.6 + 91.0 ≈ 147 m

t = 2.5 s, a = 3.4 m/s². A downgrade lengthens it; an upgrade shortens it.

Crest K @ 80 km/h (AASHTO)

K_crest = 147² / (200 · (√1.08 + √0.60)²) ≈ 21610 / 638 ≈ 34

h_eye = 1.08 m, h_obj = 0.60 m. A crest curve joining a +2% / −2% grade break (A = 4) then needs L = K·|A| ≈ 136 m.

Keep this honest When a formula or constant changes in the source, update this page in the same commit (see help/README.md). The source paths above are the single source of truth; the worked numbers are illustrative and rounded.