Sage Journals: Discover world-class research

Abstract

The contact stiffness and damping play an important role in determining the dynamic characteristics at various mechanical joint surfaces in an assembly body. In this article, a bi-fractal model is proposed to overcome the drawbacks of the G-W and M-B models, which only considered the contacts that have the same fractal parameters for the joint surfaces. A surface impact coefficient is also introduced to calculate the real contact area from the microscopic and macroscopic perspectives, which can be used to obtain the contact stiffness and damping of rough surfaces by applying the contact formula between rough planes. In addition, the contact body is approximated to a number of asperities stacked on a base, similar to a series structure of two springs. The series structure can reasonably account for the limitations of the contact stiffness and damping at joint surfaces. With these improvements, this study establishes a more precise joint surface model determining the stiffness–load and damping–load relationships. An experiment is performed to compare the contact stiffness and damping determined from this theory with the experiment results. The theoretical values correspond well with the experimental values. The contact stiffness increases non-linearly with the load; the increase is initially quick and then slows down and stabilises. A value of D = 1.42 is obtained as an estimated critical value which determined the contact damping increase or decrease with load. The model can predict the trend in the contact stiffness and damping based on the loading to better prepare for a dynamic assembly analysis.

Keywords

Joint surface Hertz contact G-W fractal contact M-B fractal contact

Introduction

The contact interfaces between parts, components, and assemblies are called mechanical joint surfaces (joint surfaces hereafter). These interfaces belong to the flexible bodies and withstand contact pressures. Furthermore, geometry errors, microscopic roughness, and a variety of oil mediums may exist between these interfaces. Any dynamic load will induce a slight relative linear displacement between the joint surfaces, storing and consuming energy both elastically and plastically via the contact stiffness and damping.¹ The contact stiffness and damping of a joint surface is often a significant contributor to the overall stiffness and damping of an assembly, and dynamic studies have demonstrated that approximately 60%–80% of the flexibility and more than 90% of the damping in an assembled body were caused by the joint surfaces.^2,3 Researches on the contact stiffness and damping at a joint surface are crucial for prediction dynamics, guidance alignment, and optimisation decisions for the assembly. The parts are assembled according to the process sequence. The different contact states of each joint surface are determined by different assembly processes, which form the isomorphic and heterogeneous contact stiffness and damping matrices (K, C). These matrices were used to establish the kinetic equation ( $M \overset{\cdot\cdot}{X} + C \overset{\cdot}{X} + K X = F$ , where M, C, K, and F are mass matrix, damping matrix, stiffness matrix, and applied force vector, respectively), which would affect the dynamics characteristics of the assembly bodies.⁴

Most of the existing researches on the contact characteristics of joint surfaces includes the Hertz contact model,^5–7 G-W contact model,⁸ and M-B contact model.⁹

Hertz theory is the first to solve the contact problem of frictionless elastic homogeneous with small deformation due to normal loading, which was mainly used to deal with the contact surface shape, contact area, and contact stress between smooth surfaces in the macroscopic domain without the consideration of microscopic characteristics, such as roughness, processing methods, and surface morphology. It presented a index relationship of load–area $P ~ A_{r}^{3 / 2}$ ( $P$ is the load between joint surfaces and $A_{r}$ is the real contact area) which reduced another index relationship of contact stiffness–area $K ~ A_{r}^{1 / 2}$ ( $K$ is the contact stiffness between joint surfaces). Hertz contact theory did not consider the transition process from elastic to plastic deformation. No energy was lost by plastic deformation during the contact process. Therefore, the contact damping was unknown. Furthermore, in reality, there are no absolutely smooth surfaces.

In the following, Hertz theory was extended by Greenwood and Williamson (G-W model), who approximated the surface topography by asperities possessing spherical summits with uniform radius of curvature and Gaussian height distribution, and assumed purely elastic deformation of the interacting asperities which are sufficiently apart from each other in order for mechanical interaction to be secondary. The G-W model examined the contact attributes of the joint surfaces statistically by relating two material properties (the elastic modulus $E$ and hardness $H$ ) and three microscopic surface properties (the surface density of asperity $η$ , asperity height standard deviation $σ$ , and asperity average curvature radius $β$ ). The plasticity index ( $ψ = (E / H) \sqrt{σ / β}$ ) was introduced to determine the deformation nature of the joint surface qualitatively. However, it cannot be utilised to calculate the energy stored via elastic deformation and the energy consumed via plastic deformation quantificationally. Therefore, the contact damping remained unknown. Many studies on the extension and variation of G-W model have also been proposed by Pullen and Williamson,¹⁰ Powierza et al.,¹¹ Greenwood and Tripp,¹² Chang et al.,¹³ McCool and Gassel,¹⁴ Halling and Nuri,¹⁵ and Stanley et al.¹⁶ whereas they still maintained the relative simplicity of the G-W statistical representation of a surface. Since the above statistical approaches based on the G-W model were limited by the dependence of statistical roughness parameters on the sampling length and resolution of the measuring instrument, they could not provide objective in its uniqueness and steady information for the surface topography.

