Temperature in the Friction Couple Consisting of Functionally Graded and Homogeneous Materials

Abstract

An analytical model was developed to determine the temperature of friction coupling, in which one element was made of a functionally graded material (FGM) and the other was homogeneous. First, for such a system, the boundary–value problem of heat conduction was formulated with consideration of the heat generation due to friction. Then, using the Laplace integral transform, an exact solution to this problem was obtained for uniform sliding, and braking with constant deceleration. A numerical analysis was performed for the selected friction pair consisting of the FGM (zircon dioxide + titanium alloy) and cast iron. It was established that the use of elements made of a FGM consisting of ZrO₂ and Ti-6Al-4V can significantly reduce the maximum temperature achieved in the friction system.

1. Introduction

Reviews of investigations on methods for establishing the temperature of systems containing friction elements made of functionally gradient materials (FGMs) can be found in previous articles [1,2,3]. In these studies, the methodology of determining the temperature in such friction couples under uniform sliding [1], during braking with time-dependent contact pressure [2], and considering the thermal sensitivity of component materials of FGMs was investigated [3]. The main factor in this methodology is an exact solution to the boundary–value heat conduction problem, taking into account the frictional heating of two semi-infinite bodies made of FGMs. It should be noted, however, that the obtained solutions did not allow determining automatically, with the help of limit transformations, solutions to the problems in the case when one of the friction pair elements is made of FGM and the other is homogeneous. Moreover, this type of friction pair is one of the most common [4]. Therefore, in this study, an attempt was made to develop a mathematical model for determining the temperature of a friction pair consisting of a body made of a two-component FGM, sliding on the surface of a homogeneous body. An exponential change in the thermal conductivity of the FGM with distance from the friction surface was assumed. Two modes of changing the sliding velocity over time were considered: uniform and linearly decreasing.

2. Statement of the Problem

The object of study is the transient temperature field, initiated in the process of frictional heating of the friction pair elements of a braking system, corresponding to the brake pad and disc. Taking into account the fact that the heat generated as a result of friction during braking is mainly directed along the normal from the friction surface to the inside of both elements [5,6], for the description of the heating process of the system, a contact scheme of two semi-infinite bodies was adopted, related to the Cartesian coordinate system (Figure 1).

The pad (body 1) is made of a two-component functionally graded material (FGM), in such a way that the friction surface is a material of low thermal conductivity and high wear resistance (ceramics etc.), while the core material has high thermal conductivity (metal alloys, copper, iron etc.). The increase in thermal conductivity of the pad material in the distance from the friction surface is exponential. On the other hand, the disc (body 2) is made of a homogeneous material (cast iron etc.). A more detailed description of the adopted model assumptions is presented in our previous articles [1,2].

The analytical model presented in the manuscript concerns the frictional system of two semi-infinite bodies, in which it is not possible to take into consideration the heat exchange between the heated elements and the surrounding environment. It is known, however, that consideration of convection cooling, would lead to a lower maximum temperature; the most important parameters in the design process of frictional systems. It should be ensured that the theoretical value of the permissible temperature for a given material (i.e., the melting point) is not exceeded. For this reason, at the design stage, calculations should be performed for the maximum temperature achieved for adiabatic conditions on the free surfaces of the friction system.

The braking process with constant deceleration was considered when the contact pressure achieved its nominal value $p_0$ immediately at the beginning of the braking, with simultaneously reduction of velocity from the initial value $V_0$ to zero at the stopping moment $t=t_s$. For such braking, the specific friction power was written in the form:

$$q(t)=q_0{q}^{\ast }(t),q_0=f_0{p}_0{V}_0,q^{\ast }(t)=1-\frac{t}{t_s},0\le t\le t_s,t_s=\frac{W_0}{q_0{A}_a},$$

(1)

where $f_0$—friction coefficient, $W_0$—initial kinetic energy of the system, and $A_a$— area of nominal contact between one brake pad and the disc.

The temperature field $T(z,t)$ in the system consisting of two sliding semi-spaces was sought based on the solution to the following thermal problem of friction:

$$\frac{\partial }{\partial z}\left[K_1(z)\frac{\partial T(z,t)}{\partial z}\right]=c_1{\rho }_1\frac{\partial T(z,t)}{\partial t},z>0,\hspace{0.17em}0<t\le t_s,$$

(2)

$$K_2\frac{{\partial }^{2}T(z,t)}{\partial z^2}=c_2{\rho }_2\frac{\partial T(z,t)}{\partial t},z<0,\hspace{0.17em}0<t\le t_s,$$

(3)

$$T(0^{+},t)=T(0^{-},t)\equiv T(t),0<t\le t_s,$$

(4)

$${K_2\frac{\partial T(z,t)}{\partial z}|}_{z=0^{-}}-{K_1(z)\frac{\partial T(z,t)}{\partial z}|}_{z=0^{+}}=q(t),0<t\le t_s,$$

(5)

$$T(z,t)\to T_0,\left|z\right|\to \infty ,0<t\le t_s,$$

(6)

$$T(z,0)=T_0,\left|z\right|<\infty .$$

(7)

where

$$K_1(z)=K_{1,1}{e}^{\gamma z},z\ge 0,\gamma \ge 0,$$

(8)

$$c_1=c_{1,1}(1-v)+c_{1,2}v,{\rho }_1={\rho }_{1,1}(1-v)+{\rho }_{1,2}v,0\le v\le 1,$$

(9)

temporal profile of specific friction power $q(t)$ was determined from Equation (1), $K_{1,m}$, $c_{1,m}$, and ${\rho }_{1,m}$—thermal conductivity, specific heat, and density of the first ($m=1$) and the second ($m=2$) component of pad material, respectively, and parameters $K_2$, $c_2$, and ${\rho }_2$—correspond to the disc material, $v$—the relative volumetric fraction of the first component of the pad material, and $T_0$—temperature of the system at the initial time moment $t=0$.

The dimensionless variables and parameters were introduced:

$$\zeta =\frac{z}{a},\tau =\frac{k_{1}t}{a^2},{\tau }_s=\frac{k_1{t}_s}{a^2},K^{\ast }=\frac{K_2}{K_{1,1}},k^{\ast }=\frac{k_2}{k_1},{\Theta }^{\ast }=\frac{T-T_0}{{\Theta }_0},{\Theta }_0=\frac{q_{0}a}{K_{1,1}^{(0)}},$$

