Beam / Slab on Elastic Foundation — MADEP CPT Interpreter

Scope and intent

1. Problem class and engineering purpose

The beam route treats a strip of finite bending stiffness EI_b resting on a soil foundation modelled as a distributed elastic reaction (Winkler) or as a spring + shear-layer pair (Pasternak). It is the appropriate screen for strip footings, narrow raft cuts, slab strips between movement joints, pipe bedding, and linear structures where bending along the strip governs the load-distribution pattern.

Transparency note. Support parameters are not read from a catalogue. They are reconstructed from the CPT-derived stiffness profile and are therefore footing-dependent, not universal soil constants. A different strip width or influence depth yields a different k_s.

Thickness sensitivity. The member thickness h is not only a section-capacity input. It changes EI_b, the characteristic length λ, and therefore the bending-moment envelope from the elastic-foundation solve. Increasing h normally reduces deflection, but it can increase M_Ed and sometimes A_s,req because the stiffer strip bridges load differently over the soil springs.

Coordinate convention. The model is one-dimensional. The x-axis is the analysed strip direction and is where L, patch start/end, point-load position, and the plotted response live. The y-axis is the model strip width b, perpendicular to x in plan. The z-axis is the vertical section depth h. For a transverse strip-footing check, x should run across the footing width and b is normally a one-meter slice along the wall. For a longitudinal beam check, x runs along the member. The app's Analysis direction selector changes these labels and the preview drawing; it does not change the underlying one-dimensional equations.

Primary quantities

w(x): Strip deflection [m], positive downwards.
M(x), V(x): Bending moment and shear for the modelled strip width; these are kN·m/m and kN/m only when b = 1.0 m.
k_s: Modulus of subgrade reaction [kN/m³ = kPa/m].
λ, β: Characteristic length and its inverse β = 1/λ [m, 1/m].

Derivation

2. Winkler equilibrium from strip kinematics

Winkler's hypothesis is that the reaction pressure at any point depends only on the local deflection: p_soil(x) = k_s w(x). Combined with Euler–Bernoulli beam kinematics (plane sections remain plane, rotations small) and equilibrium on an infinitesimal strip element of width b, the fourth-order ODE follows from sectional equilibrium on an infinitesimal slice.

M(x) = −EI_b w''(x), V(x) = M'(x) = −EI_b w'''(x)

Sectional equilibrium: dV/dx + q(x) − b k_s w(x) = 0

⇒ EI_b w''''(x) + b k_s w(x) = q(x) (Winkler ODE)

Non-dimensional form: w'''' + 4 β⁴ w = q(x)/EI_b, 4β⁴ = b k_s/EI_b

Notation

EI_b: Bending stiffness of the modelled strip width [kN·m²], or per meter only when b = 1.0 m.
q(x): Distributed applied load [kN/m].
b: Strip width used to convert modulus of subgrade reaction to a line reaction [m].