A more sophisticated approach to handling the contact between rough surfaces was offered by Majumdar and Bhushan (M-B model), who developed a fractal geometry (the fractal dimension D and fractal roughness G was introduced) model with intrinsic advantages of scale independent and self-affinity. One criticism of the M-B fractal model is that it predicted that the lighter the load is, the greater will be the percentage of contacts that are plastically deformed. This counterintuitive result arises from the fact that, for a given surface, the radius of curvature of an asperity is determined solely by the contact diameter. The M-B fractal model viewed the joint surfaces as a number of spherical contact asperities and defined distribution function of asperities^17,18 as $n (a) = (D / 2) \cdot a^{D / 2} / a^{(D / 2 + 1)}$ (where $a$ and $a_{l}$ are the actual and maximum contact areas of an asperity, respectively.). In addition, the critical contact area $a_{c} = G^{2} / {(H / (2 \cdot E))}^{2 / (D - 1)}$ was introduced to determine the deformation nature of asperities. When $a$ is larger than $a_{c}$ , the asperities are elastically deformed, while $a$ is smaller than $a_{c}$ , the asperities are plastically deformed. The fact that the relationship of stiffness–load of single spherical asperity (which was called smooth surface in Hertz model) combined with distribution function of asperities and the critical contact area could deduce the objective and stable contact stiffness $K_{n_M - B}$ of the M-B model. Furthermore, after the critical contact area was introduced, the plastic strain energy $ω_{p}$ and elastic strain energy $ω_{e}$ could be calculated, which could deduce the contact damping of joint surfaces¹⁹ $C_{n} = (ω_{p} / ω_{e}) {(m K_{n_M_B})}^{1 / 2}$ ( $m$ is the mass of joint surface). Other fractal-based models have also been developed by a number of researchers including Majumdar and Tien,²⁰ Yan and Komvopoulos,²¹ and Ciavarella et al.²² A recent discussion of Archard’s results in view of modern fractal theory was given by Ciavarella and Demelio.²³ Fractal contact models taking elastoplastic deformation into account were presented by Majumdar and Bhushan as an extension to their approach based on the Weierstrass function and by Warren and Krajcinovic,²⁴ who used a Cantor-set representation of the rough surface in conjunction with a Winkler-type foundation model. The pioneering work²⁵ has been done to estimate load–compliance relation for the contact of fractal surface with a rigid plane by using Westergaard’s solution. Komvopoulos and Ye²⁶ developed a comprehensive contact analysis that incorporates the successive occurrence of elastic, elastic–plastic, and fully plastic deformation of contacting asperities on rough surfaces with a fractal topography description. It has been shown that surface profiles formed by the electric discharge machining,²⁷ cutting, milling or grinding,^28,29 and worn surfaces³⁰ have fractal structures, and fractal parameters can reflect the intrinsic properties of surface. The M-B fractal contact model broke the restriction of Hertz contact model, namely, that the contact bodies must be smooth surfaces with different curvatures. Furthermore, the fractal dimension D and fractal roughness G of the rough surface are independent of the sampling length and instrument resolution with objective certainty and scale independence relative to the G-W model. Therefore, the fractal theory has a wider range of application than Hertz contact theory and the G-W model.

But the M-B model should be actually called the single-fractal contact model between rough plane which required contacting bodies to be with same fractal parameters and be rough planes. This article develops specific calculations and comprehensive discussions for the contact stiffness and damping of joint surfaces based on the M-B model. The bi-fractal model is introduced to overcome the limitation of requiring same fractal parameters, and the surface impact coefficient is proposed to deal with surface contact by applying the formula of planar contact. Besides, the series structure of asperities stacking on the base is put forward which can reasonably account for the finiteness of contact stiffness and damping at joint surfaces. The results obtained from the proposed model correspond well with the experimental values. The involved parameters and meaning are presented in Appendix 1.

Methods

Algorithm improvements

Bi-fractal contact model

The M-B model established a mechanical model of the joint surface via fractal geometry, which could avoid data instabilities. In the M-B model, the fractal parameters of the two rough surfaces were assumed to be equal (D = D₁ = D₂ and G = G₁ = G₂). The fractal parameters could also be equivalent to an elastic rough surface (D and G) and rigid plane in contact, so that the M-B fractal model is also called a single-fractal model. In fact, the fractal parameters of the two rough surfaces are often different (D₁≠D₂ and G₁≠G₂). In this case, the two single-fractal surfaces that are in contact are equivalent to a bi-fractal rough surface and rigid plane in contact, and the bi-fractal rough surface has an equivalent elastic modulus E and equivalent fractal parameters (D and G) (see Figure 1).

Figure 1.

Two single-fractal surfaces are simplified into a bi-fractal surface and rigid plane.

The equivalent structure function of the two single-fractal surfaces (D₁, G₁, D₂, and G₂) in contact^31–33 is

\begin{matrix} S (τ) = S_{1} (τ) + S_{2} (τ) = C_{1} G_{1}^{2 (D_{1} - 1)} τ^{(4 - 2 D_{1})} \\ + C_{2} G_{2}^{2 (D_{2} - 1)} τ^{(4 - 2 D_{2})} \end{matrix}

(1)

where $τ$ is the displacement along a certain axis, and $S_{1} (τ) and S_{2} (τ)$ are the structure function of contact interfaces 1 and 2, respectively. $C_{1} and C_{2}$ are the constants determined by the fractal parameters of contact interfaces.

When $D_{1} = D_{2}$ , the equivalent fractal dimension D and equivalent fractal roughness G are defined as

\begin{matrix} D = D_{1} = D_{2} \\ G^{2 (D - 1)} = G_{1}^{2 (D_{1} - 1)} + G_{2}^{2 (D_{2} - 1)} \end{matrix}

(2)

When $D_{1} \neq D_{2}$ , the equivalent structure function with double logarithmic is shown in Figure 2.

Figure 2.

Logarithmic curve of the structure function of a bi-fractal surface.

Where $a_{12} = τ_{12}^{2}$ can be obtained from $C_{1} G_{1}^{2 (D_{1} - 1)} τ^{(4 - 2 D_{1})} = C_{2} G_{2}^{2 (D_{2} - 1)} τ^{(4 - 2 D_{2})}$ . Based on the M-B single-fractal contact theory, the bi-fractal contact can be divided into one of the following four cases:

Case A ( $a_{c 1} < a_{12}$ ): In the case of $a < a_{c 1}$ , the asperities deform plastically, and in the case of $a_{c 1} < a < a_{12}$ , the asperities deform elastically.

Case B ( $a_{c 1} > a_{12}$ ): All asperities in the fractal region I deform plastically.

Case C ( $a_{c 2} < a_{12}$ ): All asperities in the fractal region II deform elastically.

Case D ( $a_{c 2} > a_{12}$ ): In the case of $a_{12} < a < a_{c 2}$ , asperities deform plastically, and in the case of $a < a_{c 2}$ , the asperities deform elastically.

The normal characteristics of a single asperity are as follows:

The plastic load is defined as $f_{p} = Ha$ .

The elastic load is defined as $f_{e} = (4 / 3) {E R}^{1 / 2} δ^{3 / 2} = 4 / 3 π^{1 / 2} E G^{(D - 1)} a^{(3 - D) / 2}$ .

The normal contact stiffness is defined as $k_{n} = (4 E / 3 \sqrt{π}) a^{1 / 2}$ .

The normal contact damping is defined as $c_{n} = η {(m k_{n})}^{1 / 2}$ .