(10)

where

$$a=\sqrt{3k_1{t}_s},$$

(11)

$$k_1=\frac{K_{1,1}}{c_1{\rho }_1},k_2=\frac{K_2}{c_2{\rho }_2}.$$

(12)

Taking into account the designations (10)–(12), the problem (2)–(9) was written in the form:

$$\frac{{\partial }^2{\Theta }^{\ast }(\zeta ,\tau )}{\partial {\zeta }^2}+{\gamma }^{\ast }\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \zeta }-e^{-{\gamma }^{\ast }\zeta }\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \tau }=0,\zeta >0,0<\tau \le {\tau }_s,$$

(13)

$$\frac{{\partial }^2{\Theta }^{\ast }(\zeta ,\tau )}{\partial {\zeta }^2}-\frac{1}{k^{\ast }}\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \tau }=0,\zeta <0,0<\tau \le {\tau }_s,$$

(14)

$${\Theta }^{\ast }(0^{+},\tau )={\Theta }^{\ast }(0^{-},\tau )\equiv {\Theta }^{\ast }(\tau ),0<\tau \le {\tau }_s,$$

(15)

$${K^{\ast }\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \zeta }|}_{\zeta =0^{-}}-{\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \zeta }|}_{\zeta =0^{+}}=q^{\ast }(\tau ),0<\tau \le {\tau }_s,$$

(16)

$${\Theta }^{\ast }(\zeta ,\tau )\to 0,\left|\zeta \right|\to \infty ,0<\tau \le {\tau }_s,$$

(17)

$${\Theta }^{\ast }(\zeta ,0)=0,\left|\zeta \right|<\infty ,$$

(18)

where

$$q^{\ast }(\tau )=1-\frac{\tau }{{\tau }_s},0<\tau \le {\tau }_s,$$

(19)

$${\gamma }^{\ast }\equiv a\gamma =\ln \left(\frac{K_{1,2}}{K_{1,1}}\right).$$

(20)

3. Frictional Heating under Uniform Sliding

First, the case of frictional heating process during sliding of the pad on the disc surface with constant velocity $V_0$ was considered. Then for ${\tau }_s\to \infty$ from the Equation (19), it follows that $q^{\ast }(\tau )=1$. For the boundary–value heat conduction problem (13)–(20) with a constant temporal profile of specific friction power $q^{\ast }(\tau )=1$, the integral Laplace transform was applied [7]:

$${\overline{\Theta }}^{\ast }(\zeta ,p)\equiv L[{\Theta }^{\ast }(\zeta ,\tau );\hspace{0.17em}p]={\displaystyle \underset{0}{\overset{\infty }{\int }}{\Theta }^{\ast }(\zeta ,\tau )e^{-p\tau }d\tau },$$

(21)

it was obtained:

$$\frac{d^2{\overline{\Theta }}^{\ast }(\zeta ,p)}{d{\zeta }^2}+{\gamma }^{\ast }\frac{d{\overline{\Theta }}^{\ast }(\zeta ,p)}{d\zeta }-p{e}^{-{\gamma }^{\ast }\zeta }{\overline{\Theta }}^{\ast }(\zeta ,p)=0,\zeta >0,$$

(22)

$$\frac{d^2{\overline{\Theta }}^{\ast }(\zeta ,p)}{d{\zeta }^2}-\frac{p}{k^{\ast }}{\overline{\Theta }}^{\ast }(\zeta ,p)=0,\zeta <0,$$

(23)

$${\overline{\Theta }}^{\ast }(0^{+},p)={\overline{\Theta }}^{\ast }(0^{-},p)\equiv {\overline{\Theta }}^{\ast }(p),$$

(24)

$${K^{\ast }\frac{d{\overline{\Theta }}^{\ast }(\zeta ,p)}{d\zeta }|}_{\zeta =0^{-}}-{\frac{d{\overline{\Theta }}^{\ast }(\zeta ,p)}{d\zeta }|}_{\zeta =0^{+}}=\frac{1}{p},$$

(25)

$${\overline{\Theta }}^{\ast }(\zeta ,p)\to 0,\left|\zeta \right|\to \infty .$$

(26)

An exact solution to the ordinary differential Equations (22) and (23), which meet the boundary conditions (24)–(26) has the form:

$${\overline{\Theta }}^{\ast }(\zeta ,p)=\frac{{\Delta }_1(\zeta ,p)}{p\sqrt{p\hspace{0.17em}}\Delta (p)},\zeta \ge 0,{\overline{\Theta }}^{\ast }(\zeta ,p)=\frac{{\Delta }_2(\zeta ,p)}{p\sqrt{p\hspace{0.17em}}\Delta (p)},\zeta \le 0,$$

(27)

where

$${\Delta }_1(\zeta ,p)=e^{-0.5{\gamma }^{\ast }\zeta }{I}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{p}{e}^{-0.5{\gamma }^{\ast }\zeta }\right),{\Delta }_2(\zeta ,p)=e^{\sqrt{\frac{p}{k^{\ast }}}\hspace{0.17em}\zeta }{I}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{p}\right),$$

(28)

$$\Delta (p)=I_0\left(\frac{2}{{\gamma }^{\ast }}\sqrt{p}\right)+K_{\epsilon }{I}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{p}\right),$$

(29)

$I_k(x)$—modified Bessel functions of the first kind of the kth order $k=0,\hspace{0.17em}1$ [8].

Using the inverse Laplace transform to the solution (27)–(29), the dimensional temperature rise was found in the form:

$${\Theta }^{\ast }(\zeta ,\tau )\equiv L^{-1}[{\overline{\Theta }}^{\ast }(\zeta ,p);\hspace{0.17em}\tau ]=\frac{1}{2\pi i}{\displaystyle \underset{\omega -i\infty }{\overset{\omega +i\infty }{\int }}{\overline{\Theta }}^{\ast }(\zeta ,p)e^{p\tau }dp},\tau \ge 0,\omega \equiv \mathrm{Re}p>0,i\equiv \sqrt{-1}.$$

(30)

