Hertzian contact problems

Two-dimensional plane-strain problem

The two-dimensional Hertzian contact problem validates the normal-contact algorithms in numgeo against the analytical Hertz solution. Both the element-based mortar (EBM) and segment-based mortar (SBM) discretisations are considered.¹

Figure 1 shows the mesh and the material properties—Young's modulus \(E\) and Poisson's ratio \(\nu\)—of the two deformable bodies. Each circular boundary has a radius of \(R=8\) m. A uniform vertical traction is applied to the upper body; integration over its diameter gives a resultant load of 10 kN per unit out-of-plane thickness. The nonmatching quadrilateral meshes deliberately use nonaligned surface nodes so that differences between the contact-discretisation methods become visible. At the symmetry axis, the surface edges of the lower body are twice as long as those of the upper body.

Figure 1: Model used for the two-dimensional plane-strain analyses.

Figure 2 compares the normal contact stress along the interface, measured from the symmetry axis, with the analytical Hertz solution. The quadrilateral models use EBM and SBM contact. Two additional EBM models use linear triangular elements with either conforming or nonconforming interface meshes. All four numerical solutions agree closely with the analytical pressure distribution; the largest local differences occur near the edge of the contact zone.

Figure 2: Normal contact stress for the two-dimensional plane-strain analyses.

The nonconforming triangular mesh is shown in Figure 3.

Figure 3: Nonconforming triangular interface mesh.

Three-dimensional problem with small deformations

The three-dimensional benchmark uses a spherical indenter and compares the full 3D solutions with an axisymmetric reference model. Only the EBM method is used because SBM contact is available only in two dimensions.

The upper sphere has a radius of \(R=1.25\) m and is displaced downward by 0.02 m into a much softer elastic body. The sphere is sufficiently stiff to behave almost as a rigid body. All simulations are static and geometrically linear. The benchmark includes the following meshes:

The corner nodes of the coarse u20-solid mesh coincide with those of the coarse U8-solid-3d-red mesh. Four representative 3D meshes are shown in Figure 4.

Figure 4: Representative meshes used for the three-dimensional analyses.

Figure 5 shows the normal contact stress along a radial line from the symmetry axis. With EBM contact, the contact constraints and tractions are evaluated at the surface integration points; the values written for visualisation are mapped to adjacent output locations. The plotted curves can therefore appear piecewise linear or locally nonsmooth and should not be interpreted as the interpolation used during assembly. Further details are provided in the contact-discretisation theory.

The axisymmetric and quadratic-hexahedral solutions agree closely. The linear 3D meshes show somewhat larger discretisation errors, particularly the coarse tetrahedral mesh. Refining the tetrahedral mesh substantially improves the agreement, demonstrating the expected mesh sensitivity of the contact-pressure distribution.

Figure 5: Normal contact stress for the three-dimensional small-deformation analyses.

Three-dimensional problem with large deformations

The preceding simulations are repeated with five times the indenter displacement, that is, 0.10 m. Geometric nonlinearity is enabled in this set of analyses.

Figure 6: Deformed configurations for the three-dimensional large-deformation analyses.

Figure 7 shows the normal contact stress along both in-plane symmetry directions. The unstructured 3D meshes are not exactly symmetric in the \(x\)–\(y\) plane, so the distributions along the two directions are not identical. The differences between element types and mesh densities are larger than in the small-deformation benchmark. The refined hexahedral and tetrahedral meshes reduce the directional differences, whereas the coarse tetrahedral mesh produces the largest deviation.

Figure 7: Normal contact stress along the two in-plane symmetry directions for the large-deformation analyses.

Comparison with other finite-element software

For context, the following official examples from other finite-element programs also study Hertzian contact. Their geometries, material parameters, meshes, loading procedures, and contact formulations differ from the numgeo benchmarks above, so the figures are illustrative rather than direct quantitative comparisons.

• Abaqus

Hertzian contact-pressure distribution reported in the Abaqus benchmark documentation.

• ANSYS

Hertzian contact-pressure distribution reported in the ANSYS example.

• Kratos Multiphysics

Hertzian contact result and displacement error reported in the Kratos Multiphysics example.

• COMSOL Multiphysics

Hertzian contact-pressure distribution reported in the COMSOL Multiphysics example.

• CalculiX

Hertzian contact-pressure distribution reported in the CalculiX example.