Open Access
 Issue Mechanics & Industry Volume 19, Number 4, 2018 405 14 https://doi.org/10.1051/meca/2018033 01 November 2018

## 1 Introduction

The dissipative mechanical energy due to the frictional work during the sliding contact process is transformed into frictional heating, which leads to temperature rise in the contact bodies. This thermal phenomenon exists widely in engineering applications and can significantly affect the mechanical properties of contact bodies. The knowledge of the thermoelastic stress field in sliding solid bodies is essential for accurate fatigue life and strength analysis because the total stress is comprised of the thermal stress and elastic stress. The thermal and mechanical analyses should be fully coupled and conducted simultaneously in sliding contact.

An overview of thermoelastic mechanics has been summarized in Johnson's book [1] and Yan's book [2]. The instantaneous point heat source solution and the heat source method proposed by Jaeger [3] are commonly used to solve various complicated heat transfer problems in tribology and manufacturing processes [4,5]. According to the heat source method, Hou [5] developed the comprehensive solutions for stationary/moving plane heat sources of various shaped. Levin [6] suggested that Jaeger's solution for the quasi-steady-state problem of moving heat source is inconsistent with the character of source intensity, he gave the improved fundamental solutions of various heat source (such as point, line and circular) by exponential transformation of moving coordinates, and the comparison between Jaeger's solution and fundamental solution was conducted to explain the issue he pointed out. Bosman [7] developed a semi-analytical transient temperature model by discretizing arbitrary shaped heat source into uniform square heat source, the method adopted here actually applied the idea of calculus. This solution is appropriate for microcontact and macrocontact. However, it can be only used to solve the surface temperature of bodies and cannot be used to analyze the temperature below the surface. Even so, his work is quite a useful reference in terms of temperature analysis. The separation of variables is also an effective feasible method for exploring the analytical solution of heat conduction equation [810]. Kuo et al. [9] investigated the 3-D heat conduction problem with a rectangular-shaped moving plane heat source by applying the method of separation of variables, the temperature rise expressions derived are similar to the form of solutions obtained through the heat source method. Nevertheless, the solutions of moving heat source presented by Araga et al. [9] and Belghazi et al. [10] are not convenient to use due to the complicated construction. Hamraoui et al. [11] presented an exact analytical solution of temperature in the form of series, but each term of expression is computed impossibly in actual application, as a result, his solution is still approximate. In addition, the solution is valid for the case of uniform heat source, steady- state and cylindrical coordinate, which leads to a limited application scope.

Kennedy et al. [12], a brief review of the earlier work in terms of contact temperature, pointed out that although analytical solutions of contact temperature are readily available, most of them do not consider the potentially important factors such as surface roughness, convective cooling or finite body geometry. To select the simple, accurate, and easy-to-use expressions that can be used to predict contact temperatures, Kennedy et al. [13] considered a variety of models for temperature rise, and the predictions of those models were compared with finite element thermal models of sliding contacts, the evaluations show that the models suggested by Alilat et al. [14] and Laraqi et al. [15] gave a more accurate representation of nominal steady-state temperature rise. Bansal et al. [16] studied the accuracy of approaches presented by Blok [17], Jaeger [3] and Tian and Kennedy [18] for the estimation of maximum and average temperatures at the interface for the case of sliding contacts, it was found that the accuracy of formulas developed by Tian and Kennedy is the highest, and the typical errors of the method presented by Blok less than 3%, whereas Jaeger's hypothesis leads to the accuracy that has typical errors of around 6%. Based on this research, Tian and Kennedy's formulas are the first choice, followed by Blok's method, and Jaeger's method at last for calculating the maximum and average temperature rise. Wu et al. [19] used the formulas developed by Tian and Kennedy when analyzing the characteristics of contact temperature under three-body dry friction conditions.

In addition to the analytical and semi-analytical models mentioned above, the numerical method is prevalent method for studying the contact temperatures recently, of which the finite element method is the most widely used, such as the thermal analysis of automotive brake, gear drives, etc. [2023]. However, the foregoing literatures involving thermal do not consider the effects of thermomechanical coupling. Research manifested that the thermomechanical effect plays a significant role in the variation of contact temperature and pressure in some cases [24].