The presence of $\sqrt{p}$, as well as the lack of the roots of function $\Delta (p)$, testifies that the solution (36)–(39) has a branch point for $p=0$. Therefore, to perform the integration on the complex plane ($\mathrm{Re}p,\hspace{0.17em}\mathrm{Im}\hspace{0.17em}p$), the closed contour $\Gamma$ was chosen, as demonstrated in Figure 2. The contour $\Gamma$ consists of the straight line ${\Gamma }_{\omega }$ $\mathrm{Re}\hspace{0.17em}p=\omega$, the circles ${\Gamma }_R$ and ${\Gamma }_{\delta }$ with the radii $R$ and $\delta$, respectively, with the center $p=0$, and a cut of a complex p–plane along negative real axis $\mathrm{Re}\hspace{0.17em}p<0$ and two boundaries ${\Gamma }_{\pm }$. Within the contour $\Gamma$, the integral function ${\overline{\Theta }}^{\ast }(\zeta ,p)$ in the Equation (30) is unambiguous and analytical.

Then, based on Cauchy’s theorem we obtained [9]:

$$\frac{1}{2\pi i}{\displaystyle \underset{\Gamma }{\oint }{\overline{\Theta }}^{\ast }}(\zeta ,p)e^{p\tau }dp=0.$$

(31)

Since the transform ${\overline{\Theta }}^{\ast }(\zeta ,\tau )$ carries out the conditions of Jordan’s lemma [7]:

$$\left|\frac{{\Delta }_l(\zeta ,p)}{p\sqrt{p}\Delta (p)}\right|\le \frac{const.}{p\sqrt{p}},l=1,2,$$

(32)

integrands on arcs ${\Gamma }_R$ in the Equation (31) tend to zero for $R\to \infty$; therefore, on the basis of the relations (30) and (31), the dimensional temperature rise was written in the form:

$${\Theta }^{\ast }(\zeta ,\tau )+{\Theta }_{+}^{\ast }(\zeta ,\tau )+{\Theta }_{-}^{\ast }(\zeta ,\tau )+{\Theta }_{\delta }^{\ast }(\zeta ,\tau )=0,\left|\zeta \right|<\infty ,\tau \ge 0,$$

(33)

where

$${\Theta }_{\pm }^{\ast }(\zeta ,\tau )=\frac{1}{2\pi i}{\displaystyle \underset{{\Gamma }_{\pm }}{\int }{\overline{\Theta }}^{\ast }(\zeta ,p)}{e}^{p\tau }dp,{\Theta }_{\delta }^{\ast }(\zeta ,\tau )=\frac{1}{2\pi i}{\displaystyle \underset{{\Gamma }_{\delta }}{\int }{\overline{\Theta }}^{\ast }(\zeta ,p)}{e}^{p\tau }dp.$$

(34)

In the polar coordinate system $(r,\phi )$ with center in the point $p=0$, parameter of the Laplace transform $p=r{e}^{i\phi }$, $r\ge 0$, and $|\phi |\le \pi$. Then on the boundary ${\Gamma }_{+}$ we obtained $p=r{e}^{i\pi }=-r$, $\sqrt{p}=i\sqrt{r}$, and on the edge ${\Gamma }_{-}$, respectively, $p=r{e}^{-i\pi }=-r$, $\sqrt{p}=-i\sqrt{r}$ and the first two integrals (34) took the form:

$${\Theta }_{\pm }^{\ast }(\zeta ,\tau )=\pm \frac{1}{2\pi i}{\displaystyle \underset{0}{\overset{\infty }{\int }}{\overline{\Theta }}_{\pm }^{\ast }(\zeta ,r)e^{-r\tau }dr},\left|\zeta \right|<\infty ,\tau \ge 0,$$

(35)

where ${\overline{\Theta }}_{\pm }^{\ast }(\zeta ,r)\equiv {\overline{\Theta }}^{\ast }(\zeta ,r{e}^{\pm i\pi })$.

Taking into account the dependencies [8]:

$$I_0(x)=J_0(ix),I_1(x)=-i{J}_1(ix),$$

(36)

(where $J_k(x)$ are the Bessel functions of the first kind of the kth order $k=0,1$), from Equations (27)–(29) was obtained:

$${\overline{\Theta }}_{\pm }^{\ast }(\zeta ,r)=\frac{{\Delta }_1^{\pm }(\zeta ,r)}{r\sqrt{r}\hspace{0.17em}{\Delta }^{\mp }(r)},\zeta \ge 0,{\overline{\Theta }}_{\pm }^{\ast }(\zeta ,r)=\frac{{\Delta }_2^{\pm }(\zeta ,r)}{r\sqrt{r}\hspace{0.17em}{\Delta }^{\mp }(r)},\zeta \le 0,$$

(37)

where:

$${\Delta }_1^{\pm }(\zeta ,r)=\pm i{e}^{-0.5{\gamma }^{\ast }\zeta }{J}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}{e}^{-0.5{\gamma }^{\ast }\zeta }\right),{\Delta }_2^{\pm }(\zeta ,r)=\pm i{e}^{\pm \hspace{0.17em}i\sqrt{\frac{r}{k^{\ast }}}\zeta }{J}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}\right),$$

(38)

$${\Delta }^{\pm }(\zeta ,r)=K_{\epsilon }{J}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}\right)\pm i{J}_0\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}\right).$$

(39)

On the circle ${\Gamma }_{\delta }$ it is $p=\delta e^{i\phi }$, $\sqrt{p}=\sqrt{\delta }{e}^{0.5i\phi }$, $|\phi |\le \pi$. Approaching the limit $\delta \to 0$ with consideration of the solutions forms (27)–(29), the third integral (34) was written as:

$${\Theta }_{\delta }^{\ast }(\zeta ,\tau )=\underset{\delta \to 0}{\lim }\left(-\frac{1}{2\pi i}{\displaystyle \underset{-\pi }{\overset{\pi }{\int }}{\overline{\Theta }}_{\delta }^{\ast }}(\zeta ,\delta e^{i\phi })e^{\delta e^{i\phi }\tau }i\delta e^{i\phi }d\phi \right),\tau \ge 0,$$

(40)

where

$${\overline{\Theta }}_{\delta }^{\ast }(\zeta ,\delta e^{i\phi })=\frac{{\Delta }_1(\zeta ,\delta e^{i\phi })}{\delta \sqrt{\delta }{e}^{1.5i\phi }\Delta (\delta e^{i\phi })},\zeta \ge 0,{\overline{\Theta }}_{\delta }^{\ast }(\zeta ,\delta e^{i\phi })=\frac{{\Delta }_2(\zeta ,\delta e^{i\phi })}{\delta \sqrt{\delta }{e}^{1.5i\phi }\Delta (\delta e^{i\phi })},\zeta \le 0,$$

