I was of course most interested in a paper with the above name [1] which uses an idea invented by a late colleague and myself, the Double-Hertz Model, [2] and relates closely to another of our contributions, An approximate JKR solution for an elliptical contact [3]. But that paper, as well as the numerical solutions by Wu [4], and the beautiful experimental confirmation of the predictions by Sümer at al [5], all agree that when adhesion is added to a Hertzian contact, the contact area, while remaining almost elliptical, is not the ellipse predicted by Hertz. It is always much less elliptical than the Hertzian ellipse, and its shape changes continually as the applied tension increases, becoming (modestly) more circular at pull-off. As an example, for \({R_1}/{R_2}=5\), Hertz gives \(b/a=0.346\), but at pull-off \(b/a=0.452\). Can predictions which ignore this behaviour be trusted? But more worryingly, can consideration only of the ends of the major axis be satisfactory?

In fact the authors [1] give an answer themselves. According to their Fig. 9, which gives results for almost the JKR region, the error, (using the Hertz geometry) appears to be \(0.5(1 - e)\) where \(e\) is the ellipticity \((b/a)\). So for the “roundest” ellipse studied by Johnson & Greenwood, where \(e \approx 0.4\), we may estimate an error of 30%, while for the highly elongated ellipses discussed by Zini et al. [1] in their introduction (but not studied by them) it could be of the order of 50%.

1 Why is this Analysis Wrong?

In the “approximate solution”, our aim was to obtain the same stress intensity factor all around the periphery of the contact. The “approximation” was our failure to achieve this, for our method of requiring equal values at the ends of the major and minor axes led to lower values everywhere else. But the paper by Zini et al has no such ambition, believing the behaviour at the ends of the major axis to be all that matters. There are, of course, examples in the fracture mechanics world where this is the case: with a pressurised elliptical crack no-one would look anywhere but at the ends of the major axis. But here we have a contact problem: the periphery is not determined in advance. To find it, we need the Barenblatt concept of an equilibrium SIF [6], where a crack will propagate for \(K>{K_{\text{c}}}\), and heal for \(K<{K_{\text{c}}}\). This concept becomes dominant in Schapery’s application to visco-elastic solids [7], where the crack-opening velocity increases (dramatically) when \((K - {K_{\text{c}}})\) increases, while the crack-closing (healing) velocity increases when \(({K_{\text{c}}} - K)\) increases. (See also [8]). For elastic solids, the effect is simpler: the crack tip either moves into the solid (\(K>{K_{\text{c}}}\)), or outwards (\(K<{K_{\text{c}}}\)). Thus in a study of contact, if \(K>{K_{\text{c}}}\), the contact edge will move inwards (as the gap between the two bodies extends), but if \(K<{K_{\text{c}}}\), the contact edge will move outwards. If these occur at the ends of the major and minor axes, respectively, the contact shape becomes less elliptical (Fig. 1).

Fig. 1
figure 1

Crack-opening and crack-closing at the boundary makes the contact less elliptical

2 Behaviour at the Ends of the Major Axis