The operator w'''' + 4β⁴w is biharmonic plus reaction. Its Green function decays exponentially with a spatial scale set by 1/β, so loads are localized in their structural effect: a point load influences the strip essentially over ±3/β.

Support derivation

3. Modulus of subgrade reaction from the interpreted CPT profile

For each numerical sublayer beneath the foundation, the app evaluates the stress-dependent stiffness from the interpreted Hardening Soil style layer model. The sublayer values are averaged over the chosen influence depth, and the resulting equivalent soil modulus is converted to a subgrade reaction through the Vesić expression for an infinite strip.

E_oed,i = E_oed,ref[(c′ cot φ′ + σ′_v,i) / (c′ cot φ′ + p_ref)]^m

E_s,i = E_oed,i (1 + ν_s)(1 − 2ν_s) / (1 − ν_s)

E_s,avg = (Σ E_s,i Δz_i) / (Σ Δz_i)

k_s = [0.65 E_s,avg / (B (1 − ν_s²))] · (E_s,avg B⁴ / (E_b I_b))^1/12

Notation

E_oed,ref, m, p_ref: Hardening Soil reference oedometric stiffness [kPa], stress exponent [-], and reference pressure [kPa].
σ′_v,i: Effective vertical stress at the centre of sublayer i [kPa].
c′, φ′: Effective cohesion [kPa] and effective friction angle [°] of the governing soil.
E_s,i, E_s,avg: Plane-strain secant Young modulus in sublayer i and its thickness-weighted average [kPa].
k_s: Modulus of subgrade reaction applied to the strip [kN/m³].
B: Bearing/contact width used to derive k_s [m]. It is the Vesić support width, not automatically the same as the model strip width b.
ν_s: Poisson ratio adopted for the soil support conversion [-].
E_b I_b: Bending stiffness of the modelled structural strip. It is per meter only when b = 1.0 m.

The twelfth root of (E_sB⁴/E_bI_b) captures the Vesić soil–beam stiffness coupling: stiffer strips see less uniform contact pressure.
Because k_s depends on B, E_bI_b, and the chosen influence zone, two footings on the same soil can carry materially different k_s.
The averaging depth is exposed directly in the app. Larger values smooth the interpreted CPT stiffness over more soil; smaller values make k_s more sensitive to the layer just below founding level.

Characteristic length

4. Characteristic length λ and response regimes

Non-dimensionalising the ODE identifies a single intrinsic length scale λ = (4 EI_b / (b k_s))^1/4. Deflection, moment envelope, and contact-pressure distribution all scale with λ. The strip length L relative to this scale classifies the structural behaviour into short, intermediate, or long.

λ = (4 EI_b / (b k_s))^1/4, β = 1/λ

βL < π/4 (L ≲ 0.8 λ): short strip, nearly rigid body

π/4 ≤ βL ≤ π (L ≈ 0.8–3 λ): intermediate strip, elastic bending

βL > π (L ≳ 3 λ): long strip, response localized near loads

For short strips the moment envelope approaches the rigid-beam result; subgrade stiffness scatter matters less. For long strips a point load produces a decaying response beyond roughly ±π/β from the load — the strip isolates distant regions and far-field contact pressure becomes insensitive to local detail.

Closed-form solutions

5. Hetényi solutions for canonical load cases

For an infinitely long beam on Winkler support, the fourth-order ODE admits closed-form solutions in terms of damped oscillatory kernels (Hetényi 1946). The deflection under a concentrated load at the origin is the Green function of the operator; the response to a distributed load is the convolution of that Green function with the load.

w(x) = (P β / (2 b k_s)) e^−β|x|[cos(β|x|) + sin(β|x|)] (deflection under P at x=0)

M(x) = (P / (4 β)) e^−β|x|[cos(β|x|) − sin(β|x|)]

V(x) = −(P/2) e^−β|x| cos(β|x|) sgn(x)

w_max = P β / (2 b k_s) at x = 0

M_max = P / (4 β) at x = 0

Hetényi kernels

A(βx) = e^−βx(cos βx + sin βx): Deflection kernel.
B(βx) = e^−βx sin βx: Rotation kernel.
C(βx) = e^−βx(cos βx − sin βx): Moment kernel.
D(βx) = e^−βx cos βx: Shear kernel.

The kernels decay as e^−βx, so the response to a load localized near x = 0 is numerically negligible beyond about 3/β. That bound justifies truncating the numerical strip when inputs would otherwise require an effectively infinite domain.

Sensitivity to member depth

6. Why a stiffer strip can attract a larger bending moment

On a rigid support a thicker member is unambiguously beneficial: the moment is fixed by statics and the larger lever arm reduces stress. On a Winkler foundation that intuition fails. The moment is not fixed by statics; it is set by the competition between EI and k_s. As h grows the strip stiffens faster than the support, λ grows, and the strip bridges further across the soil column. The peak moment grows with it.

EI = E_b b h³ / 12 ⇒ EI ∝ h³

k_s(EI) ∝ (E_s B⁴ / EI)^1/12 ⇒ k_s ∝ h^−1/4 (Vesić soft coupling)

λ⁴ = 4 EI / (b k_s) ⇒ λ ∝ h^13/16 ≈ h^0.81

The Hetényi closed forms then dictate how M_max scales with λ for each canonical load case:

Concentrated load P (infinite strip): M_max = P/(4β) = P λ / 4 ⇒ M_max ∝ λ ∝ h^0.81

Patch UDL q over a span shorter than λ: M_max ∝ q λ² ∝ h^1.63

UDL covering the full free-ended strip: M_max → 0 (rigid translation w ≈ q/(b k_s))

The relation that controls the reinforcement output is then μ ∝ M_max/d² with d ≈ h. Carrying the exponents through:

Asymptotic h-scalings for the small-μ (linear elastic-RC) regime where ω ≈ μ.

Load case	M_max ∝ h^α	μ ∝ h^α−2	A_s,req ∝ ω·d ∝ h^α−1	Direction
Concentrated load	h^0.81	h^−1.19	h^−0.19	A_s falls slowly with h
Patch UDL (L_patch < λ)	h^1.63	h^−0.37	h^0.63	A_s rises with h
Full-strip UDL on free-ended strip	≈ 0	≈ 0	governed by A_s,min	insensitive to h

Exponents follow from λ ∝ h^13/16 together with the Hetényi closed forms in §5. The crossover is at α = 1: any load case where M_max grows faster than h¹ will increase A_s,req.

Worked illustration. Strip 1 m wide, b = 1 m, f_ck = 30 MPa, c_nom + φ_bar/2 = 36 mm. Patch UDL load case, M₁ = 100 kN·m at h₁ = 0.30 m. Going to h₂ = 0.50 m gives: M₂/M₁ = (0.5/0.3)^1.63 = 2.15, so M₂ ≈ 215 kN·m; d₁ = 264 mm, d₂ = 464 mm; μ₁ = 0.072 → A_s,1 ≈ 904 mm²/m, μ₂ = 0.050 → A_s,2 ≈ 1094 mm²/m. The moment growth out-runs the lever-arm gain by ~21%.

The crossover thickness at which adding depth stops helping depends on the load case and on the soil column. Engineers used to rigid-support beam intuition will find this surprising; the analytical reason is that statics alone no longer sets M.
The result is exact in the small-μ regime where ω ≈ μ. For very large μ (close to ductility limit μ_lim) the closed-form ω ≈ μ approximation no longer holds, but the qualitative direction is preserved.
If the ULS bending moment comes from an external structural model rather than the soil-supported solve here, the conventional intuition (more h → less A_s) does apply. The two checks are different problems and should be kept separate.
The same logic applies to k_s calibration: doubling k_s reduces λ by a factor 2^−1/4, lowers the patch-UDL moment by 2^−0.5 ≈ 0.71×, and lowers A_s,req by ≈ 0.71×.

Boundary and loading cases

7. Point, patch, and distributed loads with finite ends

For finite strips the app superposes particular solutions and applies end boundary conditions. Free-end strip footings carry M = V = 0 at the ends; fixed-end walls carry w = w′ = 0; hinged ends fix w and M. Patch loads of intensity q over [x₁, x₂] are handled by integrating the point-load Green function, giving closed-form end moments and deflections in terms of A, B, C, D kernels.

Load units. Uniform and patch inputs are line loads q(x) in kN/m along the x direction, not pressures. A pressure p in kPa becomes q = p b. In the transverse strip-footing interpretation, a wall line load N in kN/m along the wall spread over contact width t becomes q = N/t over the patch interval.

Uniform load q over [x₁, x₂], infinite strip:

w(x) = (q / (2 b k_s))[2 − A(β|x − x₁|) − A(β|x − x₂|)]

M(x) = (q / (4 β²))[C(β|x − x₁|) + C(β|x − x₂|)]

End-moment correction (finite strip): superpose fictitious end moments enforcing M = 0 (free) or w′ = 0 (fixed).

For the common "point load on infinite strip" case, the decay factor e^−β|x| gives a fast mental estimate of the affected zone: meaningful bending extends to ±π/(2β) on each side.
Patch-load sensitivity is strongest when the patch width approaches 1/β; very wide patches approach a rigid pressure distribution, very narrow patches approach the point-load limit.
End moments under uniform load on a finite strip can reverse sign relative to the infinite-beam field; the app reports the corrected envelope, not the superposition piece alone.

Pasternak extension

8. Pasternak shear-layer coupling

Winkler independence is the main limitation of the single-parameter model: adjacent springs do not communicate, so the predicted contact-pressure field is discontinuous at load edges. The Pasternak extension restores shear coupling through a distributed shear layer on top of the springs. The governing equation acquires an extra w'' term.

p_soil(x) = k_s w(x) − G_p w''(x)

EI_b w''''(x) − b G_p w''(x) + b k_s w(x) = q(x)