(41)

$${\Delta }_1(\zeta ,\delta e^{i\phi })=e^{-0.5{\gamma }^{\ast }\zeta }{I}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{\delta }{e}^{0.5i\phi }{e}^{-0.5{\gamma }^{\ast }\zeta }\right),$$

(42)

$${\Delta }_2(\zeta ,\delta e^{i\phi })=e^{\sqrt{\frac{\delta }{k^{\ast }}}\hspace{0.17em}\zeta \hspace{0.17em}{e}^{0.5i\phi }}{I}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{\delta }{e}^{0.5i\phi }\right),$$

(43)

$${\Delta }^{\pm }(\zeta ,r)=K_{\epsilon }{J}_1\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}\right)\pm i{J}_0\left(\frac{2}{{\gamma }^{\ast }}\sqrt{r}\right).$$

(44)

Substituting the functions (41)–(44) into Equation (40), it was found:

$${\Theta }_{\delta }^{\ast }(\zeta ,\tau )=\underset{\delta \to 0}{\lim }\left(-\frac{1}{2\pi }{\displaystyle \underset{-\pi }{\overset{\pi }{\int }}\frac{{\Delta }_1(\zeta ,\delta e^{i\phi })}{\sqrt{\delta }{e}^{0.5i\phi }\Delta (\delta e^{i\phi })}}{e}^{\delta e^{i\phi }\tau }d\phi \right),\zeta \ge 0,\tau \ge 0,$$

(45)

$${\Theta }_{\delta }^{\ast }(\zeta ,\tau )=\underset{\delta \to 0}{\lim }\left(-\frac{1}{2\pi }{\displaystyle \underset{-\pi }{\overset{\pi }{\int }}\frac{{\Delta }_2(\zeta ,\delta e^{i\phi })}{\sqrt{\delta }{e}^{0.5i\phi }\Delta (\delta e^{i\phi })}}{e}^{\delta e^{i\phi }\tau }d\phi \right),\zeta \le 0,\tau \ge 0.$$

(46)

For small values of the argument [8]:

$$I_0(x)\cong 1,I_1(x)\cong 0.5x,$$

(47)

from Equations (45) and (46), the following was obtained:

$${\Theta }_{\delta }^{\ast }(\zeta ,\tau )=-\frac{1}{{\gamma }^{\ast }}{e}^{-0.5{\gamma }^{\ast }\zeta },\zeta \ge 0,{\Theta }_{\delta }^{\ast }(\zeta ,\tau )=-\frac{1}{{\gamma }^{\ast }},\zeta \le 0,\tau \ge 0.$$

(48)

Applying the function ${\Theta }_{\pm }^{\ast }(\zeta ,\tau )$ (35), (37)–(39), and ${\Theta }_{\delta }^{\ast }(\zeta ,\tau )$ (48) into the Equation (33) and introducing the notation: $\sqrt{r}=x$, $r=x^2$, the dimensional rise of temperature was found in the form:

$${\Theta }^{\ast }(\zeta ,\tau )=\frac{1}{{\gamma }^{\ast }}\left[e^{-0.5{\gamma }^{\ast }\zeta }-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}F(x)G_1(\zeta ,x)e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }dx}\right],\zeta \ge 0,\tau \ge 0,$$

(49)

$${\Theta }^{\ast }(\zeta ,\tau )=\frac{1}{{\gamma }^{\ast }}\left[1-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}F(x)G_2(\zeta ,x)e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }dx}\right],\zeta \le 0,\tau \ge 0,$$

(50)

where

$$F(x)=\frac{J_1(x)}{x^2\{{[J_0(x)]}^2+{[K_{\epsilon }{J}_1(x)]}^2\}},$$

(51)

$$G_1(\zeta ,x)=K_{\epsilon }{e}^{-0.5{\gamma }^{\ast }\zeta }{J}_1(x{e}^{-0.5{\gamma }^{\ast }\zeta }),$$

(52)

$$G_2(\zeta ,x)=K_{\epsilon }{J}_1(x)\cos \left(\frac{{\gamma }^{\ast }\zeta }{2\sqrt{k^{\ast }}}x\right)-J_0(x)\sin \left(\frac{{\gamma }^{\ast }\zeta }{2\sqrt{k^{\ast }}}x\right).$$

(53)

Substituting $\zeta =0$ into Equations (49)–(53) it was established that the temperature rise on the contact surface included in the boundary condition (24) has the form:

$${\Theta }^{\ast }(\tau )=\frac{1}{{\gamma }^{\ast }}\left[1-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}G(x)e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }dx}\right],\tau \ge 0,$$

(54)

where

$$G(x)=\frac{K_{\epsilon }{[J_1(x)]}^2}{x^2\{{[J_0(x)]}^2+{[K_{\epsilon }{J}_1(x)]}^2\}}.$$

(55)

On the basis of the Fourier’s law, the intensities of heat fluxes directed from the contact surface towards the inside of the friction pair elements were defined:

$$q_1(t)={-K_{1,1}\frac{\partial \Theta (z,t)}{\partial z}|}_{z=0^{+}},q_2(t)={K_2\frac{\partial \Theta (z,t)}{\partial z}|}_{z=0^{-}},t\ge 0,$$

(56)

The dimensionless form of dependencies (56) can be found as:

$$q_l^{\ast }(\tau )=\frac{q_l(t)}{q_0},l=1,2,$$

(57)

and taking account of (8) and (18), it was obtained:

$$q_1^{\ast }(\tau )={-\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \zeta }|}_{\zeta =0^{+}},q_2^{\ast }(\tau )={K^{\ast }\frac{\partial {\Theta }^{\ast }(\zeta ,\tau )}{\partial \zeta }|}_{\zeta =0^{-}},\tau \ge 0.$$

(58)

Substituting the dimensionless temperature rise (49)–(53) into Equation (58) and differentiating, it was found:

$$q_1^{\ast }(\tau )=1+\frac{2}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}Q(x)e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }dx},q_2^{\ast }(\tau )=-\frac{2}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}Q(x)e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }dx},\tau \ge 0,$$

(59)

where