In view of the mathematical complexity in thermomechanical problems, the finite element method has been widely used to investigate the contacts which involve frictional heating [2528]. The finite element thermomechanical models developed by Hwang [25] and Balci [26] are a one-way coupling, which only consider the effects of thermal on mechanical, and Balci's model is suitable for 2-D thermomechanical analysis. Comparing to the numerical methods, analytical methods have contributed to understand the essence of thermomechanical coupled problems if the analytical solutions of problems could be obtained. Kennedy [29] reviewed the significant developments in the research of thermomechanical phenomena. Liu and Wang [30] investigated the 3-D thermomechanical contact of non-conforming rough surfaces and derived the analytical solutions using the integral forms of known Green's functions; however, the solutions presented are adjusted to small Peclet numbers and steady-state heat transfer. They also illustrated and discussed the thermoelastic stress field inside a translating half-space with constant velocities by using the discrete convolution and fast Fourier transform method when a parabolic type or an irregularly distributed heat source is applied [31]. Chen and Wang [32] developed a three-dimensional thermoelastoplastic contact model for counterformal bodies, which considers steady state heat flux, temperature-dependent strain hardening behavior, and interaction of mechanical and thermal loads. Wang and Ai [33] established an extended precise integration model to analyze the thermomechanical behavior of multilayered transversely isotropic materials, and the actual solutions were obtained by adopting numerical quadrature schemes for the Laplace-Fourier transform inversion. Altan et al. [34] studied the elastic-plastic thermal stresses of the composite disc under known temperature rise by using the Fortran computer program. Zhang and Xiang [35] proposed a semi-analytical method for calculating the displacement and stress induced by heat transfer using Goodier's thermoelastic displacement potential and Laplace transform, but the thermal response was not influenced by the mechanical response. Kulchytsky-Zhyhailo et al. [36] studied the thermal stresses for two extreme cases of an axially symmetric contact problem with frictional heat, i.e. the unloaded bounding surfaces of bodies are thermally insulated and the temperature of outside of contact area is zero. It was found that the thermal stresses in both cases differ considerably.

Although the previous works have contained very useful information and insight into the temperature and stress fields in solids due to the stationary/moving heat source or the sliding friction, the comprehensive study of the semi-analytical transient solutions of fully coupled thermomechanical problem for the sliding contact are inadequate. Therefore, the main objective of present research was to deduce the systematic semi-analytical transient solutions of fully coupled temperature and stress fields for the sliding contact based on the heat source method and the thermal-elastic displacement potential. Results of temperature rise because of different variables and the effects of thermal behavior on stress fields are discussed and analyzed. The phenomenon of sliding friction exists widely in the field of machinery, such as gears, bearings and so on. The work done in this paper can provide some theoretical support for the research of thermomechanical coupled problem in the mechanical field.

## 2 Hypotheses and model description

The following necessary hypotheses are made in the contact model development for obtaining the semi-analytical transient solutions of temperature rise and stress:

• the material of contact bodies is same and they are perfectly elastic solids, the deformation is small;

• the contact is dry friction or boundary lubrication and the sliding friction coefficient is constant. Without considering the effects of the surface curvature micro-change and roughness;

• the area outside of the contact patches is assumed to be adiabatic, all of the dissipative mechanical energy due to the frictional work is converted into heat, and it is completely absorbed by two contact bodies with no loss of heat to the surroundings. The bodies have an initial temperature of uniform distribution.