G_p = η G_s,avg H_p

G_s,avg = E_s,avg / [2(1 + ν_s)]

H_p = z_influence

Notation

G_p: Pasternak shear-layer stiffness per unit area [kN/m²]. Physically the integrated shear stiffness of the soil column above the Winkler base.
η: Calibration factor linking averaged shear stiffness to the Pasternak layer (default 1/3 for a linearly decaying shear profile).
G_s,avg, H_p: Averaged soil shear modulus [kPa] and adopted shear-layer thickness [m].

Physically G_p smooths the contact-pressure field and introduces a correction to the characteristic length. The Pasternak ODE has two nested length scales: one dominated by flexure and one by the shear layer. The app reports both so the engineer can see whether shear-layer coupling materially changes the moment envelope.

Kerr extension

9. Kerr three-parameter model (diagnostic)

Kerr (1964) places a second spring layer in series with the shear layer. The extra parameter matches observed wheel-load response on pavements and softens the contact-pressure kinks that Pasternak can leave near load edges. The app exposes Kerr as a comparator, not a default design mode.

p = k₁ w₁, (k₂(w − w₂)) = p − G_p w₂''

Consolidating gives a sixth-order ODE in w whose lowest-order terms recover Pasternak and whose highest-order term scales with G_p/k₁.