$$Q(x)=\frac{K_{\epsilon }{J}_0(x)J_1(x)}{x\{{[J_0(x)]}^2+{[K_{\epsilon }{J}_1(x)]}^2\}}.$$

(60)

From Equations (59) and (60) it follows that $q_1^{\ast }(\tau )+q_2^{\ast }(\tau )=1$, which confirms the fulfillment of the boundary condition (16) for $q^{\ast }(\tau )=1,\hspace{0.17em}\tau \ge 0$.

4. Asymptotic Solutions

It should be noted that solutions (49)–(55) have the form of a quadrature; thus, using them, numerical integration should be performed each time on the range of bounded fields. However, in the case of small and large values of dimensionless time $\tau$ (Fourier number), the corresponding asymptotic solution will be obtained in the analytical form, not requiring numerical integration.

Small values of the Fourier number$0\le \tau <<1$ (large values of the parameter $p$ of the Laplace integral transform (30)). At large values of arguments, the modified Bessel functions behave as follows [8]:

$$I_0(x)\cong \frac{e^x}{\sqrt{2\pi x}}\left(1+\frac{1}{8x}+\frac{9}{128x^2}+\dots \right),I_1(x)\cong \frac{e^x}{\sqrt{2\pi x}}\left(1-\frac{3}{8x}-\frac{15}{128x^2}-\dots \right).$$

(61)

Limiting only to the first two components in the formula (61), the transforms of the dimensionless temperature rise (27)–(29) were written in the form:

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{-0.25{\gamma }^{\ast }\zeta -\alpha \sqrt{p}}}{(1+K_{\epsilon })p\sqrt{p}}\left(1-\frac{3{\gamma }^{\ast }{e}^{0.5{\gamma }^{\ast }\zeta }}{16\sqrt{p}}\right){\left[1+\frac{{\gamma }^{\ast }(1-3K_{\epsilon })}{16(1+K_{\epsilon })\sqrt{p}}\right]}^{-1},\zeta \ge 0,$$

(62)

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{\sqrt{\frac{p}{k^{\ast }}\hspace{0.17em}}\zeta }}{(1+K_{\epsilon })p\sqrt{p}}\left(1-\frac{3{\gamma }^{\ast }}{16\sqrt{p}}\right){\left[1+\frac{{\gamma }^{\ast }(1-3K_{\epsilon }}{16(1+K_{\epsilon })\sqrt{p}}\right]}^{-1},\zeta \le 0,$$

(63)

where

$$\alpha =\frac{2}{{\gamma }^{\ast }}(1-e^{-0.5{\gamma }^{\ast }\zeta }),\zeta \ge 0.$$

(64)

Taking into consideration that:

$$\left(1-\frac{3{\gamma }^{\ast }{e}^{0.5{\gamma }^{\ast }\zeta }}{16\sqrt{p}}\right){\left[1+\frac{{\gamma }^{\ast }(1-3K_{\epsilon }}{16(1+K_{\epsilon })\sqrt{p}}\right]}^{-1}\approx 1-\frac{{\gamma }^{\ast }}{16\sqrt{p}}\left(3e^{0.5{\gamma }^{\ast }\zeta }+\frac{1-3K_{\epsilon }}{1+K_{\epsilon }}\right),$$

(65)

$$\left(1-\frac{3{\gamma }^{\ast }}{16\sqrt{p}}\right){\left[1+\frac{{\gamma }^{\ast }(1-3K_{\epsilon }}{16(1+K_{\epsilon })\sqrt{p}}\right]}^{-1}\approx 1-\frac{{\gamma }^{\ast }}{4(1+K_{\epsilon })\sqrt{p}},$$

(66)

the transforms (62)–(64) were obtained in the form:

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{-0.25{\gamma }^{\ast }\zeta -\alpha \sqrt{p}}}{(1+K_{\epsilon })p\sqrt{p}}\left[1-\frac{{\gamma }^{\ast }}{16\sqrt{p}}\left(3e^{0.5{\gamma }^{\ast }\zeta }+\frac{1-3K_{\epsilon }}{1+K_{\epsilon }}\right)\right],\zeta \ge 0,$$

(67)

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{\sqrt{\frac{p}{k^{\ast }}\hspace{0.17em}}\zeta }}{(1+K_{\epsilon })p\sqrt{p}}\left(1-\frac{{\gamma }^{\ast }}{4(1+K_{\epsilon })\sqrt{p}}\right),\zeta \le 0.$$

(68)

Taking account of the relations [10]:

$$L^{-1}\left[\frac{e^{-a\sqrt{p}}}{p\sqrt{p^n}};\tau \right]={(4\tau )}^{\frac{n}{2}}{\mathrm{i}}^{\mathrm{n}}\mathrm{erfc}\left(\frac{a}{2\sqrt{\tau }}\right),n=1,2,a\ge 0,$$

(69)

from the transforms of solutions (67) and (68), the dimensionless temperature rises were found:

$$\begin{array}{c}{\Theta }^{\ast }(\zeta ,\tau )\cong \frac{2e^{-0.25{\gamma }^{\ast }\zeta }\sqrt{\tau }}{(1+K_{\epsilon })}\left[\mathrm{ierfc}\left(\frac{\alpha }{2\sqrt{\tau }}\right)-\frac{{\gamma }^{\ast }\sqrt{\tau }}{8}\left(3e^{0.5{\gamma }^{\ast }\zeta }+\frac{1-3K_{\epsilon }}{1+K_{\epsilon }}\right){\mathrm{i}}^2\mathrm{erfc}\left(\frac{\alpha }{2\sqrt{\tau }}\right)\right],\\ \zeta \ge 0,\end{array}$$

(70)

$$\begin{array}{c}{\Theta }^{\ast }(\zeta ,\tau )\cong \frac{2\sqrt{\tau }}{(1+K_{\epsilon })}\left[\mathrm{ierfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}\right)-\frac{{\gamma }^{\ast }\sqrt{\tau }}{2(1+K_{\epsilon })}{\mathrm{i}}^2\mathrm{erfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}\right)\right],\zeta \le 0,\\ 0\le \tau <<1,\end{array}$$

(71)

where

