3D Subsurface Stress Calculator

Use at your own risk. Verify inputs, assumptions and results independently before applying them.

3D Subsurface Stress

Inspect elastic half-space stresses beneath a Hertz point contact or a discretized normal and tangential traction map.

Analysis inputs

Material pair




Geometry and loading





Stress analysis

Maximum von Mises stress in section——
Contact radius—
Peak pressure—
Mutual approach—
Maximum shear——
Centerline stress profileEquivalent and maximum shear stress versus normalized depth

Off-axis von Mises sectionNormalized x/a and z/a; complete stress tensor

Stress tensor at section maximum MPa

σx—τxy—τxz—σy—τyz—σz—
Input definitions
Material elastic moduli E₁, E₂
Young’s moduli describing each body’s linear-elastic stiffness, in GPa; used together to calculate E*.
Material Poisson ratios ν₁, ν₂
Dimensionless lateral-to-axial strain ratios entering E*. The Hertz-mode stress field is evaluated inside body 1 using ν₁; ν₂ contributes to the combined contact modulus only.
Curvature radii R₁, R₂
Positive local spherical radii of the contacting bodies, in mm, combined by reciprocal addition. This mode assumes circular point contact.
Normal load W
Total positive compressive force carried by the contact, in N.
Tangential traction ratio μ
Dimensionless multiplier in qx = μp for fully sliding traction along x; zero gives frictionless loading.
Maximum plotted depth z/a
Dimensionless depth limit measured below the surface in multiples of the computed contact radius a.
Traction-map CSV
Each row supplies cell-centre x and y in mm, normal pressure and signed x/y shear tractions in MPa, and positive square-cell area in mm². Cell forces equal traction times area.
Half-space Poisson ratio ν
Dimensionless elastic Poisson ratio for the material receiving the specified traction map.
Section half-width
Positive lateral extent in mm: the plotted y = 0 section spans x from minus this value to plus this value.
Minimum and maximum evaluation depths
Positive lower and upper depths below the surface, in mm. Minimum must be below maximum; the surface itself is excluded because centroid point-load kernels are singular there.
How to use

Hertz point contact

  1. Select Hertz point contact.
  2. Enter both materials’ elastic moduli and Poisson ratios, the two positive radii, normal load, tangential traction ratio μ, and maximum plotted depth z/a.
  3. Select Calculate. Review contact radius, peak pressure, mutual approach, maximum von Mises and shear stresses, their locations, plots, and the stress tensor.
  4. Use μ = 0 for frictionless contact. Positive μ represents full sliding with qx = μp, not Cattaneo–Mindlin partial slip.

Traction map

  1. Select Traction map.
  2. Paste up to 2,500 square cells using the Traction map CSV columns x [mm], y [mm], normal [MPa], tangent_x [MPa], tangent_y [MPa], area [mm2].
  3. Use cell-centre x and y coordinates, non-negative normal traction, signed tangential tractions, and positive cell area.
  4. Enter Poisson ratio, section half-width, minimum positive depth, and maximum depth, then select Calculate.
  5. Check that the supplied cell tractions and areas represent the intended loads. Review the tangential resultant, maximum, contour and tensor. Repeat with a finer cell map around steep gradients to assess resolution dependence.

Reset restores the validated example values. Zero-depth evaluation is excluded because point-load kernels are singular at their centroids.

Worked example: frictionless Hertz contact

Enter 210 GPa and Poisson ratio 0.30 for both materials, 10 mm for both radii, normal load 1000 N, tangential traction ratio 0 and maximum depth 3 contact radii. The contact radius is approximately 0.3191 mm, peak pressure 4688.3 MPa and mutual approach 20.37 µm. The sampled section reports maximum von Mises stress about 2881.4 MPa and maximum shear 1440.7 MPa near x/a = 0 and z/a = 0.52.

The reported location is a sampled-section maximum, so resolution affects the last digits and location. These high elastic stresses are a mathematical benchmark, not evidence that a real material remains elastic; compare with suitable material limits before application.

Equations used

Reduced properties and circular Hertz contact

1/E* = (1 − ν1²)/E1 + (1 − ν2²)/E2
1/R* = 1/R1 + 1/R2
a = [3 W R* / (4 E*)]1/3
p0 = 3W/(2πa²),   δ = a²/R*,   p(r) = p0√(1 − r²/a²)

For ζ = z/a, compression negative:

σz/p0 = −1/(1 + ζ²)
σx/p0 = σy/p0 = −[(1 + ν)(1 − ζ atan(1/ζ)) − 1/(2(1 + ζ²))]

Discretized normal and tangential tractions

For cell force Pi = piAi, offsets x, y, positive depth z, and R = √(x² + y² + z²), the Boussinesq terms include:

σz = −3Piz³/(2πR⁵),   τxz = −3Pixz²/(2πR⁵),   τyz = −3Piyz²/(2πR⁵)

The complete Cartesian Boussinesq tensor is combined with differentiated Cerruti displacement fields for signed Qxi = qxiAi and Qyi = qyiAi. All cell tensors superpose linearly.

Stress measures

σv = √{[(σx−σy)² + (σy−σz)² + (σz−σx)² + 6(τxy² + τyz² + τxz²)]/2}
τmax = (σ1 − σ3)/2,   σh = (σx + σy + σz)/3

For traction maps, repeat the calculation with progressively finer cells and compare the stress field and peak location. A coarse centroid quadrature can misrepresent near-surface peaks; no absolute error bound is supplied.

Interpretation limits and references

The Hertz mode calculates circular point-contact geometry and uses complete Boussinesq–Cerruti elastic-half-space tensors for the off-axis section. The traction-map mode superposes user-entered normal and signed tangential tractions from square cells. Compression is negative.

Map evaluation is limited to 2,500 cells and a bounded 33 × 21 x-z section at y = 0. Evaluation begins at positive depth to avoid point-load singularities.

The model assumes a homogeneous isotropic linear-elastic half-space. It does not model yielding, coatings, adhesion, conformal geometry, thermal effects, or Cattaneo–Mindlin partial slip. Johnson, Contact Mechanics (1985) provides the Hertz and half-space formulation. Compare the related Hertzian line-contact stress calculator.