To apply the Double-Hertz approach, Zini et al. necessarily assume that an approximate solution may be obtained within the framework of the basic Hertz theory: that the contact ellipse has the Hertzian ellipticity. [We shall use the ellipticity \(e=b/a\) in preference to the eccentricity \(k=\sqrt {1 - {e^2}}\) (but \(k\) is the argument of the Legendre elliptic integrals). The use of the Double-Hertz approach requires detailed analysis using elliptic integrals, and is difficult to follow. But a direct JKR solution for their limiting case \((\mu \;\;{\text{large)}}\)provides a useful warning of the consequence of their procedure of using the Hertzian geometry and examining only the ends of the major axis, while avoiding the detailed analysis.

The Hertz pressure distribution may be written \(p(x,y)=\frac{{{p_0}}}{a}\sqrt {{a^2} - {x^2} - {y^2}/{e^2}\;}\)with \({p_0}/a\) directly related to the macro-geometry, and so fixed: we shall write\(\frac{{{p_0}}}{{\sqrt {ab} }}=\frac{{E^{\prime}}}{M}\) i.e., \(\frac{{{p_0}}}{a}=\frac{{E^{\prime}}}{M}\sqrt e\). Accordingly the change in \(p(x,y)\) due to an increase in \(a\) (as is done in applying the Double Hertz method) will be \(\frac{{{p_0}}}{a}\frac{\partial }{{\partial a}}\left\{ {\sqrt {{a^2} - {x^2} - {y^2}/{e^2}\;} } \right\}=\frac{{{p_0}}}{a}\frac{1}{{\sqrt {1 - {x^2}/{a^2} - {y^2}/{b^2}\;} }}\).

Thus, the limit of the Double Hertz analysis is the standard jkr procedure of superposing a Boussinesq stamp pressure distribution on the non-adhesive pressure distribution.

Accordingly, we apply a Boussineq distribution \(\frac{{{p_1}}}{{\sqrt {1 - {x^2}/{a^2} - {y^2}/{b^2}\;} }}\)whose magnitude gives the desired stress intensity factor at \(x= \pm \,a\): this requires\({p_1}=\sqrt {2E^{\prime}\Delta \gamma /\pi \,a}\) so the load will be \(F=(2\pi /3){p_0}\;ab - 2\pi \,{p_1}\;ab=(2\pi /3)\frac{{E^{\prime}}}{M}\;{e^{3/2}}{a^3} - \sqrt {8\pi \,\Delta \gamma } \;e\,{a^{3/2}}\).

Writing this as a quadratic in\({a^{3/2}}\) it is readily found that the minimum is \(T= - {F_{\hbox{min} }}=3{e^{1/2}}M\Delta \gamma\). Hence the non-dimensional load \({T^ * } \equiv \frac{T}{{2\pi \sqrt {{R_1}{R_2}} \;\Delta \gamma }}\)will be \({T^ * }=\frac{3}{{2\pi }}{e^{1/2}}\frac{M}{{\sqrt {{R_1}{R_2}} }}\).

What is \(M\)? From Johnson Contact Mechanics (§ 4.26a,b) we have

\(\frac{1}{{2{R_1}}}=\frac{{{p_0}}}{{E^{\prime}}}\frac{b}{{{a^2}}}D(k);\quad \frac{1}{{2{R_2}}}=\frac{{{p_0}}}{{E^{\prime}}}\frac{1}{b}B(k)\) where \(k=\sqrt {1 - {e^2}}\).

[We have simplified Johnson’s equations by the use of Emde’s Elliptic Integrals \(B(k),\;\;D(k)\) where \(B(k)=\int_{0}^{{\pi /2}} {\frac{{{{\cos }^2}\varphi \;d\varphi }}{{\sqrt {1 - {k^2}{{\sin }^2}\varphi \;} }}}\); \(D(k)=\int_{0}^{{\pi /2}} {\frac{{{{\sin }^2}\varphi \;d\varphi }}{{\sqrt {1 - {k^2}{{\sin }^2}\varphi \;} }}}\) (see Jahnke & Emde Tables of Functions) …and what a disservice to the engineering world the compilers of the NBS Tables have done by ignoring these!]

Multiplying, \(\frac{1}{{4{R_1}{R_2}}}={\left[ {\frac{{{p_0}}}{{E^{\prime}}}} \right]^2}\frac{b}{{{a^2}}}\frac{1}{b}\;D(k) \cdot B(k)\) and \(\frac{{E^{\prime}}}{{\sqrt {{R_1}{R_2}} }} \cdot \frac{1}{{2\sqrt {D(k) \cdot B(k)} }}=\frac{{{p_0}}}{a}=\frac{{E^{\prime}}}{M}\sqrt e\) so that \(\frac{M}{{\sqrt {{R_1}{R_2}} }} \cdot =2\sqrt {e\;D(k) \cdot B(k)}\).

[For a circular contact \((e=1,\;k=0)\) we have \(B(0)=D(0)=\pi /4\) and \(M/R=\pi /2\), recovering the standard \(\frac{{{p_0}}}{a}=\frac{2}{\pi }\frac{{E^{\prime}}}{R}\)].

Table 1 shows the results, and the corresponding answers from Johnson & Greenwood’s “approximate elliptical JKR theory”. For comparison, I have calculated the values using the Elliptical Double-Hertz theory, again considering only the ends of the major axis (see Supplementary Material). Zini et al study e = 0.99 and e = 0.8, but it is not clear what their answers are.

Table 1 Comparison of this analysis with approximate jkr theory

3 Discussion

Predicting the pull-off force by considering only the behaviour at the ends of the major axis gives the wrong answer. But does this analysis predict the pull-off force? Indeed, does it predict anything? I believe that it does, and that the combination of the non-adhesive Hertz solution with the Boussinesq uniform (unloading) displacement is a true initial pressure distribution, and so correctly predicts the force at which the local SIF at the ends of the major axes exceeds\(K_{1}^{C}\). This need not lead to failure: it could well indicate the initiation of the continuous geometry change found by Johnson & Greenwood [3] in their “approximate jkr solution”, and found experimentally by Sümer et al. [5]. And indeed, two recent papers [9, 10] find that in numerical solutions of contacts of variously shaped flat punches, including an elliptical one [10], they do obtain continuous peeling, leading to a more circular contact area as the load increases to the much larger value at which pull-off occurs. Does the same occur with a paraboloidal indenter?

4 Conclusion

It is impossible to find the pull-off force of an elliptical contact by analysing only the behaviour at the ends of the major axis: presumably because at pull-off the SIF must be uniform all round the periphery.

There exists the possibility that such an analysis might correctly predict the much lower load at which local peeling begins.