$$\begin{array}{c}{\mathrm{i}}^2\mathrm{erfc}(x)=0.25[\mathrm{erfc}(x)-2x\hspace{0.17em}\mathrm{ierfc}(x)],\mathrm{ierfc}(x)={\pi }^{-0.5}{e}^{-x^2}-x\hspace{0.17em}\mathrm{erfc}(x),\\ \mathrm{erfc}(x)=1-\mathrm{erf}(x),\end{array}$$

(72)

$\mathrm{erf}(x)$—Gauss error function [8]. On the contact surface $\zeta =0$ from Equations, (70) and (71) it was obtained:

$${\Theta }^{\ast }(\tau )\cong \frac{2\sqrt{\tau }}{(1+K_{\epsilon })}\left[\frac{1}{\sqrt{\pi }}-\frac{{\gamma }^{\ast }\sqrt{\tau }}{8(1+K_{\epsilon })}\right],0\le \tau <<1.$$

(73)

Approaching in Equations (70)–(73) the limit ${\gamma }^{\ast }\to 0$ ($\alpha \to \zeta$), the solution for homogeneous materials was obtained [11]:

$${\Theta }^{\ast }(\zeta ,\tau )\cong \frac{2\sqrt{\tau }}{(1+K_{\epsilon })}\mathrm{ierfc}\left(\frac{\zeta }{2\sqrt{\tau }}\right),\zeta \ge 0,$$

(74)

$${\Theta }^{\ast }(\zeta ,\tau )\cong \frac{2\sqrt{\tau }}{(1+K_{\epsilon })}\mathrm{ierfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}\right),\zeta \le 0,$$

(75)

$${\Theta }^{\ast }(\tau )\cong \frac{2}{(1+K_{\epsilon })}\sqrt{\frac{\tau }{\pi }\hspace{0.17em}},0\le \tau <<1.$$

(76)

Large values of Fourier number$\tau >>1$(small values of the parameter $p$ of the Laplace integral transform (30)). Distributions of the modified Bessel functions for small values of argument in the power series have the form [8]:

$$I_0(x)\cong 1+\frac{x^2}{4}+\dots ,I_1(x)\cong \frac{x}{2}\left(1+\frac{x^2}{8}+\dots \right).$$

(77)

Taking into account the relations (77), the Laplace transforms of dimensionless temperature rise (27)–(29) were written as:

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{-{\gamma }^{\ast }\zeta }}{{\gamma }^{\ast }}\left[\frac{\beta }{p(\beta +\sqrt{p})}+\frac{e^{-{\gamma }^{\ast }\zeta }}{2K_{\epsilon }{\gamma }^{\ast }(\beta +\sqrt{p}}\right],\zeta \ge 0,$$

(78)

$${\overline{\Theta }}^{\ast }(\zeta ,p)\cong \frac{e^{\sqrt{\frac{p}{k^{\ast }}\hspace{0.17em}}\zeta }}{{\gamma }^{\ast }}\left[\frac{\beta }{p(\beta +\sqrt{p})}+\frac{1}{2K_{\epsilon }{\gamma }^{\ast }(\beta +\sqrt{p}}\right],\zeta \le 0,$$

(79)

where

$$\beta =\frac{{\gamma }^{\ast }}{K_{\epsilon }}.$$

(80)

Using the dependencies [10]:

$$L^{-1}\left[\frac{e^{-a\sqrt{p}}}{(\beta +\sqrt{p})};\tau \right]=\frac{e^{-\frac{a^2}{4\tau }}}{\sqrt{\pi \tau }}-\beta \hspace{0.17em}{e}^{a\beta +{\beta }^2\tau }\mathrm{erfc}\left(\frac{a}{2\sqrt{\tau }}+\beta \sqrt{\tau }\right),$$

(81)

$$L^{-1}\left[\frac{\beta \hspace{0.17em}{e}^{-a\sqrt{p}}}{p(\beta +\sqrt{p})};\tau \right]=\mathrm{erfc}\left(\frac{a}{2\sqrt{\tau }}\right)-e^{a\beta +{\beta }^2\tau }\mathrm{erfc}\left(\frac{a}{2\sqrt{\tau }}+\beta \sqrt{\tau }\right),a\ge 0,$$

(82)

from the transform solutions (78) and (79), the dimensionless temperature rises were obtained in the form:

$${\Theta }^{\ast }(\zeta ,\tau )\cong \frac{e^{-{\gamma }^{\ast }\zeta }}{{\gamma }^{\ast }}\left\{1-e^{{\beta }^2\tau }\mathrm{erfc}(\beta \sqrt{\tau })+\frac{e^{-{\gamma }^{\ast }\zeta }}{2K_{\epsilon }{\gamma }^{\ast }}\left[\frac{1}{\sqrt{\pi \tau }}-\beta \hspace{0.17em}{e}^{{\beta }^2\tau }\mathrm{erfc}(\beta \sqrt{\tau })\right]\right\},\zeta \ge 0,\tau >>1,$$

(83)

$$\begin{array}{l}{\Theta }^{\ast }(\zeta ,\tau )\cong \frac{1}{{\gamma }^{\ast }}\{\mathrm{erfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}\right)-e^{\frac{\beta \left|\zeta \right|}{\sqrt{k^{\ast }}}+{\beta }^2\tau }\mathrm{erfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}+\beta \sqrt{\tau }\right)+\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\frac{1}{2K_{\epsilon }{\gamma }^{\ast }}\left[\frac{e^{-\frac{{\zeta }^2}{4k^{\ast }\tau }}}{\sqrt{\pi \tau }}-\beta e^{\frac{\beta \left|\zeta \right|}{\sqrt{k^{\ast }}}+{\beta }^2\tau }\mathrm{erfc}\left(\frac{\left|\zeta \right|}{2\sqrt{k^{\ast }\tau }}+\beta \sqrt{\tau }\right)\right]\},\hspace{0.17em}\hspace{0.17em}\zeta \le 0,\hspace{0.17em}\hspace{0.17em}\tau >>1.\end{array}$$

(84)

Substituting $\zeta =0$ into Equations (83) and (84), it was found:

$${\Theta }^{\ast }(\tau )\cong \frac{1}{{\gamma }^{\ast }}\left\{1-e^{{\beta }^2\tau }\mathrm{erfc}(\beta \sqrt{\tau })+\frac{1}{2K_{\epsilon }{\gamma }^{\ast }}\left[\frac{1}{\sqrt{\pi \tau }}-\beta \hspace{0.17em}{e}^{{\beta }^2\tau }\mathrm{erfc}(\beta \sqrt{\tau })\right]\right\},\tau >>1.$$