In order to illustrate the universality of the problem, the sliding contact of non-conforming surfaces is the object of the fully coupled thermomechanical investigation, as shown in Figure 1a. The dimensions of the interface Ω between contact bodies, which are subject to the normal load Fn, are generally small when compared to the characteristic size of the bodies. Hence, each contact body can be treated as a half-space elastic solid, as shown in Figure 1b. In the figure, the symbol V is the sliding velocity of body 2 relative to body 1. The symbol Vp is the moving velocity of contact patch, Vp = V for the body 1, and Vp = 0 for the body 2. The coordinate system o-xyz is fixed with the contact body, the initial contact point is the coordinate origin and the contact surface is the xy-plane; the coordinate system o'-x'y'z' is connected with the contact patch Ω and the center of patch is the coordinate origin. The positive of x-axis (or x'-axis) is regarded as the moving direction of contact patch. The two coordinate systems coincide under the initial state.

The interface Ω has to be subjected to the distributed normal pressure p(x',y'), tangential traction s(x',y') and heat flux q(x',y') due to the friction of the contact bodies. If the two bodies are same, the tangential traction does not disturb the distribution of normal pressure, and contact characteristics are determined by the profiles of the two surfaces and the normal force [1], so the law between them can be expressed as (1a) (1b) (1c) where p0 is the centeral pressure of the Hertzian contact ellipse; a and b are the semi-axes of the ellipse in the x' and y' direction respectively; f is the sliding friction coefficient on the interface.

 Fig. 1The sliding contact between non-conforming surfaces.

## 3 Decouple of governing equations

Based on the above assumptions and simplified model, ignoring the influence of the volume force and the inertia term, the thermomechanical coupling governing equations can be denoted as (2) where λ and G are Lame's constants. β is thermal stress coefficient, β = /(1−2µ). E, µ and α are the Young's modulus, Poisson ratio and thermal expansion coefficient, respectively. ρ is the density of material. kd is the thermal diffusivity, kd = k/ρc. k and c are the thermal conductivity and specific heat capacity respectively. uζ (ζ = x, y, z) is the displacement components. ε is the volumetric strain. θ is the temperature rise, which equals to the difference of the temperature at time t and the initial temperature, θ = TtT0.

The temperature field and the stress field (or strain field) are fully coupled, they interact with each other. Therefore, the governing equations should be decoupled firstly for obtaining the transient solutions of the coupled thermomechanical contact if we want to use the analytical method. The approach that introduces the thermoelastic displacement potential Φ is used to do it, and then the decoupled governing equations can be obtained (3) (4)

## 4 Semi-analytical transient solutions for contact bodies

### 4.1 Temperature rise

The basis for the heat source method is the solution for the instantaneous point heat source [3], i.e. (5) where Q is the amount of heat liberated by the instantaneous point heat source; R2 = (x−x0)2 + (y−y0)2 + (z−z0)2, (x0, y0, z0) is the coordinate of the point heat source. The following classical Fourier heat conduction equation can be satisfied with the equation (5) (6)

It is obvious that the equations (4) and (6) are identical in form, and the only difference is the coefficients in each equation. Therefore, the temperature rise considering the thermomechanical coupling due to the plane heat source of half-space elastic solid is deduced based on the solution of the instantaneous point heat source.

For the body 1, it can be equivalent to a half-space elastic solid which is subjected to an elliptical moving plane heat source on the surface, as shown in Figure 2. In order to analyse the temperature rise of the arbitratry point in contact body at the observation time t, a new time parameter τ (known as the historical time) is required. The elliptic heat source is moved from the initial position to the dotted-line position in the graph at the τ moment. The infinitesimal A-area on the moving elliptic heat source is taken, and the heat acting on the A-area in the dτ time is expressed as (7) where γ is the heating partition coefficient (the amount of heat absorbed by the body 1 accounts for the proportion of the total frictional heating). It is inconvenient to determine the complete heating partition coefficient by doing a specific numerical analysis for each sliding cantact case, most investigators have used an approximation, such as the one first made by Blok [17]. Blok's criterion determined the heating partition coefficient can be expressed as (8)

The temperature rise of the arbitrary point (x, y, z) in contact body 1 at the observation time t, which is induced by the heat dQA, can be written as (9)

On the level of theory, the integration of the equation (9) over the entire contact ellipse and the entire time obtains the temperature rise at any point as well as at any observation time in the body 1. The integral operation may be easier if the shape of the heat source is rectangular, but for the general shape of the heat source, even with the use of numerical integration, this integral operation may also be very difficult. The method of direct integration is not adopted in this paper. The semi-analytical transient solution of the temperature rise for the contact body is deduced by using an equivalent moving point heat source instead of the moving plane heat source. As shown in Figure 2, an equivalent moving point heat source dQe is assumed at the center of the contact ellipse, and the associated temperature rise is equal to the result generating by the dQA at arbitrary point and arbitrary time in the contact body. Then the problem of moving plane heat source can be equivalent to the issue of moving point heat source, and the works may become relatively simple.

Due to the hypothetical equivalent moving point heat source dQe, the temperature rise at arbitrary point and arbitrary time in the contact body 1 is expressed as (10)

As shown in Figure 3, the three-dimensional region near the interface is selected as the computational target zone for the correlation fields. The region is discretized into M × N × L elements. The pressure on the area outside the contact ellipse is zero. It is assumed that the pressure acting on each element on the surface is uniform distribution [37]. The assumed equivalent moving point heat source Qe corresponding to the plane heat source at the τ moment can be expressed as the following semi-analytical form by uniting the equations (7), (9)(10) and integrating over the interface (11)

To make according to the solution of instantaneous point heat source, the temperature rise of the arbitrary point in the contact body 1 at the t moment is expressed as (12)

The equation (12) can be regarded as the operation of discrete convolution with considering the transformation relation between coordinate systems o-xyz and o'-x'y'z', therefore the discrete convolution and fast Fourier transform technique can be applied to improve the computational efficiency. To facilitate the presentation, the pi,j is used to denote the pressure p(x'i, y'j) of the element (i, j) on the surface, and the is used to denote the influence of the pressure pi,j on the temperature rise of the specified point (xm, yn, zl), i.e. pi,j = p(xi, yj), . The discrete convolution form of the equation (12) is expressed as (13) where m = 1, 2,…, M; n = 1, 2,…, N; l = 1, 2,…, L; The sign “**” denotes the two-dimensional discrete convolution operation; “≈” denotes the two-dimensional discrete fast Fourier transform; “IFFT” represents the inverse fast Fourier transform. For Matlab, the two-dimensional discrete convolution operation is only a directive, which is very convenient to use. Of course, we can also calculate the temperature rise using the fast Fourier transform based on the discrete convolution theorem.

For the body 2, it can be treated as a half-space elastic solid which is subjected to the elliptic stationary plane heat source on the surface due to Vp = 0. Then the temperature rise of the arbitrary point and arbitrary time in the body 2 is expressed as (14)

And the discrete convolution form (15)

Although the above semi-analytical transient solutions are developed based on the elliptical heat source and the Hertzian pressure distribution, they are equally effective for other shaped heat sources, (e.g. rectangular, circular or even irregular shapes) and pressure distribution patterns (parabolic, homogeneous, etc.). For various heat sources, the corresponding temperature rise can be obtained by altering the pressure distribution in the above-mentioned formulas.

 Fig. 2The elliptic moving heat source on the surface of the body 1.
 Fig. 3The computational target region of the temperature and stress.

### 4.2 Stress

Total stress of contact body can be decomposed into two parts: elastic stress produced by the surface distributed tractions p(x', y') and s(x', y'), and the thermal stress caused by temperature variation θ, (16)

The elastic stress near the contact interface is given by [3739] (17) where and are the influence coefficients for the normal pressure and the tangential traction, respectively, the detailed expressions are listed in Appendix.

Because of the similar expressions of the temperature rise for contact bodies, the identical approach and principle is employed to analyze the thermal stress field of the contact body 1 or body 2 due to the temperature variation. For the body 1, the alternative expression of the temperature rise generated by the moving plane heat source is given as (18)

Substituting the equation (18) into the (3), then (19) where , ▽2 is Laplace operator.

The following transformations are conducted for the integral variables

And then the particular solution of the equation (19) is expressed as (20) here, or, (21)

The thermal stress associated with the special solution of thermoelastic displacement potential can be calculated using the constitutive equation of elasticity, (22) the corresponding second-order numerical partial derivatives are calculated by using the central differences, for example (23a) (23b)

The rest of second-order partial derivatives can be calculated in a similar fashion.

The area inside the contact interface can not expand freely under the action of frictional heat and the pressure on the contact patch will increase due to the mutual constraints; while the area outside the interface is a free surface where the distribution pressure is zero. Therefore, the boundary conditions on the surface of the contact body can be expressed as (24)

The stress formed by the superposition of and can only satisfy the first boundary condition of equation (24), the following assumed normal distribution pressure should be applied to the surface of the contact body in order to satisfy the second boundary condition (25)

The stress produced by the assumed pressure in the contact body can also be calculated by the influence coefficient method. The superposition of and is the required thermal stress, i.e. (26)

The total stress considering the thermomechanical coupling for the contact body 1 is obtained by substituting equations (26) and (17) into (16).

For the contact body 2, the relevant mathematical formulas are as follows (27) (28)

## 5 Results and discussion

Numerical simulation results of transient temperature rise and stress of the contact interface which is simulated by the sliding contact of an ellipsoid and a half-space solid is presented in this section. The materials of two contact bodies are assumed to be bearing steel (GCr15), and the mechanical and thermal properties of GCr15 [40] and the simulation parameters are summarized in Table 1. The Cartesian coordinate equation of the ellipsoid is , and the coordinates of the contact point with the half-space body on the surface of the elliposid is (0,0,4), then the corresponding main curvature radiuses at the contact point are R1 = 9 mm, R2 = 4 mm. The Hertzian solutions [1] of contact elliptical semi-axes and peak pressure are a = 0.1918 mm, b = 0.1279 mm and p0 = 1946.8 MPa when the load Fn = 100 N.

Table 1

Parameters and material (GCr15) properties in simulations.

### 5.1 Verification

In order to verify the current method, the steady-state maximum temperature rise (Fig. 4) and the pressure distribution (Fig. 5) on the contact interface are calculated based on different approaches [4144], and the results obtained are compared. The comparison results indicate that the solutions obtained by present method have an approximate agreement with those calculated by other investigators' methods.

 Fig. 4Comparison of the steady-state maximum temperature rise.
 Fig. 5Comparison of surface pressure along the x-axis for body 1.

### 5.2 Temperature rise

Owing to the heat allocation scheme from Blok's work is adopted in this paper, the maximum surface temperature rises of the two contacting bodies can only be guaranteed to be equal, as shown in Figure 6, and it is impossible to ensure that the temperatures at all points within the contact interface are equal. It can be seen from the data in the figure that although the maximum surface temperature rises of two contacting bodies are equal, the gradient of surface temperature rise for the body 1 is greater than the temperature rise gradient of the body 2, this indicats that the movement of heat source will aggravate the uneven degree of temperature distribution on the contact surface.

For the sliding contact problems with elliptical contact patch, the dimensionless Péclet number Pe is usually defined as Pe = Va/2kd, it is an extremely useful parameter for the relative velocity of motion of the heat source considering the thermal properties of the conduction medium which determines the speed of dissipation of heat in the medium. According to Blok's heat allocation scheme, Figure 7a presents that the variation of heating partition coefficient with sliding time under various Péclet numbers, and the variation of corresponding maximum contact temperature rise is shown in Figure 7b. As shown in Figure 7a, it is known that the heating partition coefficient is approximately equal to 0.5 when the Péclet number Pe ≤ 0.1; for Pe ≥ 0.1, the greater the relative sliding velocity and the longer the sliding time, the higher the proportion of heat distributed to the contact body 1. Figure 7b manifests that the contact interface temperature rise can reach the quasi-steady state in a very short time during the sliding process under the current constant velocity and constant force, the current quasi-steady state time is about 0.1 s. However, there are few cases of constant force and constant velocity in practical engineering. In most cases, the force and speed are time-varying, such as alternating load in bearing and frequent variations of velocity in drive system, and it will take a long time to reach quasi-steady state in these cases. The transient models developed here can solve the thermal problems of these cases and contribute to better understand the changing process of interface temperature. After the quasi-steady state time, the change of the contact temperature rise is very small and the increase of corresponding heating partition coefficient has become smaller.

The difference of surface temperature rises between the coupled and uncoupled models for body 1 is depicted in Figure 8, which can reflect the effect of thermomechanical coupled term on the temperature rise. When Pe ≤ 1, it means that the relative sliding velocity between contact bodies is small, the difference between coupled and uncoupled models on the temperature rise is small, hence the coupled term can be neglected for this situation, but there are few cases of small sliding velocity in actual sliding contact process, the sliding velocity is larger in most cases of sliding contact; for the case with the Pe ≥ 5, the effect of coupled term must be considered, it is obvious that the effect of the thermomechanical coupled term is more significant with the increase of the sliding velocity. According to the trend shown in the figure, for the cases with large sliding velocity (e.g. automotive braking, grinding process), the coupled model can more accurately calculate the temperature rises. Besides, the figure also shows that the higher the temperature rise on the contact patch, the greater the effect of thermomechanical coupled term.

The effects of sliding velocity, friction coefficient, initial temperature on the temperature rise are described in Figure 9. The curves in the figure demonstrate that the variations of sliding velocity, frictional coefficient can significantly affect the temperature rise, and the temperature rise is not sensitive to the variation of the initial temperature. Also, the surface temperature rise of the contact body is rapidly reduced from the peak position to the surrounding area, this phenomenon indicates that the heat affected zone of the contact bodies is concentrated only in a small region near the contact patches.

For the case with a sliding velocity which belongs to the Péclet number Pe ≤ 1, as shown in Figure 9a1, the curves reveal that the characteristics of surface temperature rise distributions of the body 1 for the low sliding velocity are still similar to the symmetrical distribution of the contact center, and the peak temperature rises are located on the contact center. However, the maximum temperature rises markedly shift from the center to the rear edge of the contact region when Pe ≥ 5, and the trend of temperature with the sliding velocity is enhanced.

 Fig. 6Surface temperature rise distribution. (a): body 1; (b): body 2.
 Fig. 7Variation of the heating partition coefficient and the maximum temperature rise.
 Fig. 8Difference of surface temperature rise along the x-axis between coupled and uncoupled models for body 1.
 Fig. 9The effects of sliding velocity, friction coefficient, initial temperature on temperature rise, (a1)∼(a3): body 1; (b1)∼(b3): body 2.

### 5.3 Stress

The surface pressure distribution of the body 1 calculated by the different models is shown in Figure 10. Figure 10a is obtained by the Hertzian model and Figure 10b is obtained by the model discribed herein. Figure 10 manifests that the boundary conditions on the surface of the contact body can be satisfied by applying the assumed normal distribution pressure listed in equation (24), i.e. the surface outside the contact zone is free, the pressure distribution on it is zero, and the pressure on the contact patch is increased due to the surface inside the contact patch can not expand freely. Figure 11 shows that the surface pressure considering the thermal effect of the body 1 under various sliding velocities. For the case of Pe ≤ 1, due to the contact temperature rise caused by friction heat is small, the effect of thermal behavior on the surface pressure is weak; but for Pe ≥ 5, the friction heat has a significant effect on the surface contact pressure.

On the surface of the contact body, the effect of the thermal stress produced by the uneven temperature distribution on the each stress component along the x-axis is illustrated in Figure 12. Negative values in the figures mean compressive stress and positive values represent the tensile stress. For the three normal stress components, as shown in Figure 12ac, the thermal stress increases the compressive stress on the contact patch, and the effect of thermal stress on the normal stresses σ1x and σ1y is greater than that on the normal stress σ1z; as the free expansion of the surface outside the contact patch, the thermal stress only increases the value of the normal stresses σ1x and σ1y, and the normal stress σ1z remains at zero.

Figure 12df show the effect of thermal stress on the three shear stress components on the surface of contact body, the shear stress component τ1xy has both tensile shear stress and compressive shear stress on the contact patch, and the thermal stress reduces the tensile shear stress and increases the compressive shear stress, the value of τ1xy is small by comparing with the normal stress; the surface pressure without considering the thermal effect is symmetric about the x-axis, so the elastic shear stress along the x-axis is zero, but the slight compressive shear stress along the x-axis will be generated when the thermal effect is considered; for the shear stress component τ1zx, the thermal stress increases the compressive shear stress on most contact areas and slightly decreases the value of τ1zx on the small part of the contact area, the thermal stress causes the compressive shear stress on the area ahead of the patch and causes the tensile shear stress on the rear area of the patch.

The above-mentioned analyses indicate that the thermal stress has different effects on the various stress components, and it can not be ascertained whether the stress state of a specified point in the elastic body is enhanced or weakened by the thermal stress, also the stress state of a specified point can not be determined only by a certain stress component. The Mises equivalent stress is defined as it takes into account the influence of each stress component, and its value is independent of the selection of the coordinate system. Therefore, the best choice to measure the stress state of the elastic body is Mises equivalent stress. The effect of thermal stress on the Mises equivalent stress on the surface of contact body is shown in Figure 13. By comparing three curves in the figure, it can be seen that the thermal stress will decrease the Mises equivalent stress on the contact patch and increase the value of equivalent stress on the surface outside the patch. The reason why Mises equivalent stress shows this feature can be explained by the definition of Mises equivalent stress and Figure 14. The total Mises equivalent stress is not a simple superposition of the elastic Mises equivalent stress and the thermal Mises equivalent stress when considering the thermal effect, it should be determined by the value of the difference between the normal stress components and the value of the shear stress components. As shown in Figure 12, the shear stress is small compared to normal stress and thus the normal stress is considered as a major factor affecting Mises equivalent stress; according to the information reflected in Figure 14, it is known that the difference value between the normal stress components on the contact patch is decreased and it is increased for the surface outside the contact patch when considering the thermal effect, therefore the Mises equivalent stress on the surface of contact body shows the feature shown in Figure 13.

 Fig. 10Distribution of surface pressure for body 1.
 Fig. 11Surface pressure along the x-axis under various sliding velocities for body 1.
 Fig. 12The effect of thermal stress on each stress component along the x-axis on the surface of body 1.
 Fig. 13The effect of thermal stress on the Mises equivalent stress on the surface of body 1 along the x-axis.
 Fig. 14The difference between the normal stress components on the surface of body 1 along the x-axis.

## 6 Conclusions

In this paper, the plane-point heat source equivalent method, the thermoelastic displacement potential and the central difference were used to derive the semi-analytical solutions of the transient temperature rise and stress of the fully thermomechanical coupling for the sliding contact. The presented method had be verified by using comparative analysis approach, the effects of thermomechanical coupled term, sliding velocity, frictional coefficient and initial temperature on temperature rise and the effect of thermal stress on stress fields were discussed, the results indicate that

during the sliding contact process, for the case of the Péclet number Pe ≤ 1, the effect of thermomechanical coupled term can not be considered when calculating the contact temperature rise, but it will have a significant effect on the temperature rise for the Pe ≥ 5;

• the contact temperature rise of bodies is greatly affected by the sliding velocity and the friction coefficient, the effect of the initial temperature on the temperature rise is small;

• when Pe ≤ 1, the effect of thermal behavior on the surface pressure is weak, and increasing sliding velocity leads to a noticeable increase in surface contact pressure for the case of Pe ≥ 5;

• in general, the thermal stress increases the various stress components on the surface of contacting bodies. However, it has different effects on the Mises equivalent stress on the surface inside and outside the contact area: the Mises equivalent stress outside the contact area is increased and the value on the contact area is decreased.

## Nomenclature

a, b: Semi-axes of the contact ellipse (mm)

c: Specific heat capacity (J kg−1 K−1)

E: Young's moduls (MPa)

f: Sliding friction coefficient

G: Lame's constant

(i, j): Element number on the surface, i = 1, 2,…, M; j = 1, 2,…, N

k: Thermal conductivity (W mm−1 K−1)

kd: Thermal diffusivity, kd = k/ρc (mm2 K−1)

m, n, l: m = 1, 2,…, M; n = 1, 2,…, N; l = 1, 2,…, L

M, N, L: Amount of elements in the , and directions

p: Normal pressure (MPa)

p0: The centeral pressure of contact ellipse (MPa)

p*t: Assumed pressure to satisfy the boundary condition (MPa)

q: Heat flux (W mm−2)

Q: Amount of heat liberated by the instantaneous point heat source (J)

Qe: Hypothetical equivalent point heat source (J)

R1, R2: Main curvature radiuses (mm)

s: Tangential traction (MPa)

t: Time (s)

Tt, T0: Temperature at time t and the initial temperature (K)

uζ: Displacement component

V: Sliding velocity of body 2 relative to body 1 (mm s−1)

Vp: Moving velocity of contact patch (mm s−1)

x, y, z, , , : Space coordinates

Δxʹ, Δyʹ, Δzʹ: Intervals in the , and directions

α: Thermal expansion coefficient (K−1)

β: Thermal stress coefficient, β = /(1-2µ)

λ: Lame's constant

µ: Poisson ratio

ρ: Density of material (kg mm−3)

τ: Time parameter (s)

γ: Heating partition coefficient

ε: Volumetric strain

θ: Temperature rise (K)

η, ζ: x, y or z

: Total stress, elastic stress, thermal stress (MPa)

σx, σy, σz: Normal stress components (MPa)

τxy, τyz, τzx: Shear stress components (MPa)

σvm: Mises equivalent stress (MPa)

Ω: Contact patch

Φ, Φ*: Thermoelastic displacement potential and its a particular solution

2: Laplace operator

Pe: Peclet numbers, Pe = Va/2kd

IFFT: Inverse fast Fourier transform

: Influence coefficients for normal pressure and tangential traction

**: Two-dimensional discrete convolution operation

≈: Two-dimensional discrete fast Fourier transform

## Acknowledgements

The authors would like to gratefully acknowledge the financial support from the National Natural Science Foundation of China (No.51275440).

## Appendix

The influence coefficients and can be expressed as [39] (A.1) (A.2) where x+ = xmxi+Δx/2, x = xmxiΔx/2, y+ = ynyj + Δy/2, y = ynyjΔy/2. (A.3)and, (A.4)where R = [x2 + y2 + z2]1/2; .

## References

1. K.L. Johnson, Contact mechanics, Cambridge University Press, London, 1985 [Google Scholar]
2. Y. Zongda, W. Hongli, Thermal stress, Higher Education Press, Beijing, 1993 [Google Scholar]
3. J.C. Jaeger, Moving sources of heat and the temperature at sliding contacts, J. Proc. R. Soc. New South Wales, 76 (1942) 203–224 [Google Scholar]
4. H. Chen, Y. Hu, H. Wang, et al., Calculation of temperature fields of bodies in sliding contact without lubrication, J. Tsinghua Univ. Sci. Technol. 47 (2007) 1962–1964 [Google Scholar]
5. Z.B. Hou, R. Komanduri, General solutions for stationary/moving plane heat source problems in manufacturing and tribology, Int. J. Heat Mass Transfer 43 (2000) 1679–1698 [CrossRef] [Google Scholar]
6. P. Levin, A general solution of 3-D quasi-steady-state problem of a moving heat source on a semi-infinite solid, Mech. Res. Commun. 35 (2008) 151–157 [CrossRef] [Google Scholar]
7. R. Bosman, M.B. de Rooij, Transient thermal effects and heat partition in sliding contacts, J. Tribol. 132 (2010) 021401–021409 [CrossRef] [Google Scholar]
8. W.L. Kuo, J.F. Lin, General temperature rise solution for a moving plane heat source problem in surface grinding, Int. J. Adv. Manuf. Technol. 31 (2006) 268–277 [CrossRef] [Google Scholar]
9. G. Araya, G. Gutierrez, Analytical solution for a transient, three-dimensional temperature distribution due to a moving laser beam, Int. J. Heat Mass Transfer 49 (2006) 4124–4131 [CrossRef] [Google Scholar]
10. H. Belghazi, M. El Ganaoui, J.C. Labbe, Analytical solution of unsteady heat conduction in a two-layered material in imperfect contact subjected to a moving heat source, Int. J. Therm. Sci. 49 (2010) 311–318 [CrossRef] [Google Scholar]
11. M. Hamraoui, T. Osman, A. Boucheffa, et al., Analytical modelling of the three-dimensional steady-state temperature in a bearing ring, Mech. Ind. 12 (2011) 1–4 [Google Scholar]
12. F.E. Kennedy, X. Tian, Modeling sliding contact temperatures, including effects of surface roughness and convection, J. Tribol. 138 (2016) 042101 [CrossRef] [Google Scholar]
13. F.E. Kennedy, Y. Lu, I. Baker, Contact temperatures and their influence on wear during pin-on disk tribotesting, Tribol. Int. 82 (2015) 534–542 [CrossRef] [Google Scholar]
14. N. Alilat, A. Bairi, N. Laraqi, Three-dimensional calculation of temperature in a rotating disk subjected to an eccentric circular heat source and surface cooling, Numer. Heat Transfer, 46 (2004) 167–180 [CrossRef] [Google Scholar]
15. N. Laraqi, N. Alilat, J.M. Garcia de Maria, A. Bairi, Temperature and division of heat in a pin-on-disc frictional derive−exact analytical solution, Wear, 266 (2009) 765–770 [CrossRef] [Google Scholar]
16. D.G. Bansal, J.L. Streator, On estimation of maximum and average interfacial temperature rise in sliding elliptical contacts, Wear, 278–279 (2012) 18–27 [CrossRef] [Google Scholar]
17. H. Blok, Theoretical study of temperature rise at surfaces of actual contact under oiliness lubricating condition, in: Proceeding of the Institute of Mechanical Engineers General Discussion of Lubrication, Institute of Mechanical Engineers, London, 1937, pp. 222–235 [Google Scholar]
18. X. Tian, F.E. Kennedy, Maximum and average flash temperatures in sliding contacts, J. Tribol. 116 (1994) 167–174 [CrossRef] [Google Scholar]
19. H.W. Wu, Y.Y. Chen, J.H. Horng, Contact temperature under three-body dry friction conditions, Wear, 330–331 (2015) 85–92 [CrossRef] [Google Scholar]
20. A. Belhocine, M. Bouchetara, Transient thermal Ansys analysis of dry contacts-Application to automotive braking, Mech. Ind. 13 (2012) 45–57 [CrossRef] [Google Scholar]
21. Y. Wang, W. Tang, Y. Chen, et al., Investigation into the meshing friction heat generation and transient thermal characteristics of spiral beral gears, Appl. Therm. Eng. 119 (2017) 245–253 [CrossRef] [Google Scholar]
22. H. Xiao, D. Tang, Z. Deng, et al., Thermal analysis and experimental verification of the transmission in a lunar drilling system, Appl. Therm. Eng. 113 (2017) 765–773 [CrossRef] [Google Scholar]
23. W. Li, J. Tian, Unsteady-state temperature field and sensitivity analysis of gear transmission, Tribol. Int. 116 (2017) 229–243 [CrossRef] [Google Scholar]
24. J. Wen, M.M. Khonsari, Thermomechanical effects on transient temperature in non-conformal contacts experiencing reciprocating sliding motion, Int. J. Heat Mass Transfer, 52 (2009) 4390–4399 [CrossRef] [Google Scholar]
25. P. Hwang, X. Wu, Investigation of temperature and thermal stress in ventilated disc brake based on 3D thermo-mechanical coupling model, J. Mech. Sci. Technol. 24 (2010) 81–84 [Google Scholar]
26. M.N. Balci, S. Dag, B. Yildirim, Subsurface stresses in graded coatings subjected to frictional contact with heat generation, J. Therm. Stresses, 40 (2017) 517–534 [CrossRef] [Google Scholar]
27. G. Zhou, L. Hua, D. Qian, et al., Effects of axial rolls motions on radial-axial rolling process for large-scale alloy steel ring with 3D coupled thermo-mechanical FEA, Int. J. Mech. Sci. 59 (2012) 1–7 [CrossRef] [Google Scholar]
28. A. Zahedi, M.R. Movahhedy, Thermo-mechanical modeling of high speed spindles, Sci. Iran. B, 19 (2012) 282–293 [CrossRef] [Google Scholar]
29. F.E. Kennedy Jr, Thermal and thermomechanical effects in dry sliding, Wear, 100 (1984) 453–476 [CrossRef] [Google Scholar]
30. S. Liu, Q. Wang, A three-dimensional thermomechanical model of contact between non-conforming rough surfaces, J. Tribol. 123 (2001) 17–26 [CrossRef] [Google Scholar]
31. S. Liu, Q. Wang, Transient thermoelastic stress field in a half-space, J. Tribol. 125 (2003) 33–43 [CrossRef] [Google Scholar]
32. W.W. Chen, Q.J. Wang, Thermomechanical analysis of elastoplastic bodies in a sliding spherical contact and the effects of sliding speed, heat partition, and thermal softening, J. Tribol. 130 (2008) 041402–041410 [CrossRef] [Google Scholar]
33. L.J. Wang, Z.Y. Ai, Plane strain and three-dimensional analysis for thermo-mechanical behavior of multilayered transversely isotropic materials, Int. J. Mech. Sci. 103 (2015) 199–221 [CrossRef] [Google Scholar]
34. G. Altan, M. Topcu, N.B. Bektas, et al., Elastic-plastic thermal stress analysis of an aluminum composite disc under parabolic thermal load distribution, J. Mech. Sci. Technol. 22 (2008) 2318–2327 [CrossRef] [Google Scholar]
35. Y. Zhang, Y.Y. Xiang, A semi-analytical method and its application for calculating the thermal stress and displacement of sparsely fractured rocks with water flow and heat transfer, J. Zhejiang Univ. Sci. A Appl. Phys. Eng. 16 (2015) 922–934 [CrossRef] [Google Scholar]
36. R.D. Kulchytsky-Zhyhailo, Z.S. Olesiak, Stress distribution in rotating solids with frictional heat excited over contact region, J. Therm. Stresses, 29 (2006) 957–972 [CrossRef] [Google Scholar]
37. J.J. Kalker, Numerical calculation of the elastic field in a half-space, Commun. Appl. Numer. Methods, 2 (1986) 401–410 [CrossRef] [Google Scholar]
38. S. Liu, Q. Wang, G. Liu, A versatile method of discrete convolution and FFT (DC-FFT) for contact analyses, Wear, 243 (2000) 101–111 [CrossRef] [Google Scholar]
39. S. Liu, Q. Wang, Studying contact stress fields caused by surface tractions with a discrete convolution and fast Fourier transform algorithm, J. Tribol. 124 (2002) 36–45 [CrossRef] [Google Scholar]
40. X. Li, H. Ding, Z. Tang, Study of themo-physical properties for GCr15 bearing steel in continuous casting, J. Mater. Metall. 9 (2010) 241–244 [Google Scholar]
41. K. Knothe, S. Liebelt, Determination of temperatures for sliding contact with application for wheel-rail systems, Wear, 189 (1995) 91–99 [CrossRef] [Google Scholar]
42. F.E. Kennedy, Contact temperature of a moving solid surface, Springer US, 2013, pp. 544–548 [Google Scholar]
43. M.A. Tanvir, Temperature rise due to slip between wheel and rail-an analytical solution for Hertzian contact, Wear, 61 (1980) 295–308 [CrossRef] [Google Scholar]
44. J.R. Barber, Distortion of the semi-infinite solid due to transient surface heating, Int. J. Mech. Sci. 14 (1972) 377–393 [CrossRef] [Google Scholar]

Cite this article as: Z. An, P. Sun, Semi-analytical transient solutions of fully thermomechanical coupled problem for sliding contact, Mechanics & Industry 19, 405 (2018)

## All Tables

Table 1

Parameters and material (GCr15) properties in simulations.

## All Figures

 Fig. 1The sliding contact between non-conforming surfaces. In the text
 Fig. 2The elliptic moving heat source on the surface of the body 1. In the text
 Fig. 3The computational target region of the temperature and stress. In the text
 Fig. 4Comparison of the steady-state maximum temperature rise. In the text
 Fig. 5Comparison of surface pressure along the x-axis for body 1. In the text
 Fig. 6Surface temperature rise distribution. (a): body 1; (b): body 2. In the text
 Fig. 7Variation of the heating partition coefficient and the maximum temperature rise. In the text
 Fig. 8Difference of surface temperature rise along the x-axis between coupled and uncoupled models for body 1. In the text
 Fig. 9The effects of sliding velocity, friction coefficient, initial temperature on temperature rise, (a1)∼(a3): body 1; (b1)∼(b3): body 2. In the text
 Fig. 10Distribution of surface pressure for body 1. In the text
 Fig. 11Surface pressure along the x-axis under various sliding velocities for body 1. In the text
 Fig. 12The effect of thermal stress on each stress component along the x-axis on the surface of body 1. In the text
 Fig. 13The effect of thermal stress on the Mises equivalent stress on the surface of body 1 along the x-axis. In the text
 Fig. 14The difference between the normal stress components on the surface of body 1 along the x-axis. In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.