1
Introduction to Contact Problems in Structural Mechanics

Often, a body is not isolated in space. Whether motionless or moving, it is generally brought into contact with other bodies, over an essentially long period of time. In the case involving contact between several solid bodies (fluids are excluded here), physical observation reveals that the bodies cannot interpenetrate, and that their interaction takes place along their external boundaries. New mechanical connections are then formed; stresses develop on the contact surface and modify their shapes more or less strongly, as well as the trajectory of each body. Therefore, mechanics naturally include the study of contact with the aim of modeling this type of interaction between the bodies and calculating their resulting movement.

The first mathematical relation to express consists of interpreting the non-interpenetrability condition that can be observed, which appears as an inequality to be satisfied. Other relations are then to be added, the number and complexity of which vary according to the physical nature of the contact (with or without friction, with or without stick, static or dynamic, etc.). In general, a contact problem is always difficult to solve, insofar as it involves not only equations but also inequations.

In 1881, the frictionless contact problem between two spheres was solved analytically by Hertz, thus giving rise to the contact mechanics of structures. Hertz assumed that the contact portion is elliptical in shape, and that the bodies in contact behave in a linearly elastic way under small perturbations. Since then, other work has been carried out, allowing Hertz’s results to be extended to more complex configurations such as the solutions of Mindlin [94] and Spence [127, 128]. The reader can consult the works of Johnson [67] and Goldsmith [45], which provide a collection of analytical solutions of contact problems with or without friction, in statics or dynamics, and for elastic or elastoplastic bodies. More recently, Storakers [130] succeeded in solving Hertz’s problem with friction for various contact profiles by combining both analytical and numerical methods. In general, theoretical solutions of contact problems are difficult to construct; they are rare, and have often been obtained with simplifying assumptions, in particular that of small perturbations.

The development of computers since the 1970s has enabled researchers to solve more complicated and realistic contact problems numerically. Besides the typical nonlinearities – whether it be material or geometric – a numerical contact method must be able to deal with the non-linearity of the boundary interactions between two or more bodies. Difficulties stem from the fact that the contact laws reflecting these interactions are often irregular and multi-valued.

In the following, we present a brief review of the key numerical formulations that are used to solve contact problems. Some underlying mathematical results will be recalled succinctly and without proof. In order to understand how the numerical solving methods operate, we will rely on the following simple example of a mass-spring system:

EXAMPLE.– Consider a particle of mass m connected to a spring of stiffness k. In the reference state, there is no gravity, the spring is at rest, and the particle is at a distance h from a fixed plane (Figure 1.1(a)). Under the effect of gravity g, the particle moves by an amount u measured in the direction of gravity and calculated relative to the reference position of the particle (Figure 1.1(b)). The proximity between the particle and the rigid plane, denoted g(u) (no possible confusion with the acceleration of gravity g), is defined by g(u) ≡ u − h. This quantity must always be either negative or zero.

The problem therefore consists of determining the equilibrium position of the particle, that is to say the displacement u, under the force of gravity in the presence of an obstacle, which here is the rigid plane, and thus imposes the inequality constraint g(u) = u − h ≤ 0.

Figure 1.1. Mass-spring system. (a) Reference position (no gravity); (b) equilibrium position under gravity.

We easily obtain the exact solution of this problem:

1. If mg < kh, then the displacement is the particle is displaced without touching the rigid plane. The reaction from the rigid plane on the particle, denoted λ, is zero. The condition mg < kh signifies that the particle weight must be low enough so that no contact occurs.
2. If mg ≥ kh, then u = h: the particle enters into contact with the rigid plane. The contact reaction equals λ = mg − kh = mg − ku ≥ 0.

The total potential energy of the system is

In the h > mg/k case (Figure 1.2(a)), the solution u corresponds to the absolute minimum of Π. In the other case (Figure 1.2(b)), the solution u = h does not correspond to the absolute minimum of Π, but instead to the minimum of Π on the space of admissible solutions defined by g(u) = u − h ≤ 0.

The example above shows that the proximity g(u) and the contact reaction λ must satisfy the inequalities g ≤ 0 and λ ≥ 0, such that:

Figure 1.2. The displacement u which solves the mass-spring problem, represented by the thick black point

• if contact does not take place, that is to say if g(u) < 0, then λ = 0;
• if contact does take place, that is to say if g(u) = 0, then λ ≥ 0.

We can group these two cases together using the following mathematical relations, called the Hertz–Signorini–Moreau conditions:

Figure 1.3 shows the special relationship between g and λ which is imposed by the Hertz–Signorini–Moreau conditions.

Figure 1.3. Hertz–Signorini–Moreau conditions

In the general case of the contact problem between two 3D solids, it becomes difficult to solve it, due to the fact that the portions of the surfaces, which are likely to come into contact, are not a priori known, and can have complicated shapes, especially when the bodies undergo large deformations. The relations [1.2] remain valid, however, depending on a suitably defined proximity g and the normal contact reaction that we denote by λ, and must be written for any point of the surface at which contact may occur.

1.1. Solving a contact problem numerically via the penalty method

A classical method to solve the contact problem, and which is easy to implement, is the penalty method. This method can be presented within the framework of optimization theory, but it is easier to explain in terms of mechanics.

EXAMPLE.– Let us go back to the mass-spring system considered above, where we had to calculate the displacement u of the particle while respecting the inequality constraint g(u) = u − h ≤ 0. The idea of the penalty method is to abandon the constraint g ≤ 0 and, at the level of the rigid plane, to add an artificial spring of stiffness , in such a way that:

• if the particle does not make contact with the rigid plane, then this spring does not come into play (Figure 1.4(a));
• if the particle exceeds the rigid plane, then the artificial spring is compressed and opposes the particle with a resistive force proportional to the penetration g(u) = u − h (Figure 1.4(b)).

Thus, the artificial spring introduced makes it possible to oppose the penetration, and the larger its stiffness ϵ, the stronger the artificial spring. The condition g ≤ 0 is not completely satisfied, although intuitively, it should be expected that it should become more so as ϵ tends to infinity.

By employing the penalty method, we have to solve a classical problem found in structural mechanics, instead of having to deal with an inequality constraint.

Figure 1.4. Mass-spring system solution via the penalty method.

To solve the mass-spring problem, we therefore seek a u which minimizes the following function, and which is equal to the sum of the total potential energy [1.1] and the strain energy of the artificial spring:

where 〈g〉 stands for the positive part of 〈g〉: g = g if g ≥ 0, = 0 otherwise. In optimization theory, the additional term is called the penalty term, and expresses the fact that the artificial spring will only work when compressed. The ϵ term itself is called the penalty parameter.

The extremum of Πpen is obtained by setting the derivative with respect to u of this function equal to zero. Using (where H is the Heaviside function), we...