(85)

5. Temperature Field during Braking with Constant Deceleration

Based on Duhamel’s theorem [12], the dimensionless temperature rise during braking with constant deceleration was sought in the form:

$${\widehat{\Theta }}^{\ast }(\zeta ,\tau )=\frac{\partial }{\partial \tau }{\displaystyle \underset{0}{\overset{\tau }{\int }}{q}^{\ast }(\tau -s){\Theta }^{\ast }(\zeta ,s)ds},\left|\zeta \right|<\infty ,0\le \tau \le {\tau }_s,$$

(86)

where the temporal profiles of the specific friction power $q^{\ast }(\tau )$ and function ${\Theta }^{\ast }(\zeta ,\tau )$ were determined from Equations (19), (49), and (50), respectively. Performing the integration first, and then differentiating, from the Equation (86) we obtained:

$${\widehat{\Theta }}^{\ast }(\zeta ,\tau )=\frac{1}{{\gamma }^{\ast }}\left[e^{-0.5{\gamma }^{\ast }\zeta }{q}^{\ast }(\tau )-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}F(x)G_1(\zeta ,x)P(\tau ,x)dx}\right],\zeta \ge 0,0\le \tau \le {\tau }_s,$$

(87)

$${\widehat{\Theta }}^{\ast }(\zeta ,\tau )=\frac{1}{{\gamma }^{\ast }}\left[1-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}F(x)G_2(\zeta ,x)P(\tau ,x)dx}\right],\zeta \le 0,0\le \tau \le {\tau }_s,$$

(88)

where

$$P(\tau ,x)=e^{-{(0.5{\gamma }^{\ast }x)}^2\tau }-\frac{(1-e^{-{(0.5{\gamma }^{\ast }x)}^2\tau })}{{(0.5{\gamma }^{\ast }x)}^2{\tau }_s},$$

(89)

and functions $F(x)$, $G_1(\zeta ,x)$, and $G_2(\zeta ,x)$ can be found from the formulas (51)–(53).

The temperature change on the friction surface was found, substituting $\zeta =0$ into the Equations (87) and (88), in the form:

$${\widehat{\Theta }}^{\ast }(\tau )=\frac{1}{{\gamma }^{\ast }}\left[1-\frac{4}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}G(x)P(\tau ,x)dx}\right],\tau \ge 0,$$

(90)

where functions $G(x)$ and $P(\tau ,x)$ were determined from relations (55) and (89), respectively.

Knowing the dimensionless temperature rise (87), (88), from formulas (58) the dimensionless intensities of frictional heat fluxes were found:

$${\widehat{q}}_1^{\ast }(\tau )=q^{\ast }(\tau )+\frac{2}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}Q(x)P(\tau ,x)dx},{\widehat{q}}_2^{\ast }(\tau )=-\frac{2}{\pi }{\displaystyle \underset{0}{\overset{\infty }{\int }}Q(x)P(\tau ,x)dx},0\le \tau \le {\tau }_s,$$

(91)

where functions $Q(x)$ and $P(\tau ,x)$ have the forms (60) and (89), respectively. From Equation (91) it follows that ${\widehat{q}}_1^{\ast }(\tau )+{\widehat{q}}_2^{\ast }(\tau )=q^{\ast }(\tau )$, which confirms the fulfillment of the boundary condition (16) with the dimensionless specific friction power $q^{\ast }(\tau )$ in the form (19).

6. Numerical Analysis

Calculations were performed for a friction pair, where the first element (pad) is made of two-component FGM: zircon dioxide ZrO₂ (friction surface) and titanium alloy $\mathrm{Ti}-6\mathrm{Al}-4\mathrm{V}$(core). While the second material (brake disc) is homogeneous: cast iron ChNMKh. The properties of the materials are included in

Table 1. Material properties at the initial temperature $T_0$ [3,13].

Material	$\begin{array}{c}\mathbf{Thermal}\mathbf{Conductivity}\\ {\mathbf{Wm}}^{-1}{\mathbf{K}}^{-1}\end{array}$	$\begin{array}{c}\mathbf{Specific}\mathbf{Heat}\mathbf{Capacity}\\ \mathbf{J}\hspace{0.17em}{\mathbf{kg}}^{-1}{\mathbf{K}}^{-1}\end{array}$	$\mathbf{Density}\hspace{0.17em}\mathbf{kg}\hspace{0.17em}{\mathbf{m}}^{-3}$
ZrO2	1.94	452.83	6102.16
Ti-6Al-4V	6.87	538.08	4431.79
ChNMKh	52.17	444.6	7100

.

The values of the remaining input parameters used to perform the calculations are listed in

Table 2. Input parameters [14].

$$\begin{gathered} \mathbf{Friction} \\ \mathbf{Coefficient}{\mathit{f}}_0 \end{gathered}$$	$$\mathbf{Nominal}\mathbf{Pressure}{\mathit{p}}_0^{},\hspace{0.17em}\mathbf{MPa}$$	$$\mathbf{Initial}\mathbf{Sliding}\mathbf{Speed}{\mathit{V}}_0^{},\hspace{0.17em}{\mathbf{ms}}^{-1}$$	$$\mathbf{Initial}\mathbf{Kinetic}\mathbf{Energy}{\mathit{W}}_0^{},\hspace{0.17em}\mathbf{kJ}$$	$$\mathbf{Nominal}\mathbf{Contact}\mathbf{Area}{\mathit{A}}_{\mathit{a}},\hspace{0.17em}{\mathbf{m}}^2$$	$$\begin{gathered} \mathbf{Initial} \\ \mathbf{Temperature}{\mathit{T}}_0^{},°\mathbf{C} \end{gathered}$$
0.27	0.602	23.8	103.54	0.00221	20

.