The distribution function for the contact area is divided into the following two cases

Case 1 (a_{l} < a_{12}) n (a) = \frac{D_{1}}{2} (\frac{a_{l}^{\frac{D_{1}}{2}}}{a^{(\frac{D_{1}}{2} + 1)}})

(3)

Case 2 (a_{1} > a_{12}) {\begin{matrix} n_{1} (a) = \frac{D_{1}}{2} {(\frac{a_{1}}{a_{12}})}^{\frac{D_{2}}{2}} (\frac{a_{12}^{\frac{D_{1}}{2}}}{a^{(\frac{D_{1}}{2} + 1)}}) a < a_{12} \\ n_{2} (a) = \frac{D_{2}}{2} (\frac{a_{1}^{\frac{D_{2}}{2}}}{a^{(\frac{D_{2}}{2} + 1)}}) a_{12} < a < a_{1} \end{matrix}

(4)

In Case 1, the asperity deformation does not enter fractal region II, so the stiffness and damping is determined in accordance with the traditional M-B single-fractal model using the equivalent fractal dimension D and fractal roughness G.

In Case 2, the asperity deformation has entered fractal region II and can be divided into the following four situations to determine the loading, contact stiffness, and contact damping (subscript $n_M_B_D$ indicates the normal direction of bi-fractal M-B model):

Cases A and C: $a_{c 1} < a_{12} and a_{c 2} < a_{12}$

\begin{matrix} F = \int_{0}^{a_{c 1}} f_{p} (a) n_{1} (a) da + \int_{a_{c 1}}^{a_{12}} f_{e} (a) n_{1} (a) da \\ + \int_{a_{12}}^{a_{l}} f_{e} (a) n_{2} (a) da \\ K_{n_M_B_D} = \int_{a_{c 1}}^{a_{12}} k_{n} (a) n_{1} (a) da + \int_{a_{12}}^{a_{l}} k_{n} (a) n_{2} (a) da \\ C_{n_M_B_D} = \frac{\frac{1}{2} \int_{0}^{a_{c 1}} f_{p} (a) δ (a) n_{1} (a) da}{\frac{1}{2} \int_{a_{c 1}}^{a_{12}} f_{e} (a) δ (a) n_{1} (a) da + \frac{1}{2} \int_{a_{12}}^{a_{l}} f_{e} (a) δ (a) n_{2} (a) da} \\ \cdot \sqrt{m K_{n_M_B_D}} (5) \end{matrix}

(5)

Cases A and D: $a_{c 1} < a_{12} and a_{c 2} > a_{12}$

\begin{matrix} F = \int_{0}^{a_{c 1}} f_{p} (a) n_{1} (a) da + \int_{a_{c 1}}^{a_{12}} f_{e} (a) n_{1} (a) da \\ + \int_{a_{12}}^{a_{c 2}} f_{p} (a) n_{2} (a) da + \int_{a_{c 2}}^{a_{1}} f_{e} (a) n_{2} (a) da \\ K_{n_M_B_D} = \int_{a_{c 1}}^{a_{12}} k_{n} (a) n_{1} (a) da + \int_{a_{c 2}}^{a_{1}} k_{n} (a) n_{2} (a) da \\ C_{n_M_B_D} = \frac{\frac{1}{2} \int_{0}^{a_{c 1}} f_{p} (a) δ (a) n_{1} (a) da + \frac{1}{2} \int_{a_{12}}^{a_{c 2}} f_{p} (a) δ (a) n_{2} (a) da}{\frac{1}{2} \int_{a_{c 1}}^{a_{12}} f_{e} (a) δ (a) n_{1} (a) da + \frac{1}{2} \int_{a_{c 2}}^{a_{1}} f_{e} (a) δ (a) n_{2} (a) da} \\ \cdot \sqrt{m K_{n_M_B_D}} (6) \end{matrix}

(6)

Cases B and C: $a_{c 1} > a_{12} and a_{c 2} < a_{12}$

\begin{matrix} F = \int_{0}^{a_{12}} f_{p} (a) n_{1} (a) da + \int_{a_{12}}^{a_{l}} f_{e} (a) n_{2} (a) da \\ K_{n_M_B_D} = \int_{a_{12}}^{a_{l}} k_{n} (a) n_{2} (a) da \\ C_{n_M_B_D} = \frac{\frac{1}{2} \int_{0}^{a_{12}} f_{p} (a) δ (a) n_{1} (a) da}{\frac{1}{2} \int_{a_{12}}^{a_{1}} f_{e} (a) δ (a) n_{2} (a) da} \cdot \sqrt{m K_{n_M_B_D}} \end{matrix}

(7)

Cases B and D: $a_{c 1} > a_{12} and a_{c 2} > a_{12}$

\begin{matrix} F = \int_{0}^{a_{12}} f_{p} (a) n_{1} (a) da + \int_{a_{12}}^{a_{c 2}} f_{p} (a) n_{2} (a) da \\ + \int_{a_{c 2}}^{a_{l}} f_{e} (a) n_{2} (a) da \\ K_{n_M_B_D} = \int_{a_{c 2}}^{a_{l}} k_{n} (a) n_{2} (a) da \\ C_{n_M_B_D} = \frac{\frac{1}{2} \int_{0}^{a_{12}} f_{p} (a) δ (a) n_{1} (a) da + \frac{1}{2} \int_{a_{12}}^{a_{c 2}} f_{p} (a) δ (a) n_{2} (a) da}{\frac{1}{2} \int_{a_{c 2}}^{a_{l}} f_{e} (a) δ (a) n_{2} (a) da} \\ \cdot \sqrt{m K_{n_M_B_D}} (8) \end{matrix}

(8)

Surface fractal contact model

The G-W model and M-B fractal contact model were concerned only with the rough planes contact without involving rough surfaces. Specifically, these models only considered the impact of the contact surface microscopic morphology on the contact stiffness and damping without considering the impact of the macroscopic surface shape on the contact stiffness and damping.

When two rough surfaces come into contact, the actual contact area is less than the one of the same size rough plane due to the impact of the surface shapes. The surface impact coefficient was introduced to account for this phenomenon (see Figure 3)³⁴

λ = {(\frac{S_{h}}{\sum S})}^{X_{h}}

(9)

where $S_{h}$ is the theoretical contact area of the two surfaces which can be obtained according to different contact styles (see Table 1), $\sum S$ is the sum of the areas of the two surfaces, and $X_{h}$ is a synthetic curvature coefficient which is defined as

X_{h} = \frac{1}{R_{1}} + \frac{1}{R_{2}}

(10)

where $F$ is the normal load between the surfaces, $E$ is the equivalent elastic modulus, $R_{1}$ and $R_{2}$ are curvature radius of two surfaces, respectively. If one of the surfaces is plane, the curvature radius $R_{1} (R_{2}) \to \infty$

Table 1.

The theoretical contact area of two surfaces.

Contact style	Contact area
	$S_{h} = π {(\frac{3 FR}{4 E})}^{2 / 3}$ ; $1 / R = 1 / R_1 \pm 1 / R_2$ ,where ‘+’ indicates outer contact and ‘−’ indicates inner contact
	$S_{h} = 4 l {(\frac{FR}{π lE})}^{1 / 2}$ ; $1 / R = 1 / R_1 \pm 1 / R_2$ , where ‘+’ indicates outer contact and ‘−’ indicates inner contact

Figure 3.

Micro-model of the surface contact.

The variable $\sum S$ is the sum of the areas of the two surfaces

\sum S = S_{1} \pm S_{2}

(11)

where $S_{1}$ is the area of surface 1 and $S_{2}$ is the area of surface 2 (a ‘+’ indicates outer contact and a ‘−’ indicates inner contact).

The new contact area distribution function is defined as

n' (a) = λ \cdot n (a)

(12)

The real contact area is defined as $A'_{r} = \int_{0}^{a_{l}} n' (a) ada = λ A_{r}$ , where $n_{i} (a)$ in the bi-fractal M-B contact model should be replaced with $n'_{i} (a) (i = 1, 2)$ . Ultimately, the contact stiffness and contact damping of rough surfaces is defined as (subscript $n_M_B$ indicates the normal direction of curved surface M-B model)

\begin{matrix} K_{n_M_B_C} = λ K_{n_M_B_D} \\ C_{n_M_B_C} = \sqrt{λ} C_{n_M_B_D} \end{matrix}

(13)

Series fractal contact model

The G-W model and M-B fractal contact model only considered the contact stiffness and damping changes with the actual contact area. In fact, the contact body could be simplified to a smooth base with a rough surface connected in series (see Figure 4). During loading, the smooth base also incurs a slight deformation. The contact stiffness and damping of the smooth base are constant, whereas the contact stiffness and damping of a surface asperity changes non-linearly with any increase in the load.

Figure 4.

Simplified model of the contact surface in series.

The base stiffness is defined as

K_{n_base} = \frac{EA}{L}

(14)

where A is the cross-sectional area of contact body, and L is the thickness of contact body.

The asperity stiffness for a planar contact is defined as

K_{n_interface_p} = K_{n_M_B_D}

(15)

The asperity stiffness for a surface contact is defined as

K_{n_interface_c} = K_{n_M_B_C} = λ K_{n_M_B_D}

(16)

The contact stiffness can be related in parallel according to spring formula as follows

\frac{1}{K_{n_total}} = \frac{1}{K_{n_interface_p}} + \frac{1}{K_{n_base}}

(17)

The total stiffness for a planar contact is defined as

\frac{1}{K_{n_total}} = \frac{1}{K_{n_interface_p}} + \frac{1}{K_{n_base}}

(18)

The total stiffness for a surface contact is defined as

\frac{1}{K_{n_total}} = \frac{1}{K_{n_interface_c}} + \frac{1}{K_{n_base}}

(19)

Design of experiments

The principles of the experimental design in this study are

To guarantee that pressure is evenly distributed on the joint surfaces;

To minimise the measurement of physical quantities and facilitate more direct measurements;

To ensure high sensor resolution and high measurement accuracy and reliability;

To ensure that the dynamic characteristics of the joint surface can be extrapolated from the experimental system.

According to these principles, the investigated interface is formed by the square contact surfaces of $100 \times 100$ between the upper and the lower specimens without lubrication medium, which are produced by milling. The joint material is carbon steel with $E = 207 GPa$ , $μ = 0.29$ , and $HB = 200$ (HB: Brinell hardness). The static preload at the interface is applied by means of four bolts attached displacement sensor, and the exciting force is excited on the surface of the upper specimen through the exciter. A three-phase accelerometer is fixed on the upper specimen. The specimen is suspended by the rubber rope to prevent the response of the basement structure from influencing the test rig. The response signal was collected from the fixed three-phase accelerometer and displacement sensors (Figure 5). In the test process, it should be ensured as far as possible that the distribution of normal static load at the interface is uniform. So two specimens must be as small as possible in dimension and connected by multiple symmetry bolts. Four bolts must keep same preload length by controlling torque. The relative deformation of the two specimens was treated as the bolt length variation measured by displacement sensors (Figure 5(b)). Therefore, we can directly measure the amount of bolt deformation to obtain the overall relative dynamic displacement and relative speed of the joint surfaces. It is easier and more precise than a method which obtains the overall relative dynamic displacement and speed of joint surfaces by subtracting the dynamic displacement and speed of specimens 1 and 2. Bolt deformation is an average amount of four displacement sensors. The stiffness in the four bolts (M12), $4 k_{b}$ , will influence the normal overall contact stiffness and damping of joint surface as $K_{n}$ and $C_{n}$ , respectively, where $k_{b} = EA / l$ , $E = 207 GPa$ is the elastic module of bolt, $A = π d^{2} / 4$ (d is the minimum diameter of bolt) is the sectional area of bolt, and $l$ is the length of bolt. The deformation and stiffness of the preloaded bolt are shown in Figure 5(b). The mathematical model of the experiment³⁵ is (see Figure 5(c))

\begin{matrix} m_{1} {\overset{\cdot\cdot}{x}}_{1} (t) + f = f_{n} \\ f_{n} = (4 k_{b} + K_{n}) \cdot Δ x + C_{n} \cdot Δ \overset{\cdot}{x} \end{matrix}

(20)

The variable ${\overset{\cdot\cdot}{x}}_{1}$ is the normal acceleration of specimen 1, $Δ x$ is the normal dynamic displacement of the joint surface, $f$ is the applied force on specimen 1, $f_{n}$ is the normal dynamic load on the joint surface, $m_{1}$ is the mass of specimen 1, $k_{b}$ is the stiffness of a single bolt, $K_{n}$ is the overall normal stiffness of the joint surface, and $C_{n}$ is the overall normal damping on the joint surface.

Figure 5.

Experimental program: (a) the experiment of normal stiffness and damping, (b) the stiffness of bolt and (c) the model of stiffness and damping of joint surface.

In the experiment, the specimens are excited by a harmonic signal. The dynamic load on the joint surface can be obtained by measuring the acceleration signal from specimen 1

f_{n} = F_{n} \cos ω t = f - m_{1} {\overset{\cdot\cdot}{x}}_{1}

(21)

The average linear displacement of the bolts (the average relative linear displacement of the joint surface) is defined as

Δ x = Xcos (ω t + ϕ_{n})

(22)

The average linear velocity of the bolts (the average relative linear velocity of the joint surface) is defined as

Δ \overset{\cdot}{x} = - X ω \sin (ω t + ϕ_{n})

(23)

According to the principles established for this experiment

f_{n} = (4 k_{b} + K_{n}) \cdot Δ x + C_{n} \cdot Δ \overset{\cdot}{x} = f - m_{1} {\overset{\cdot\cdot}{x}}_{1} (t) = F_{n} \cos ω t

(24)

and

4 k_{b} + K_{n} = \frac{F_{n} \cos ϕ_{n}}{X_{n}}

(25)

C_{n} = \frac{F_{n} \sin ϕ_{n}}{ω X_{n}}

(26)

Therefore, the contact stiffness and damping per unit area are

k_{n} = (\frac{F_{n} \cos ϕ_{n}}{X_{n}} - 4 k_{b}) / S

(27)

C_{n} = \frac{F_{n} \sin ϕ_{n}}{ω X_{n} S}

(28)

Preparations for experiment

Before performing experiments, the specimens of different roughnesses should be milled and measured. Then, we should calculate the fractal parameters of specimen interfaces and choose the specified specimens as experimental objects. We applied the structure function method to calculate the fractal dimension D and fractal roughness G of surface profiles. The W-M function³⁶ is continuous and non-differentiable everywhere and takes on a statistical self-affinity, which can be used to simulate the rough surface morphology as follows^37,38

Z (x) = G^{(D - 1)} \sum_{n = n_{1}}^{\infty} \frac{\cos (2 π γ^{n} x)}{γ^{(2 - D) n}} 1 < D < 2, γ > 1

(29)

where $γ = 1.5$ can characterise the randomness of high spectral density and phase. The structure function $S (τ)$ of the W-M function can be expressed in the form

S (τ) = 〈 {[z (x) - z (x + τ)]}^{2} 〉 = C G^{2 (D - 1)} τ^{4 - 2 D}

(30)

where $τ$ is the displacement along the x-direction, $〈〉$ denotes averaging over the statistical ensemble of $z (x)$ , and $C$ is the scaling coefficient given by the formula

C = \frac{Γ (2 D - 3) \sin [(2 D - 3) π / 2]}{(4 - 2 D) \ln γ}

(31)

The plot of $S (τ)$ as a function of $τ$ in equation (30) is a straight line on double logarithmic coordinates

\log S (τ) = (4 - 2 D) \log τ + \log C + 2 (D - 1) \log G

(32)

Therefore, the fractal parameters D and G can be obtained from the structure function graph of the surfaces. It is seen from equation (32) that the fractal dimension D can reflect the slope of the straight line and fractal roughness G can determine the intercept of the straight line³⁹

D = \frac{(4 - k_{τ})}{2}

(33)

\log G = \frac{h_{τ} - \log C}{2 (D - 1)}

(34)

The analogue signals of specimen contour $Z (x)$ measured by contourgraph are inputted computer through analogue-to-digital (A/D) conversion. The N contour points are sampled at $Δ t$ sampling interval and the height is denoted as

Z (x_{i}) = Z_{i} (i = 0, 1, 2, 3, \dots, N - 1)

(35)

Set $τ = n Δ t$ , then, the structure function is

\begin{matrix} S (τ) = 〈 {[Z (x + n Δ t) - Z (x)]}^{2} 〉 \\ = \frac{1}{N - n} \sum_{i = 0}^{N - n} {(Z_{i + n} - Z_{i})}^{2} \end{matrix}

(36)

Plot the double logarithmic graphs of structure function of specimen interface and measure the slope and intercept of the straight line, then we can obtain the fractal dimension D and fractal roughness G of contact surfaces.

Results

Theoretical calculations

Series fractal contact model calculations

The contact stiffness of the M-B fractal model with D and G change is as shown in Figure 6. When a single-fractal surface considering the series structure of asperities with base and a rigid smooth plane are in contact (see section ‘Series fractal contact model’), the contact stiffness and damping curves determined from D and G change, are shown in Figures 7 and 8. The parameters of the contact bodies are listed in Table 2.

Figure 6.

Change of contact stiffness of the fractal model with D and G.

Figure 7.

Change of contact stiffness of the series fractal model with D and G.

Figure 8.

Change of contact damping of the series fractal model with D and G.

Table 2.

Contact body attributions.

E (N/mm²)	$μ$	H (HB)	Geometry (mm × mm × mm)
2.07E+5	0.29	200	100 ×100 × 20

For Figure 6, the contact stiffness is linearly increased as the load increases all along, and the greater the value of D is and the smaller the value of G is, the larger the contact stiffness is. For Figures 7 and 8, the contact stiffness increases non-linearly as the load increases and the increase is initially quick and then slows down and stabilises. The load becomes smaller and the contact stiffness stabilises quicker for larger values of D and lower values of G. The contact damping changes non-linearly with increasing load. The fractal dimension value D = 1.42 is a critical value. When D > 1.42, the contact damping increases with the load, and the contact damping increases more significantly for greater values of D. When D < 1.42, the contact damping decreases as the load increases, and the contact damping decreases more slowly for smaller values of D. As the value of G decreases, the contact damping value decreases.

Bi-fractal contact model calculations

The curves for the normal contact stiffness and damping for a bi-fractal rough plane contact as a function of D and G are shown in Figure 9 (see section ‘Bi-fractal contact model’). The parameters of the contact bodies are listed in Table 2. Two single-fractal rough plane contacts can be simplified to a contact between bi-fractal rough surface and a rigid plane, which can be calculated using the M-B single-fractal model when the value of D is the same. When the values of D for two single-fractal surfaces are different, four cases arise according to the bi-fractal theory.

Figure 9.

Normal contact stiffness and damping of the bi-fractal model.

For Cases B and D: $a_{c 1} > a_{12}$ and $a_{c 2} > a_{12}$ , the parameter G = 1E−10 is selected to be a constant value and D is selected to be 1.48, 1.44, 1.4, and 1.36 to calculate a mapping relationship between the contact stiffness–load and contact damping–load of the single-fractal surfaces with a rigid plane in contact. Bi-fractal theory is also used to obtain the curves for the contact stiffness–load and the contact damping–load of the two single-fractal surfaces in contact. As shown in Figure 10, when the two single-fractal surfaces with different D values contact, the contact stiffness–load curve is smaller than the stiffness–load curve for two single-fractal surfaces and a rigid plane contact, and the curve is similar to that of a smaller D value. The contact damping–load curve is located in the middle of the curve of two single-fractal surfaces and a rigid plane in contact and is similar to the contact damping–load curve for a smaller D value.

Figure 10.

Bi-fractal model simplification.

Surface fractal contact calculations

The variation curve of the normal contact stiffness and damping of rough cylinders in contact with the change in the curvature radius is shown in Figure 12 (see section ‘Surface fractal contact model’). The parameters of the contact bodies are listed in Table 2. The stiffness–load curve of the rough plane contact, rough cylindrical outer contact, and rough cylindrical inner contact is shown and compared.

Experimental test

Four bolts (M12) increase the pressure at 0.5-MPa intervals through different preloads, and every time, the bolts’ preloaded length has to be the same by controlling the torque. The contact stiffness and damping of the joint surface were measured at different pressures, and the theoretical and experimental values are compared in Table 3 and Figures 12 and 13.

Table 3.

Comparison of the experimental and theoretical results.

Pressure of the joint surface (MPa)	Theoretical value of the normal stiffness of the joint surface (N/mm³)	Experimental value of the normal stiffness of the joint surface (N/mm³)	Theoretical value of the normal damping of the joint surface (N s/mm³)	Experimental value of the normal damping of the joint surface (N s/mm³)
0	0	0	1.250E−5	1.275E−5
0.5	2358	998	0.715E−5	0.845E−5
1	3623	2920	0.685E−5	0.986E−5
1.5	4519	5324	0.678E−5	0.784E−5
2	5194	6185	0.650E−5	0.796E−5
2.5	5728	6076	0.646E−5	0.813E−5
3	6153	6425	0.638E−5	0.765E−5
3.5	6511	6988	0.633E−5	0.778E−5
4	6814	6424	0.626E−5	0.788E−5
4.5	7075	6600	0.623E−5	0.748E−5
5	7301	6420	0.618E−5	0.756E−5

The theoretical and experimental values obtained for the stiffness and damping between the two rough and milled planes are compared in Figure 12. The parameters of the contact bodies are listed in Table 2. The experimental contact stiffness increases with the load increase and remains unchanged after a certain load is reached, which is all along fluctuated around the theoretical values. The experimental results for the contact damping decrease as the load increases; the increase is initially quick and then slows down and stabilises. The experimental results are always greater than the theoretical results.

The variations and the comparison curve of the contact stiffness and damping are shown in Figure 13 by measuring a milled rough plane (D = 1.58) and comparing these values with three other milled, contacting, rough planes (D = 1.35, 1.41, and 1.44). The parameters of the contact bodies are listed in Table 2. The contact stiffness increases with the load. The contact stiffness will be lower with rougher surfaces, and the load will change slightly as the contact damping increases or decreases. The experimental and theoretical trends are the same, but there are numerical differences, which can be used to predict the trend in machine dynamics.

Discussion

Effects of D and G on the stiffness and damping of series fractal model

For Figure 6, the contact stiffness of the M-B fractal model is non-linearly increased all along as the load increases, and the greater the value of D is and the smaller the value of G is, the larger the contact stiffness is. The contact stiffness and the real contact area are subjected to the index relationship $K_{n_M_B} ~ A_{r}^{1 / 2}$ . With the load increase, the real contact area of asperities is increased non-linearly, so the contact stiffness is increased all along without stable value.

For Figures 7 and 8, the contact stiffness of series fractal model increases non-linearly as the load increases; the increase is initially quick and then slows down and stabilises. The joint surface is formed by stacking the asperities on the base. As the load increases, both the number of asperities in contact and the area of the asperities increase simultaneously. If every asperity is considered to be a spring, the joint surface can be considered a parallel structure of springs. As the number of springs in the parallel and single spring stiffness increases, the contact stiffness at the joint surface increases non-linearly with increasing load. At first, the number and area of the asperities involved at the joint contact increase rapidly. As the load continues to increase, the rate of the increase in the number of asperities and the total asperity area slows down such that the contact stiffness becomes more rigid and then slows down. When the load continues to increase, the number of asperities in contact begins to decrease, whereas the asperity contact area continues to increase. When equilibrium is established, the contact stiffness stabilises and the actual contact area is nearly equal to the nominal contact area. The joint surface is pressed flat, and the contact stiffness is approximately equivalent to the tensile stiffness of the base.

A greater value of D and a smaller value of G lead to a smaller load and smaller contact stiffness stabilisation time. A greater value of D and a smaller value of G lead to better consistency in the height and area of the asperities. As the load increases, the number of asperities in contact at the same time causes rapid deformation due to the height consistency, which results in an actual contact area reaches the nominal contact area soon after and the joint surface is pressed flat. In this case, the contact stiffness stabilises very quickly and the large number of small asperities can withstand only a small load before deforming.

When D > 1.42, the contact damping increases with the load, and the contact damping increases more quickly for greater values of D. When D < 1.42, the contact damping decreases with the load, and the contact damping decreases more slowly for lower values of D. The contact damping $C_{n} = (ω_{p} / ω_{e}) {(m K_{n_M_B})}^{1 / 2}$ is used to measure the energy consumption capacity of the structure. There are not rigorous interpretations for the increase or decrease of damping with load. In fact, the D = 1.42 is just an estimated critical value.

Comparison of single-fractal model and bi-fractal model

For Figures 9 and 10, the contact stiffness–load curve of bi-fractal model is smaller than any one of the two single-fractal models and is closer to the contact stiffness–load curve of single-fractal model with a smaller value of D. The two single-fractal surface contacts are regarded as two springs in series, and the composite stiffness of the two springs in series is smaller than any one of the two spring stiffnesses. Two rough surface contacts can be simplified as a rough surface contacting with a rigid plane, similar to the G-W and M-B fractal models. The joint surface of the smoother contact (the larger D) is seen as a rigid smooth plane, and the joint surface of the rougher contact (the smaller D) is regarded as an equivalent rough surface; thus, the curve behaves similarly to a joint contact with a smaller D value.

The contact damping–load curve of bi-fractal model is located in the middle of the curves for the two single-fractal models and is closer to the contact damping–load curve with a smaller D value. The damping is a parameter that measures energy consumption. When two single-fractal surfaces come into contact, the energy consumption should be a weighted average so that the damping curve is located between the two single-fractal surface curves. The damping curve is similar to that of a smaller value of D because two single-fractal surface contacts are seen as a rough surface (the smaller D) and a rigid plane (the larger D) in contact.

Effects of curvature radius on the stiffness of surface fractal model

For Figure 11, the curvature radius and contact style play an important role on the contact stiffness values of joint surfaces but have no influence on the change trend of them. We can explain as follows: According to the contact stiffness formula for Hertz theory, $k_{n} = 4 / 3 E R^{1 / 2} δ^{1 / 2}$ , for the same deformation $δ$ , the synthetic contact radius is $R = R_{1} R_{2} / (R_{1} \pm R_{2}) = R_{1} / (R_{1} / R_{2} \pm 1)$ .

Figure 11.

Contact stiffness of the surface fractal model.

If $R_{1} > R_{2}$ , a ‘+’ indicates outer contact and a ‘−’ indicates inner contact. The synthetic contact radius of the inner contact is clearly greater than that for the outer contact. Therefore, the inner contact stiffness is greater than the outer contact stiffness. Suppose that contact body 1 is a rough plane, namely $R_{1} \to \infty$ and that contact body 2 is a rough cylinder of radius $R_{2}$ . When $R_{2} \to \infty$ , the synthetic contact radius increases, which causes the rough plane contact stiffness to be greater than the rough cylindrical contact stiffness. For inner contact, a value of $R_{1} = 100 mm$ is assumed. The variable $R_{2}$ is then set to 40, 70, 100, and 130 mm. Because of $R = R_{1} R_{2} / (R_{1} \pm R_{2}) = R_{1} / (R_{1} / R_{2} \pm 1)$ , as $R_{2}$ increases, the synthetic contact radius, and thus the contact stiffness increase. When $R_{2} = R_{1}$ , the radii are equivalent to that of a rough plane, and the rough plane contact and contact stiffness reach a maximum value. When $R_{2} > R_{1}$ , the rough plane contact returns to a rough cylindrical inner contact, causing the synthetic contact radius and contact stiffness to decrease.

Comparison of theoretical calculation and experimental test

For Figure 12(a), the main differences between the experiment curve and theoretical curve of stiffness–load are the slope of them. The slope of theoretical curve has been reduced until zero, while one of the experimental curves first increases and then decreases to zero finally. The experimental curve has experienced a process of less than, greater than, and less than the theoretical curve. The height of asperities at interface is commonly subject to the Gauss distribution. In the beginning, the number of asperities in contact and the single asperity area are relatively small, and most of asperities are deformed plastically, so the actual contact stiffness is less than the theoretical values. With increasing load, the number of asperities in contact and asperity area rapidly increase simultaneously, which results in the actual contact stiffness and rate of stiffness change exceed the theoretical values. When a certain load is reached, the number of asperities in contact begins to reduce, while the asperity area still continues to increase. When equilibrium is established, the contact stiffness keeps stable. In addition, due to the stress hardening of asperities, the actual contact stiffness is less than the theoretical contact stiffness.

Figure 12.

For Figure 12(b), the shape and trend of the theoretical curve and the experimental curve about damping–load is basically accordant, and the experimental curve is overall greater than the theoretical curve. The main reason is that the energy loss of plastic deformation in theoretical calculation is only considered without involving the else energy loss due to contact creep and friction.

For Figure 13, the experimental curve of stiffness–load is similar to the theoretical curve, while the damping–load trend is not clear when the load is small. Mainly because the number of asperities in contact is small at present and not only exist plastic deformation and elastic deformation, but also a variety of forms of energy loss, it is unable to determine the dominant form of the damping.

Figure 13.

The measured results of stiffness and damping with different Ds.

Limitations

The limitations of the presented model need to also be mentioned: (1) The non-normalised formulas will only guarantee the consistency with the trends without stable results. (2) The mechanism for the increase or decrease of damping with load has not been clear. The critical value of D determined the increase or decrease has not been given strictly. (3) The accurate relationship between roughness and fractal parameters has not been presented rigorously. (4) The interfaces in contact joint by bolt have been viewed as experimental object to verify the presented model without involving other typical joint structure, such as gear transmission and bearing support.

Conclusion

This article summarises the characteristics and limitations of the Hertz, G-W, and M-B contact models and proposes the three improvements based on the M-B model. The conclusions can be summed up as follows:

The contact stiffness of the M-B fractal model increases non-linearly all along without stable value, while the contact stiffness of series fractal model increases non-linearly with the load, and the increase is initially quick and then slows down and stabilises. D = 1.42 is as an estimated critical value which determined the contact damping increase or decrease with load.

The contact stiffness–load curve of bi-fractal model is smaller than any one of the two single-fractal models and is closer to the contact stiffness–load curve of single-fractal model with a smaller value of D. The contact damping–load curve of bi-fractal model is located in the middle of the curves of the two single-fractal models and is closer to the contact damping–load curve with a smaller D value.

The curvature radius and contact style play an important role on the contact stiffness values of joint surfaces but have no influence on the change trend of them.

The theoretical calculations have been verified by the experimental test. The improved theoretical model can match the experimental measured data well.

Footnotes

Appendix 1 Declaration of conflicting interests

The authors declare that there is no conflict of interest.

Funding

This work was supported by Department Pre-Research Fund of China, project number 9140A18010312KG 0101.

References

Xueliang

. The dynamic characteristic and application of mechanical joint surface. Beijing: Press of Science and Technology of China, 2002.

Ren

Bears

. Identification of effective linear joints using coupling and joint identification techniques. J Vib Acoust 1998; 120: 331–338.

Majumdar

Bhushan

. Role of fractal geometry in roughness characterization and contact mechanics of surfaces. J Tribol: T ASME 1990; 112: 205–216.

Wan

Gong

. Research on dynamic alignment for sophisticated system based on ACP technology. Adv Mat Res 2013; 630: 473–478.

Hertz

. Über die Berührung fester, elastischer Körper. J Reine Angew Math 1881; 92: 156–171.

Timoshenkl

Goodier

. Theory of elasticity. New York: McGraw-Hill, 1951.

Johnson

. Contact mechanics. Cambridge: Cambridge University Press, 1987.

Greenwood

Williamson

JBP

. Contact of nominally flat surfaces. P R Soc Lond A Math 1966; 295: 300–319.

Majumdar

Bhushan

. Fractal model of elastic-plastic contact between rough surfaces. J Tribol: T ASME 1991; 113: 1–11.

10.

Pullen

Williamson

JBP

. On the plastic contact of rough surfaces. P R Soc Lond A Math 1972; 327: 157–173.

11.

Powierza

Klimczak

Polijaniuk

. On the experimental verification of the Greenwood–Williamson model for the contact of rough surfaces. Wear 1992; 154: 115–124.

12.

Greenwood

Tripp

. The elastic contact of rough spheres. J Appl Mech: T ASME 1967; 34: 153–159.

13.

Chang

Etsion

Bogy

. An elastic–plastic model for the contact of rough surfaces. J Tribol 1987; 110: 50–56.

14.

McCool

Gassel

. The contact of two surfaces having anisotropic roughness geometry. ASLE Spec Publ SP 1981; 7: 29–38.

15.

Halling

Nuri

. The elastic–plastic contact of rough surfaces and its relevant in study of wear. P I Mech Eng C: J Mec 1988; 202(C4): 269–274.

16.

Stanley

Etison

Bogy

. Adhesion of contacting rough surfaces in the presence of sub-boundary. J Tribol 1990; 112: 98–104.

17.

Burdekin

Back

Cowley

. Analysis of the local deformation in machine joints. 1978.

18.

Yamada

Kakubari

. Prediction of dynamic characteristics of beam containing junction. 1986.

19.

Jinmin

Tianning

. Fractal model for normal dynamic parameters of joints surfaces. J Xi’an Jiaotong Univ 2009; 43: 91–94.

20.

Majumdar

Tien

. Fractal characterization and simulation of rough surfaces. Wear 1990; 136: 313–327.

21.

Yan

Komvopoulos

. Contact analysis of elastic-plastic fractal surfaces. J Appl Phys 1998; 84: 3617–3624.

22.

Ciavarella

Demelio

Barber

. Linear elastic contact of the Weierstrass profile. Proc R Soc Lond: A 2000; 456: 387–405.

23.

Ciavarella

Demelio

. Elastic multiscale contact of rough surfaces: Archard’s model revisited and comparison with modern fractal models. J Appl Mech: T ASME 2001; 68: 496–498.

24.

Warren

Krajcinovic

. Random Cantor set models for the elastic-perfectly plastic contact of rough surfaces. Wear 1996; 196: 1–15.

25.

Ciavarella

Murolo

Demelio

. Elastic contact stiffness and contact resistance for the Weierstrass profile. J Mech Phys Solids 2004; 52: 1247–1265.

26.

Komvopoulos

. Three-dimensional contact analysis of elastic-plastic layered media with fractal surface topographies. J Tribol: T ASME 2001; 123: 632–640.

27.

Zhou

Leu

Blackmore

. Fractal geometry model for wear prediction. Wear 1993; 170: 1–14.

28.

Zhu

. The fractal character of processed metal surfaces. Wear 1997; 208: 17–24.

29.

Brown

Savary

. Describing ground surface texture using, contact profilometry and fractal analysis. Wear 1991; 141: 211–226.

30.

Tricot

Ferland

Baran

. Fractal analysis of worn surfaces. Wear 1994; 172: 127–133.

31.

Majumdar

Bhushan

. Elastic-plastic contact model for bifractal surface. Wear 1992; 153: 53–64.

32.

Gao

. Bifractal characteristics of ground workpiece surfaces. Key Eng Mat 2001; 202–203: 25–28.

33.

Bhushan

. Nanotribology and nanomechanics an introduction. 2nd ed. Heidelberg: Springer, 2008.

34.

Kang

Han

. Research on fractal contact model of cylinders’ surface (The 5th China international symposium on tribology). Adv Tribol 2009; 11: 90–93.

35.

Konowalski

. Experimental research and modeling of normal contact stiffness and contact damping of machined joint surface. Adv Manuf Sci Technol 2009; 33: 53–68.

36.

Mandelbrot

. The fractal geometry of nature. San Francisco, CA: Freeman, 1982.

37.

Bhushan

. Introduction to tribology[M]. New York: John Wiley & Sons, 2002.

38.

Berry

Lewis

. On the Weierstrass-Mandelbrot fractal function. P Roy Soc Lond Ser: A 1980; 370: 459–484.

39.

Jiang

Zheng

Zhu

. A contact stiffness model of machined plane joint based on fractal theory [J]. J Tribol: T ASME 2010; 132: 011401.

An improved algorithm for the normal contact stiffness and damping of a mechanical joint surface

Abstract

Keywords

Introduction

Methods

Algorithm improvements

Bi-fractal contact model

Surface fractal contact model

Series fractal contact model

Design of experiments

Preparations for experiment

Results

Theoretical calculations

Series fractal contact model calculations

Bi-fractal contact model calculations

Surface fractal contact calculations

Experimental test

Discussion

Effects of D and G on the stiffness and damping of series fractal model

Comparison of single-fractal model and bi-fractal model

Effects of curvature radius on the stiffness of surface fractal model

Comparison of theoretical calculation and experimental test

Limitations

Conclusion

Footnotes

Appendix 1

Declaration of conflicting interests

Funding

References