Kerr is supplied as a comparator for problems where Pasternak yields suspicious end moments or contact-pressure kinks.
Production reporting remains on Winkler / Pasternak.

Numerical discretization

10. Finite-element implementation of the strip

The strip is discretised with two-node cubic-Hermite beam elements. Each element carries two nodal rotations and two nodal deflections (4 DOF), and Winkler support is assembled as a consistent line spring. Pasternak shear coupling adds a consistent element contribution proportional to G_p ∫ N′_i N′_j dx.

K_e^beam = ∫ EI N''^T N'' dx

K_e^Winkler = ∫ b k_s N^T N dx

K_e^Pasternak = ∫ b G_p N'^T N' dx

K w = f, f_e = ∫ N^T q(x) dx + nodal actions

Element size defaults to min(λ/10, L/50) so the shortest retained wave of the Hetényi kernel is resolved by at least 10 elements per decay length.
Point loads map onto nodal actions exactly; patch loads use Gauss integration on the element span.
The linear system is symmetric positive definite for Winkler and for Pasternak when G_p ≥ 0.

Limitations

11. Assumptions and boundaries of validity

One-dimensional strip. Two-dimensional plate action (slab on column grid, orthotropic reinforcement) is outside this route.
Winkler and Pasternak are stiffness screens derived from the CPT column. They do not capture yield beneath local high-pressure zones — that is the domain of the Stage 2 deformation solver.
Contact uplift (tension under the footing) is not modelled in the public route. Strongly eccentric loads that would produce tensile reactions are flagged rather than auto-resolved.
Support stiffness is footing-dependent (B, E_bI_b, influence depth) — neighbouring strips on the same soil can carry different k_s.
Pasternak's η is a calibration. The engineer should compare Winkler and Pasternak runs before committing to the larger envelope.
The route is a preliminary sizing screen. A final structural design should carry the computed M, V, and contact-pressure fields into a dedicated structural model.

References

Hetényi, M. (1946). Beams on Elastic Foundation. University of Michigan Press. Origin of the A/B/C/D damped-kernel closed forms used in §5.
Vesić, A. B. (1961a). Bending of beams resting on isotropic elastic solid. J. Engng. Mech. Div. ASCE 87(EM2). Derivation of the twelfth-root k_s(B, EI) coupling in §3.
Vesić, A. B. (1961b). Beams on elastic subgrade and the Winkler hypothesis. Proc. 5th ICSMFE, Vol. 1.
Pasternak, P. L. (1954). On a new method of analysis of an elastic foundation by means of two foundation constants. Gosudarstvennoe Izdatel'stvo Literaturi po Stroitel'stvu i Arkhitekture, Moscow.
Kerr, A. D. (1964). Elastic and viscoelastic foundation models. J. Appl. Mech. 31(3), 491–498.
Selvadurai, A. P. S. (1979). Elastic Analysis of Soil–Foundation Interaction. Elsevier.