Then, from formulas (1) and (19), the nominal value of specific friction power $q_0=3.87\mathrm{MW}{\mathrm{m}}^{-2}$, braking time $t_s=12.1\mathrm{s}$, and gradient parameter ${\gamma }^{\ast }=1.26$ were determined. Based on Equation (9), for an equal volumetric fraction of FGM component $(v=0.5)$, the effective values of specific heat capacity and density of the pad material were obtained, $c_1=495.45\mathrm{J}{\mathrm{kg}}^{-1}{\mathrm{K}}^{-1}$, ${\rho }_1=5266.97\mathrm{kg}{\mathrm{m}}^{-3}$, respectively. Thereafter, the following parameters were calculated sequentially: thermal diffusivity $k_1=0.743\cdot {10}^{-6}{\mathrm{m}}^2{\mathrm{s}}^{-1}$ and $k_2=1.65\cdot {10}^{-5}{\mathrm{m}}^2{\mathrm{s}}^{-1}$, the effective depth of heat penetration of the pad $a=5.2\mathrm{mm}$, the dimensionless braking time ${\tau }_s=0.33$, and the temperature scaling factor ${\Theta }_0=\mathrm{10,373}{}^{\circ }\mathrm{C}$, based on Equations (10)–(12).

The integrals in the obtained solutions were calculated numerically using the QAGI procedure of the QUADPACK package [15]. Changes of the dimensionless temperature rise and intensities of heat fluxes during sliding with a constant velocity are presented in Figure 3, Figure 4, Figure 5 and Figure 6.

Temporal profiles of the dimensionless temperature rise ${\Theta }^{\ast }(\zeta ,\tau )$ (49)–(53) at a few distances from the friction surface are shown in Figure 3. The temperature of both elements increased monotonically over time. The highest temperature was reached on the friction surface, and decreasing moving away from it. For a given distance from this surface, the temperature of the homogeneous cast iron element was always higher than the temperature of the functionally graded element. Having a much greater thermal conductivity, the cast iron was heated to a much deeper extent than the FGM (Figure 4).

Temporal profiles of dimensionless heat flux intensities $q_l^{\ast }(\tau )$, $l=1,2$ (59), (60) are demonstrated in Figure 5. It was found that the main element that absorbs frictional heat is the cast-iron disc, especially at the initial stage of the heating process. The amount of heat directed from the friction surface towards the inside of the pad increases with time, and towards the inside of the disc it decreases. A comparison of dimensionless temperature values ${\Theta }^{\ast }(\tau )$ of the friction surface, found by means of the exact (54), (55) and asymptotic solutions (74), (75) are shown in Figure 6. In the considered range of Fourier number $0\le \tau \le {\tau }_s$, the respective temperature values were almost the same.

Relevant results, obtained in the case of a linearly decreasing velocity (so-called braking with a constant deceleration), are presented in Figure 7, Figure 8, Figure 9 and Figure 10. The temporal profile of the dimensionless temperature rise ${\widehat{\Theta }}^{\ast }(\zeta ,\tau )$ (87)–(89) during the braking process was different than during uniform sliding (Figure 7). The dimensionless time to reach the maximum temperature on the friction surface was ${\tau }_{\max }\approx 0.5{\tau }_s$ and became higher when increasing the distance from it. After reaching the maximum value, the temperature dropped. More vividly, such a concentration of high temperature near the friction surface is shown in the distribution of isotherms, as illustrated in Figure 8. Apparently, as in the case of uniform sliding, the greater part of the frictionally-generated heat is absorbed by the cast iron disc ($\approx 85\%$) (Figure 9). The intensities of heat fluxes ${\widehat{q}}_l^{\ast }(\tau )$, $l=1,2$ (91) are almost unchanged during the entire braking process.

The results presented in Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8 and Figure 9 were obtained with the same ($v=0.5$) volumetric components fractions of ZrO₂ and Ti-6Al-4V, determining the effective specific heat capacity and density using formula (9). On the other hand, the change of the maximum temperature $T_{\max }\equiv T(0,t_{\max })$ with the increase of the parameter $v$ is presented in Figure 10. The highest value $T_{\max }=1117{}^{\circ }\mathrm{C}$ was achieved in the case of the pad made of pure zirconium dioxide, and the lowest $T_{\max }=1052{}^{\circ }\mathrm{C}$, when it was made of the titanium alloy.

7. Conclusions

A mathematical model was proposed to determine the transient temperature field in a friction pair, in which one element is made of a functionally graded material and the other is made of a homogeneous material. It was assumed that the thermal conductivity of a FGM increases exponentially with the distance from the contact surface. An exact solution of the appropriate boundary–value problem of heat conduction was formulated and then solved, with consideration of frictional heat generation. Two cases of the friction power temporal profiles were analyzed in detail: constant (uniform sliding), and linearly decreasing in time (braking with constant deceleration). A numerical analysis was performed for a two-component FGM (ZrO₂ + Ti-6Al-4V) sliding on the cast-iron disc. It was found that the greater part of heat generated due to friction was absorbed by the cast iron (about 85%), which resulted in a greater depth of effective heat penetration in this element, due to the high thermal conductivity of cast iron. At a fixed distance from the friction surface, the temperature of the cast iron element is higher than that of the FGM element, in both the considered cases: uniform sliding, and during braking. Thus, in order to protect systems against such undesirable phenomena as overheating and thermal cracking etc., the use of FGM on the friction elements may be justified. It is also worth emphasizing that in the analyzed range of the Fourier number change $0\le \tau \le 0.33$, the appropriate asymptotic solution can be effectively used, giving a high accuracy of calculations, without the inconveniences related to numerical integration in an exact solution.

It should be noted that the shape of the friction pair elements, as well as their positional relationship, can be considered in some spatial problems of friction solved by numerical methods, in particular the finite element method (FEM). The temperature evolution obtained by them oscillates, as a result of the heating area moving on the surface of the brake disc. The model proposed in this paper is one-dimensional, based on a physically-justified assumption that heat, generated by friction of two elements, propagates in the direction perpendicular to the contact surface. This allows determining the mean temperature (from the above-mentioned oscillations) on the friction surfaces of both elements.

According to the current state of knowledge [16,17], the temperature of the friction surface is the sum of the volume temperature (average temperature in volume), the mean temperature, and the flash temperature. The flash temperature is the component that takes into consideration the texture of the friction surfaces. The flash temperature calculation models need appropriate experimental data as input parameters. In the case of homogeneous materials, such data can be found in the article in ref. [13]. However, we have not found such data for the considered friction pair. The development of models for determining the flash temperature of such couples is a potential direction for our research. In the future, we intend to expand the proposed mathematical model with the possibility of testing the temperature of friction systems of this type (functionally graded and homogeneous materials) made of thermally sensitive materials and the temperature-dependent friction coefficient.