*Article* **The Heat Partition Ratio during Braking in a Functionally Graded Friction Couple**

**Aleksander Yevtushenko , Katarzyna Topczewska \* and Przemysław Zamojski**

Faculty of Mechanical Engineering, Bialystok University of Technology (BUT), 45C Wiejska Street, 15-351 Bialystok, Poland; a.yevtushenko@pb.edu.pl (A.Y.); p.zamojski@pb.edu.pl (P.Z.)

**\*** Correspondence: k.topczewska@pb.edu.pl

**Abstract:** The theoretical scheme for determining the heat partition ratio (HPR) in a friction couple made of functionally graded materials (FGMs) was proposed. As a result, the formula for the calculation of the HPR was found, which depends on the thermal properties and the parameters of the material's gradient. In specific cases of these parameters, the known formulas for estimating the HPR for homogeneous materials were obtained. Calculations were carried out for the friction couple consisting of the following two-component FGMs: Al2O3–Cu (first body) and ZrO2–Ti–6Al–4V (second body), under the conditions corresponding to a single braking with a constant deceleration. It was established that the vast majority (almost 90%) of heat that was generated by friction was absorbed by the first body in the selected couple. The possibilities of using the obtained results were discussed herein.

**Keywords:** heat partition ratio; functionally graded materials; frictional heating; temperature; braking

### **1. Introduction**

Functionally graded materials (FGMs) are materials in which, along at least one specific direction, a continuous change in the functional or structural properties has been obtained in a selected technological process. The possibilities of using graded materials in technology seem to be almost unlimited. Proposals of possible FGM applications appear for many industries, such as the following: optical, energy [1,2], aviation [3,4], medical [5], etc. Moreover, they have been used in friction elements, including braking systems [6,7].

The significant influence of temperature on the friction and wear characteristics of friction materials is well known and has been studied by many researchers who are involved in determining the temperature mode of disc brakes [8–10]. The basis for establishing the temperature are the solutions (analytical, numerical, or analytical–numerical) to the thermal problems of friction, i.e., the boundary-value problems of heat conduction with two specific boundary conditions on the contact surface of the pad and the disc [11,12]. One of them, the energetic solution, determines the equality of the specific friction power *q* and the sum of the heat flux intensities *ql*,*l* = 1, 2 that are directed along the normal surface to the contact surface towards the insides of the friction elements *q*<sup>1</sup> + *q*<sup>2</sup> = *q*. Whereas the second solution concerns the type of thermal contact—with (imperfect) or without (perfect) consideration of the thermal resistance on the friction surfaces. In the latter case, the surface temperatures of the friction elements are the same. The coupling of the temperature fields of both of the friction elements through the above-mentioned boundary conditions means that obtaining the analytical [13–16] or the numerical [17–19] solutions of the thermal problems of friction requires the performance of complex mathematical transformations or long-term calculations.

In order to reduce these difficulties, there is also another approach to estimate the temperature of the braking systems, on the basis of solutions to the thermal problems of friction. It is based on a virtual separation of the friction pair elements and the subsequent

**Citation:** Yevtushenko, A.; Topczewska, K.; Zamojski, P. The Heat Partition Ratio during Braking in a Functionally Graded Friction Couple. *Materials* **2022**, *15*, 4623. https://doi.org/10.3390/ ma15134623

Academic Editor: Pawel Pawlus

Received: 1 June 2022 Accepted: 28 June 2022 Published: 30 June 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

replacement of the frictional interaction on their working surfaces with the heat fluxes of the following intensities: *ql* = *αlq*, *l* = 1, 2, *α*<sup>1</sup> = *α*, and *α*<sup>2</sup> = 1 − *α*, respectively, where *α* is the heat partition ratio [20,21]. Reviews of the experimental methods, as well as the theoretical methods, for determining the HPR in braking systems have been developed [22– 25]. The theoretical methods primarily rely on the fact that the temperature field, which is first found by means of analytical solutions, contains *a priori* unknown coefficient *α*, which is then determined from the condition of the equality of the maximum or the mean temperature on the friction surfaces of the pad and the disc. As a result, the formulas for the calculation of *α* contain the thermo-physical properties of the friction materials and, in some cases, the operating parameters of the process. Substituting thus determines the value of *α* back into the analytical or the numerical solutions, and the temperature mode of the braking system is estimated.

It should be noted that the formulas that have been obtained so far for determining the HPR, based on the analytical solutions to the thermal problems of friction, concern homogeneous materials [26–30]. The proportion of the heat distribution by means of the HPR between polymer and steel elements was established experimentally [31]. The aim of this study is to obtain the theoretical formulas for estimating the HPR during the braking of friction pair elements that are made of functionally graded materials (FGMs). One of the proposed applications of the results could be the determination of the bulk temperature of friction elements that are made of thermally sensitive FGM, during a single braking process [32] or a repetitive short-term braking mode [33].

#### **2. Heating of the FGM Semi-Space by the Heat Flux with Constant Intensity**

Consider a temperature of *T* of a semi-limited body *z* ≥ 0, which is uniform and equal to *T*<sup>0</sup> at the initial point in time of *t* = 0. The body is made of a functionally graded material (FGM) with the thermal conductivity of *K*, increasing exponentially in the positive direction of the axis *z* as follows:

$$K(z) = K\_0 e^{\gamma z}, z \ge 0,\\ K\_0 \equiv K(0), \tag{1}$$

where *γ* > 0 is the gradient of the material. Next, the semi-space is heated on the outer surface as *z* = 0 by the heat flux with a constant intensity of *q* = *q*<sup>0</sup> in time as *t* > 0 (Figure 1).

**Figure 1.** Scheme of the heating of the FGM semi-space.

The transient temperature field of the semi-space was searched in following the form:

$$T(z,t) = T\_0 + \Theta(z,t), \; z \ge 0, \; t \ge 0,\tag{2}$$

where the temperature rise of Θ(*z*, *t*) was determined from the solution to the following boundary-value heat conduction problem:

$$\frac{\partial}{\partial z}\left[K(z)\frac{\partial\Theta(z,t)}{\partial z}\right] = \rho c \frac{\partial\Theta(z,t)}{\partial t}, \; z > 0, \; t > 0,\tag{3}$$

$$K\_0 \frac{\partial \Theta(z, t)}{\partial z} \Big|\_{z=0^+} = -q\_0, \ t > 0,\tag{4}$$

Θ(*z*, *t*) → 0 , *z* → ∞, *t* > 0, (5)

$$
\Theta(z,0) = 0, \; z \ge 0,\tag{6}
$$

and *ρ*, *c* are the density and the specific heat capacity of the material, respectively.

The following dimensionless variables and parameters were introduced:

$$\zeta = \frac{z}{a}, \ \gamma = \frac{\gamma^\*}{a}, \ \tau = \frac{kt}{a^2}, \ k = \frac{K\_0}{c\rho}, \ q^\* = \frac{q}{q\_0}, \ \Theta^\* = \frac{\Theta}{\Theta\_0}, \ \Theta\_0 = \frac{q\_0 a}{K\_0} \tag{7}$$

where *a* is the effective heating depth, i.e., the distance from the heated surface, at which the temperature is equal to 5% of the maximum temperature on the heated surface [2]. Taking into account the designations (7), the problem (3)–(6) was written in following the form:

$$\frac{\partial^2 \Theta^\*(\zeta, \tau)}{\partial \zeta^2} + \gamma^\* \frac{\partial \Theta^\*(\zeta, \tau)}{\partial \zeta} - e^{-\gamma^\* \zeta} \frac{\partial \Theta^\*(\zeta, \tau)}{\partial \tau} = 0, \zeta > 0, \ \tau > 0,\tag{8}$$

$$\left. \frac{\partial \Theta^\*(\zeta\_{\prime} \tau)}{\partial \zeta} \right|\_{\zeta=0} = -1, \ \tau > 0,\tag{9}$$

$$
\Theta^\*(\zeta, \tau) \to 0 \text{ } \zeta \to \infty \text{ } \tau > 0,\tag{10}
$$

$$
\Theta^\*(\zeta, 0) = 0, \,\,\zeta \ge 0. \tag{11}
$$

Applying the following Laplace integral transform [34]:

$$\overline{\Theta}^\*(\zeta, p) \equiv \mathcal{L}[\Theta^\*(\zeta, \tau); p] = \int\_0^\infty \Theta^\*(\zeta, \tau) e^{-p\tau} d\tau, \text{Re}p \ge 0,\tag{12}$$

to the boundary-value problem (8)–(11), the following were obtained:

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

$$\left. \frac{d\overline{\Theta}^\*(\zeta\_\prime p)}{d\zeta} \right|\_{\zeta=0^+} = -\frac{1}{p'} \tag{14}$$

$$
\overline{\Theta}^\*(\zeta, p) \to 0 \; , \; \zeta \to \infty. \tag{15}
$$

The solution to the boundary problem (13)–(15) takes the following form:

$$\overline{\Theta}^\*(\zeta, p) = e^{-\frac{1}{2}\gamma\_1^\*\zeta} \frac{\varrho(\zeta, p)}{\Phi(p)}, \zeta \ge 0,\tag{16}$$

where

$$\varrho(\zeta, p) = \left. I\_1 \left( \frac{2}{\gamma^\*} \sqrt{p} \, e^{-\frac{1}{2} \gamma^\* \zeta} \right) \right| \tag{17}$$

$$\Phi(p) = p\sqrt{p}\,\mathrm{l}\_0\left(\frac{2}{\gamma^\*}\sqrt{p}\right),\tag{18}$$

*Ik*(*x*), *k* = 0, 1 are the modified Bessel functions of the first kind [35]. Differentiating the function (18) with consideration of the relation of *I* <sup>0</sup>(*x*) = *I*1(*x*), the following was found:

$$\Phi'(p) = \frac{3}{2}\sqrt{p}\operatorname{I}\_0\left(\frac{2}{\gamma^\*}\sqrt{p}\right) + \frac{p}{\gamma^\*}\operatorname{I}\_1\left(\frac{2}{\gamma^\*}\sqrt{p}\right). \tag{19}$$

The transition from the space of the transforms to the originals in the solution (16)–(19) was carried out on the basis of the Vashchenko-Zakharchenko theorem [36,37] as follows:

$$\Theta^\*(\zeta, \tau) = e^{-\frac{1}{2}\gamma^\*\zeta} \left[ \lim\_{p \to 0} \frac{\varrho(\zeta, p)p}{\Phi(p)} + \sum\_{n=1}^{\infty} \frac{\varrho(\zeta, p\_n)}{\Phi'(p\_n)} e^{p\_n \tau} \right], \zeta \ge 0, \ \tau \ge 0,\tag{20}$$

where

$$I\_0\left(\frac{2}{\gamma\_1^\*}\sqrt{p\_n}\right) \equiv 0, \ n = 1, 2, \dots \tag{21}$$

By using the expressions of the modified Bessel functions [35] as follows:

$$I\_0(\mathbf{x}) \cong 1 + \frac{\mathbf{x}^2}{4} + \frac{\mathbf{x}^4}{64} + \dots, \ I\_1(\mathbf{x}) \cong \frac{\mathbf{x}}{2} + \frac{\mathbf{x}^3}{16} + \dots \tag{22}$$

and by limiting them to the first components, the following representations of functions (17) and (18) for the small values of the parameter *p* were obtained:

$$\varphi(\zeta, p) \cong \frac{1}{\gamma^\*} e^{-\frac{1}{2}\gamma^\* \overline{\zeta}} \sqrt{p} \, \, \Phi(p) \cong p\sqrt{p}.\tag{23}$$

Taking into account the Formulas (23) and relationships [35] as follows:

$$I\_0(\mathbf{x}) = I\_0(i\mathbf{x}),\ I\_1(\mathbf{x}) = -iI\_1(i\mathbf{x}),\ i \equiv \sqrt{-1},\tag{24}$$

where *Jk*(*x*), *k* = 0, 1 are the Bessel functions of the first kind, solutions (20) and (21) were written in the following form:

$$\Theta^\*(\zeta,\tau) = \frac{1}{\gamma^\*} e^{-\frac{1}{2}\gamma^\*\zeta} \left[ e^{-\frac{1}{2}\gamma^\*\zeta} - 4 \sum\_{n=1}^\infty \frac{f\_1(\mu\_n e^{-\frac{1}{2}\gamma^\*\zeta})}{\mu\_n^2 f\_1(\mu\_n)} e^{-\lambda\_n \tau} \right], \zeta \ge 0, \ \tau \ge 0,\tag{25}$$

where

$$
\lambda\_n = \left(\frac{1}{2}\gamma^\*\mu\_n\right)^2,\tag{26}
$$

$$J\_0(\mu\_n) \equiv 0, \ n = 1, 2, \dots \tag{27}$$

By adopting *ζ* = 0 in the solution (25)–(27), the dimensionless temperature rise on the heated surface was found in the following form:

$$\Theta^\*(\tau) \equiv \Theta^\*(0, \tau) = \frac{1}{\gamma^\*} \left( 1 - 4 \sum\_{n=1}^{\infty} \frac{e^{-\lambda\_n \tau}}{\mu\_n^2} \right), \tau \ge 0. \tag{28}$$

Verification of the developed model was carried out by checking the boundaries (9), (10), and the initial (11) conditions. For this, by differentiating the solution (25) with respect to the spatial variable *ζ*, and taking into account the relationship [35] as follows:

$$J\_{\,1}^{\prime}(\mathbf{x}) = J\_0(\mathbf{x}) - \mathbf{x}^{-1} J\_1(\mathbf{x}) \tag{29}$$

the following was found:

$$\frac{\partial \Theta^\*(\zeta\_{\prime}\tau)}{\partial \zeta} = -e^{-\gamma^\*\zeta} \left[ 1 - 2 \sum\_{n=1}^{\infty} \frac{J\_0(\mu\_n e^{-\frac{1}{2}\gamma^\*\zeta})}{\mu\_n J\_1(\mu\_n)} e^{-\lambda\_n \tau} \right], \zeta \ge 0, \ \tau \ge 0. \tag{30}$$

Approaching Formula (30) to the limit of *ζ* → 0, with consideration of the fact that *μ<sup>n</sup>* are the roots of Equation (27), the following was achieved:

$$\left.\frac{\partial\Theta^\*(\zeta\_{\prime}\tau)}{\partial\zeta}\right|\_{\zeta=0^+} = -1,\tag{31}$$

which confirms that the boundary condition (9) was satisfied. However, by approaching to the limit *ζ* → ∞ in the solution (25), we have confirmed the fulfillment of the condition of the temperature fade (10). In addition, it was numerically established that the initial condition (11) was met. In particular, according to the Formula (25), on the heated surface *ζ* = 0 at the initial point in time, the dimensionless temperature rise is zero if the following occurs: <sup>∞</sup>

$$\sum\_{n=1}^{\infty} \frac{1}{\mu\_n^{\gamma}} = \frac{1}{4}.\tag{32}$$

Based on the calculations, it was found that the sum of the first 10<sup>2</sup> components of the series on the left-hand side of Equation (32) is equal to 0.248985.

In addition to the exact (25), the appropriate asymptotic solutions to the problem at small and large values of the Fourier number *τ* (dimensionless time) were also found.

*Small values τ* (*large values of the parameter p*). Taking into account that, in Formulas (17) and (18), the first components of the asymptotic of modified Bessel functions at large values of the argument [35] were as follows:

$$I\_0(\mathbf{x}) \cong \frac{e^{\mathbf{x}}}{\sqrt{2\pi\mathbf{x}}} \left( 1 + \frac{1}{8\mathbf{x}} + \dots \right), \ I\_1(\mathbf{x}) \cong \frac{e^{\mathbf{x}}}{\sqrt{2\pi\mathbf{x}}} \left( 1 + \frac{3}{8\mathbf{x}} + \dots \right), \tag{33}$$

the transformed solution (16) was written as follows:

$$\overline{\Theta}^\*(\zeta, p) = e^{-\frac{1}{4}\gamma^\*\zeta} \frac{e^{-b\sqrt{p}}}{p\sqrt{p}}, \; b = \frac{2}{\gamma^\*} \left(1 - e^{-\frac{1}{2}\gamma^\*\zeta} \right), \; \zeta \ge 0. \tag{34}$$

Using the relation [38] as follows:

$$L^{-1}\left[\frac{e^{-b\sqrt{p}}}{p\sqrt{p}};\,\tau\right] = 2\sqrt{\tau}\,\text{ierfc}\left(\frac{b}{2\sqrt{\tau}}\right),\tag{35}$$

from the Formula (34), the following form of dimensionless temperature rise at the initial moments of heating was obtained:

$$\Theta^\*(\zeta, \tau) \cong 2e^{-\frac{1}{4}\gamma^\*\zeta} \sqrt{\tau} \operatorname{ierfc}\left(\frac{b}{2\sqrt{\tau}}\right), \zeta \ge 0, \ 0 \le \tau << 1,\tag{36}$$

where ierfc(*x*) = *π*<sup>−</sup> <sup>1</sup> <sup>2</sup> *e*−*x*<sup>2</sup> − *x*erfc(*x*), erfc(*x*) = 1 − erf(*x*), erf(*x*) is the Gauss error function [35].

At *ζ* = 0 from the solution (36), a known result for the evolution of the temperature on the heated surface of a homogeneous semi-space was obtained [18] as follows:

$$\Theta^\*(\tau) \cong 2\sqrt{\frac{\tau}{\pi'}}, \ 0 \le \tau << 1,\tag{37}$$

*Large values τ (small values of the parameter p).* By including the first two components in the distributions (22), from the transformed solution (16)–(18) it follows that:

$$\overline{\Theta}^\*(\zeta, p) \cong \frac{e^{-\gamma^\* \zeta}}{\gamma^\*} \left[ \frac{1}{p} - \left( 1 - \frac{1}{2} e^{-\gamma^\* \zeta} \right) \frac{1}{(\gamma^{\*2} + p)} \right], \zeta \ge 0. \tag{38}$$

Considering the relations [38] as follows:

$$L^{-1}[p^{-1}; \pi] = 1,\\ L^{-1}[\left(\gamma^{\*2} + p\right)^{-1}; \pi] = e^{-\gamma^{\*2}\tau},\tag{39}$$

from the solution (38), the following asymptotic representation for the dimensionless temperature rise at large values of time was obtained:

$$\Theta^\*(\zeta,\tau) \cong \frac{e^{-\gamma^\*\zeta}}{\gamma^\*} \left[ 1 - \left( 1 - \frac{1}{2} e^{-\gamma^\*\zeta} \right) e^{-\gamma^\*\tau} \right], \zeta \ge 0, \ \tau >> 1. \tag{40}$$

The temperature rise of the heated surface was found, substituting to the formula (40) *ζ* = 0 as follows:

$$\Theta^\*(\tau) \cong \frac{1}{\gamma^\*} \left( 1 - \frac{1}{2} e^{-\gamma^\* \tau} \right), \text{ } \tau >> 1. \tag{41}$$

### **3. Heating of the FGM Semi-Space by Heat Flux with the Intensity Linearly Decreasing in Time**

As presented above, the exact (25), asymptotic (36), and (40) solutions were obtained at a constant intensity of heat flux of *q* = *q*0. This section concerns the heating of the surface of the semi-space by heat flux with the following time profile of intensity:

$$q(t) = q\_0 q^\*(t), \ q^\*(t) = 1 - \frac{t}{t\_s}, \ 0 \le t \le t\_{\text{s}} \tag{42}$$

where *ts* is the stop moment of heating. The evolution of the heat flux intensity (42) corresponds to the temporal profile of the specific friction power during braking with constant deceleration [39]. The corresponding dimensionless temperature rise of Θˆ <sup>∗</sup> was searched based on Duhamel's theorem [40] as follows:

$$\hat{\Theta}^\*(\zeta, \tau) = \frac{\partial}{\partial \tau} \int\_0^\tau q^\*(\tau - s) \Theta^\*(\zeta, s) ds, \zeta \ge 0, 0 \le \tau \le \tau\_\text{s}. \tag{43}$$

where Θ<sup>∗</sup> is the dimensionless temperature rise (25)–(27), *q*∗(*τ*) is the function (42), and *τ<sup>s</sup>* is the dimensionless stop time as follows:

$$
\pi\_s = \frac{kt\_s}{a^2}.\tag{44}
$$

After the integration with the next differentiation on the right-hand side of Formula (42), the following was obtained:

$$\hat{\Theta}^{\*}\left(\zeta,\tau\right) = \frac{1}{\gamma^{\*}} e^{-\frac{1}{2}\gamma^{\*}\overline{\zeta}} \left[ q^{\*}\left(\tau\right) e^{-\frac{1}{2}\gamma^{\*}\overline{\zeta}} - 4 \sum\_{n=1}^{\infty} \frac{f\_{1}\left(\mu\_{n}e^{-\frac{1}{2}\gamma^{\*}\overline{\zeta}}\right)}{\mu\_{n}^{2}f\_{1}\left(\mu\_{n}\right)} G\_{n}\left(\tau\right) \right], \zeta \ge 0, 0 \le \tau \le \tau\_{\text{s}} \tag{45}$$

where

$$G\_{\rm nl}(\tau) = e^{-\lambda\_{\rm ll}\tau} - \frac{\left(1 - e^{-\lambda\_{\rm tr}\tau}\right)}{\lambda\_{\rm nl}\tau\_{\rm s}},\tag{46}$$

the coefficients *λ<sup>n</sup>* were determined from Formula (26), and the numbers *μn*, *n* = 1, 2, ..., are the single positive roots of Equation (27).

By putting *ζ* = 0 in the solution (45), the following formula was used to determine the evolution of the dimensionless temperature rise of the heated surface that was obtained:

$$\hat{\Theta}^\*(\tau) = \frac{1}{\gamma^\*} \left[ q^\*(\tau) - 4 \sum\_{n=1}^\infty \frac{G\_{\mathcal{U}}(\tau)}{\mu\_n^2} \right], \ 0 \le \tau \le \tau\_s. \tag{47}$$

### **4. The Heat Partition Ratio**

We will generalize the results that were obtained above to the case of two (*l* = 1, 2) FGM semi-spaces that exponentially increased with the distance from the surface *z* = 0 of thermal conductivity as follows:

$$K\_l(z) = K\_{l,0} e^{\gamma z}, \; z \ge 0, \; K\_{l,0} \equiv K\_l(0), \; \gamma\_l \ge 0, \; l = 1, 2. \tag{48}$$

The surfaces of *z* = 0 of each semi-space were heated with the heat flux at the following intensities:

$$q\_l(t) = \mathfrak{a}\_l q(t), \ 0 \le t \le t\_{\mathfrak{s}}, \ l = 1, 2, \ \mathfrak{a}\_1 = \mathfrak{a}, \ \mathfrak{a}\_2 = 1 - \mathfrak{a},\tag{49}$$

where *q*(*t*) is the function (42) and *α* is the unknown heat partition ratio (HPR) (Figure 2).

**Figure 2.** Scheme of separation of the friction pair elements.

Taking into account the form of the solutions (45) and (46), the temperature field in each semi-space heated by heat fluxes of intensity (49) respectively, was written in the following form:

$$T\_l(z,t) = T\_0 + \Theta\_{l,0} \hat{\Theta}\_l^\*(\zeta\_{l\prime}\tau^\*), \; z \ge 0, \; 0 \le t \le t\_{\rm s\prime} \; l=1,2,\tag{50}$$

where

$$\hat{\Theta}\_{l}^{s}(\underline{\zeta}\_{l},\boldsymbol{\tau}^{\*}) = \frac{1}{\gamma\_{l}^{s}} e^{-\frac{1}{2}\gamma\_{l}^{\ast}\underline{\zeta}\_{l}} \left[ \left( 1 - \boldsymbol{\tau}^{\*} \right) e^{-\frac{1}{2}\gamma\_{l}^{\ast}\underline{\zeta}\_{l}} - 4 \sum\_{n=1}^{\infty} \frac{I\_{1}(\mu\_{n}e^{-\frac{1}{2}\gamma\_{l}^{\ast}\underline{\zeta}\_{l}})}{\mu\_{n}^{2}I\_{1}(\mu\_{n})} G\_{l,n}(\boldsymbol{\tau}^{\*}) \right], \; \zeta\_{l} \ge 0, \; 0 \le \boldsymbol{\tau}^{\*} \le 1,\tag{51}$$

$$G\_{l,n}(t^\*) = e^{-\lambda\_{l,n}\tau^\*} - \frac{(1 - e^{-\lambda\_{l,n}\tau^\*})}{\lambda\_{l,n}\tau\_{l,s}},\tag{52}$$

$$
\lambda\_{l,n} = \left(\frac{1}{2}\gamma\_l^\*\mu\_n\right)^2 \mathbf{r}\_{l,s\prime} \tag{53}
$$

$$\tau^\* = \frac{\mathbf{t}}{\mathbf{t}\_s}, \ \zeta\_l = \frac{z}{a\_l}, \ \tau\_{l,s} = \frac{k\_l t\_s}{a\_l^2}, \ \gamma\_l = \frac{\gamma\_l^\*}{a\_l}, \ k\_l = \frac{\mathbf{K}\_{l,0}}{c\_l \rho\_l}, \ \Theta\_l^s = \frac{\Theta\_l}{\Theta\_{l,0}}, \ \Theta\_{l,0} = \frac{a\_l q\_0 a\_l}{\mathbf{K}\_{l,0}}, l = 1, 2,\tag{54}$$

*cl*, *ρ<sup>l</sup>* are the specific heat capacity and the density of the materials, respectively, *μ<sup>n</sup>* > 0, *n* = 1, 2, . . . are the single roots of Equation (27).

Substituting *z* = 0 (*ζ<sup>l</sup>* = 0) to the solution (51)–(54), the change in the temperature on the heated surface was obtained in the following form:

$$T\_l(t) = T\_0 + \Theta\_{l,0} \hat{\Theta}\_l^\*(\tau^\*), \ 0 \le t \le t\_s, \ l = 1, 2,\tag{55}$$

where

$$\Theta\_l^\*(\tau^\*) = \frac{1}{\gamma\_l^\*} \left[ 1 - \tau^\* - 4 \sum\_{n=1}^\infty \frac{G\_{l,n}(\tau^\*)}{\mu\_n^2} \right], \ 0 \le \tau^\* \le 1. \tag{56}$$

The coefficients of Θ*l*,0 (54) in Formulas (50) and (55) contain an unknown heat partition ratio of *α*. In order to determine it, the equality condition of the mean temperature in time on the heated surfaces was used as follows:

$$
\tilde{T}\_1 = \tilde{T}\_{2\prime} \tag{57}
$$

where

$$\widetilde{T}\_{l} = \frac{1}{t\_{\rm s}} \int\_{0}^{t\_{\rm s}} T\_{l}(t)dt,\ l = 1,2,\tag{58}$$

Considering the temperature on the heated surfaces (55) and (56) in the Formula (58), after integration was found as follows:

$$
\bar{T}\_l = T\_0 + \Theta\_{l,0} \bar{\Theta}\_l^\* \; , \; l = 1, 2, \tag{59}
$$

where

$$\tilde{\Theta}\_l^\* = \frac{1}{2\gamma\_l^\*} \left( 1 - 8 \sum\_{n=1}^{\infty} \frac{\tilde{G}\_{l,n}}{\mu\_n^2} \right) . \tag{60}$$

$$\widetilde{G}\_{l,\rm{n}} = \frac{(1 - e^{-\lambda\_{l,\rm{n}}})}{\lambda\_{l,\rm{n}}} \left( 1 + \frac{1}{\lambda\_{l,\rm{n}} \tau\_{l,\rm{s}}} \right) - \frac{1}{\lambda\_{l,\rm{n}} \tau\_{l,\rm{s}}}, \; l = 1, 2, \tag{61}$$

and the coefficients of *λl*,*<sup>n</sup>* were determined from the formula (53).

Substituting the mean temperature (59)–(61) to Equation (57), the HPR was obtained in the following form:

$$\alpha = \frac{K^\*}{a^\* \widetilde{\Theta}^\* + K^\*} \,\tag{62}$$

where

$$a^\* = \frac{a\_1}{a\_2}, \ K^\* = \frac{K\_{1,0}}{K\_{2,0}}, \ \tilde{\Theta}^\* = \frac{\Theta\_1^\*}{\tilde{\Theta}\_2^\*}. \tag{63}$$

If the effective heating depths are determined from the Formula [9] as follows:

$$a\_l = \sqrt{3k\_l t\_s}, \ l = 1, 2,\tag{64}$$

then from Formula (54) it follows that the dimensionless heating time of each semi-space is equal to *τl*,*<sup>s</sup>* = 3−<sup>1</sup> ∼= 0.333, and from Formulas (62) and (63) it follows that:

$$\alpha = \frac{K\_{\varepsilon}}{\tilde{\Theta}^{\*} + K\_{\varepsilon}},\tag{65}$$

where

$$K\_{\varepsilon} = \frac{K^\*}{\sqrt{k^\*}}, \; k^\* = \frac{k\_1}{k\_2}. \tag{66}$$

Omitting the component containing a series, on the right-hand side of Formula (60), i.e., assuming the following:

$$
\widetilde{\Theta}\_l^\* \cong \frac{1}{2\gamma\_l^{\*,\prime}}, \ l = 1, 2,\tag{67}
$$

where

$$\alpha \cong \frac{K^\* \gamma^\*}{a^\* + K^\* \gamma^\*},\tag{68}$$

and Formula (65) is as follows:

$$\alpha \cong \frac{K\_{\varepsilon} \gamma^\*}{1 + K\_{\varepsilon} \gamma^\*} \tag{69}$$

where

$$
\gamma^\* = \frac{\gamma\_1^\*}{\gamma\_2^\*}.\tag{70}
$$

With the same effective heating depths (*a*∗ = 1) and dimensionless gradients of materials (*γ*∗ = 1) from Formula (68), the known Block' result [20] as follows:

$$\alpha \cong \frac{K^\*}{1 + K^\*},\tag{71}$$

and from Formula (69), the classic Charron's Formula [41] was achieved as follows:

$$\alpha \cong \frac{K\_{\varepsilon}}{1 + K\_{\varepsilon}} \, \tag{72}$$

which is often used in analytical and numerical modeling of the frictional heating of homogeneous materials during braking [42–44]. For the same materials of both semi-spaces (*K*∗ = *k*∗ = *γ*∗ = *a*∗ = 1), all of the above obtained formulas give the value of the heat partition ratio equal to *α* = 0.5.

It should be noted that the received formulas for the determination of HPR need to be verified in the future with appropriate experimental data. Due to the lack of such opportunity at that moment, the authors would be grateful for the provision of such data or for carrying out the cooperative research. This would also allow us to establish the limits of the applicability of the formulas. However, already at this stage of research, we can assume the possibility of the practical use of the obtained formulas to determine HPR, due to the usage of the classic methodology that has been approved for homogeneous materials. Furthermore, the experimentally confirmed formulas for homogeneous materials can be obtained from the proposed solution for FGMs, as a results of the proper limit approach.

An important element of the most general Formula (62) for determining the HPR are the dimensionless time-averaged temperature rises Formulas (60) and (61) found in the case of the heat flux intensity linearly decreasing in time. This case is often considered when calculating the temperature of the brake systems that are operating in the mode with a sudden increase in the contact pressure to the nominal value at the beginning of braking. A classification of the remaining heat flux intensity temporal profiles has been proposed for homogeneous materials, without and with consideration of the contact pressure rise [45,46]. Obtaining the appropriate solutions in the case of FGM is also one of the directions of our research in the future.

#### **5. Example of Calculation of the Heat Partition Ratio for an FGM Couple**

Calculations were performed for two semi-limited bodies made of two-component FGMs. The first (*l* = 1) element forms aluminum oxide Al2O3 (base, *m* = 0) and copper Cu (core, *m* = 1), and the other (*l* = 2) contains zircon dioxide ZrO2 (base, *m* = 0) and titanium alloy Ti-6Al-4V (core, *m* = 1). The thermo-physical properties of these materials are demonstrated in Table 1.


**Table 1.** Thermo-physical properties of FGM components [32].

The other operating parameters are the nominal heat flux intensity of *q*<sup>0</sup> = 3.78 MW m<sup>−</sup>2, the braking time of *ts* = 121 s, and the temperature at the initial time moment of *T*<sup>0</sup> = 20 ◦C [32]. The dimensionless material gradient parameters of each element were calculated from the relation [47] as follows:

$$\gamma\_l^\* = \ln \left( \frac{K\_{l,1}}{K\_{l,0}} \right), \ l = 1, 2. \tag{73}$$

Obtaining the values of *γ*∗ <sup>1</sup> = 2.38, *γ*<sup>∗</sup> <sup>2</sup> = 1.26. The specific heat capacity and density of materials of the heated elements were determined according to the following mixture law:

$$\mathcal{L}\_l = \mathcal{c}\_{l,1}v + (1-v)\mathcal{c}\_{l,0\prime}\,\rho\_l = \rho\_{l,1}v + (1-v)\rho\_{l,0\prime}\,l = 1,2,\tag{74}$$

where *v* is the volume fractions of the base and core components. For their equal participation (*v* = 0.5), the following properties of FGM were found:

$$c\_1 = 437.32 \text{ J/kg}^{-1} \text{K}^{-1}, \,\rho\_1 = 6469.42 \text{ kgm}^{-3}, \,\text{k}\_1 = 1.32 \cdot 10^{-5} \,\text{m}^2\text{s}^{-1}, \tag{75}$$

$$\sigma\_2 = 495.46 \text{ J/kg}^{-1} \text{K}^{-1}, \,\rho\_2 = 5266.96 \text{ kgm}^{-3}, \,\text{k}\_2 = 7.43 \cdot 10^{-7} \,\text{m}^2\text{s}^{-1}.\tag{76}$$

From Equation (64) the effective depths of the heat penetration were determined as *a*<sup>1</sup> = 21.8 mm, *a*<sup>2</sup> = 5.2 mm. Next, based on Formulas (63), (66), and (70), the values of the dimensionless parameters were calculated as follows: *a*∗ = 4.208, *K*∗ = 19.196, *k*∗ = 17.706, *K<sup>ε</sup>* = 4.562, *γ*∗ = 1.883. This allowed for an estimation by means of Equations (60), (61) dimensionless, time-averaged temperature rises of Θ<sup>∗</sup> <sup>1</sup> = 0.329, Θ<sup>∗</sup> <sup>2</sup> = 0.703, and their ratio of Θ<sup>∗</sup> = 0.468 (63).

By substituting the found parameter values successively to the right-hand sides of Formulas (65) and (69), the proper values of the HPR were calculated as *α* = 0.907 and *α* = 0.896. A slight (1.2%) difference in the obtained results allowed us to analyze the influence of the dimensionless parameters of the thermal activity *Kε* (66) and the relative gradient of the FGMs *γ*∗ (70) on the HPR value (Figure 3), on the basis of Formula (69) only. For a fixed value of the thermal activity coefficient *Kε*, a rise of the parameter *γ*∗ causes an increase in the amount of heat that is directed to the first element of the heated couple. Conversely, by increasing the thermal activity of the friction pair at a predetermined value *γ*∗, the amount of heat that is directed to the first element increases.

For the HPR values *α*<sup>1</sup> = 0.896, *α*<sup>2</sup> = 0.104 (49), from Equation (54), the scaling factors of the temperature rises were determined to be Θ1,0 = 2034 ◦C and Θ2,0 = 1080 ◦C. Then, based on the relations Formulas (55) and (56) the evolutions of the temperature on the heated surfaces of each element were found (Figure 4).

**Figure 3.** Isolines of the heat partition ratio *α* (69) in the coordinate system (*γ*∗, *Kε*).

**Figure 4.** Evolutions of temperature *Tl*, *l* = 1, 2 (55) and (56) on the heated surfaces *z* = 0 of the considered semi-spaces made of FGMs.

It has been established that the temporal profiles of the temperature of the elements are different, in particular in the final stage of the heating process. The time course of the temperature of the first element (*l* = 1, Al2O3–Cu) is typical for the evolution of the friction surface temperature during braking with a constant deceleration—a rapid increase in the temperature at the beginning of braking, reaching its maximum value in the middle of the process, followed by a temperature reduction until the standstill. However, in the second element (*l* = 2, ZrO2–Ti–6Al–4V) a rise of temperature on the heated surface is monotonic during the whole process. Such temperature behavior is decisively influenced by the thermo-physical properties of the component materials of each element. In the functionally graded friction couple under consideration, the materials of both of the components of the first element have a significantly greater ability to dissipate the heat from the heated surface than the materials of the second element (Table 1). Moreover, this is confirmed, by the values of the coefficients *α*<sup>1</sup> = 0.896, *α*<sup>2</sup> = 0.104, which prove that the first element

absorbs almost 90%, and the second only slightly more than 10% of the entire heat flux intensity *q* (42). Due to the low thermal conductivity of the component materials of the second element, in particular the zircon dioxide, the temperature of the heated surface of this element continues to rise during heating, even with a linearly decreasing intensity of the heat flux.

The results that are demonstrated in Figures 3 and 4 were obtained by using the HPR value *α* = 0.896 that was found from the dependency (69). The approximated formula for determining the *α* (69), as well as the exact Equation (65), was obtained by assuming that the effective depths of the heating of the elements *al*, *l* = 1, 2 were determined from the empirical Formula (64), so the corresponding values of the Fourier numbers *τl*,*s*, *l* = 1, 2 (54) were the same and equal to 1/3. The introduction of Formula (64) to the model is justified when the thickness of the heated element is greater than the appropriate value *al*, *l* = 1, 2 [2]. In general, for the determination of the temperature mode at the design stage of the brake, the Fourier numbers *τl*,*s*, *l* = 1, 2 are given, and the effective heating depths based on the relation (54) are calculated in the following form:

$$a\_l = \sqrt{\frac{k\_l t\_s}{\tau\_{l,s}}}, l = 1, 2. \tag{77}$$

With consideration of Equation (75) for the selected friction couple and the established time of heating (braking) *ts* = 12.1 s the ratio of the depths of heating *a*<sup>∗</sup> (63) is determined from the following relation:

$$a^\* = \sqrt{\frac{k^\*}{\tau\_s^{\*'}}}, \ \tau\_s^\* = \frac{\tau\_{1,s}}{\tau\_{2,s}}.\tag{78}$$

where *k*∗ is the ratio of the thermal diffusivities of the FGMs (66). By treating the Fourier numbers *τl*,*s*, *l* = 1, 2 (54) as independent variables in Formula (62), their influence on the HPR *α* was investigated with the previously determined values of the dimensionless parameters, *K*∗ = 19.196, *k*∗ = 17.706, *γ*∗ <sup>1</sup> = 2.38, and *γ*<sup>∗</sup> <sup>2</sup> = 1.26 for a given pair of elements. The dimensionless parameter Θ<sup>∗</sup> (63) included in Formula (62) was determined from Formulas (60) and (61), and the finding of *a*<sup>∗</sup> was made by using Formula (76). The results of the calculations are presented in Figure 5.

**Figure 5.** Dependences of the heat partition ratio *α* (62) on: (**a**) Fourier number *τ*1,*<sup>s</sup>* for different values of *τ*2,*s*; (**b**) Fourier number *τ*2,*<sup>s</sup>* for different values of *τ*1,*s*.

At a fixed value *τ*2,*s*, the HPR *α* quickly increases with increasing *τ*1,*s*, reaching its nominal value at *τ*1,*<sup>s</sup>* ∼= 1.5 (Figure 5a). The largest nominal value *α* is achieved for *τ*1,*<sup>s</sup>* = 1/3, i.e., when determining the effective heating depth using Formula (64). The results obtained

at *τ*2,*<sup>s</sup>* = 1 and *τ*2,*<sup>s</sup>* = 2 practically coincide. A different nature of the HPR *α* change occurs with a fixed value of the Fourier number *τ*1,*<sup>s</sup>* and increasing the values of *τ*2,*<sup>s</sup>* (Figure 5b). By increasing *τ*2,*s*, the amount of heat that is absorbed by the first element decreases, reaching a minimum value at *τ*2,*<sup>s</sup>* ∼= 1.5. It should be noted that the Fourier numbers *τl*,*s*, *l* = 1, 2 (54) are dimensionless input parameters that play an important role in estimating the temperature mode of the brakes. The designers calculate their values after determining the geometric parameters *al*, *l* = 1, 2 (63). In the quality of these parameters, they choose the thickness of the friction elements, or, if they are greater than the effective depth of the heat penetration (64), then the latter. Then, having previously found the maps of those, as shown in Figure 5, the designer can easily estimate the HPR of a given braking system.

#### **6. Conclusions**

A methodology was proposed in order to determine the heat partition ratio (HPR) in the braking systems with friction pair elements that were made of functionally graded materials (FGMs). The basis of this methodology is an exact solution to the initial-boundary problem of heat conduction that is formulated for the FGM semi-space that is heated on its surface by a heat flux with an intensity that decreases linearly with time and includes an unknown a priori HPR. For its determination, the condition of the equality of the timeaveraged temperature on the heated surfaces of the two different semi-spaces was used.

As an example, the heat partition ratio was found for the FGM friction couple that consisted of Al2O3–Cu and ZrO2–Ti–6Al–4V. It was established that element Al2O3–Cu absorbs most of the heat that is generated due to friction (almost 90%). The maximum temperature on the friction surface of the Al2O3–Cu element is about 820 ◦C and is achieved in the middle of the braking time. However, the highest temperature of the friction surface of the ZrO2–Ti–6Al–4V element is achieved at the stop moment and amounts to 970 ◦C. Thus, the crucial influence on the evolution of the temperature in the functionally graded friction couple have the thermo-physical properties of the component materials.

Then, the simulation of the HPR dependency on the dimensionless parameters, such as the thermal activity of the friction couple, the gradient materials, and the Fourier numbers was carried out. It was found that increasing the values of the friction couple thermal activity or the gradient materials ratio causes an increase in the HPR. The scheme of the map development for a given functionally graded friction pair was proposed, which allows for a quick estimation of the HPR, depending on such input parameters as their Fourier numbers. The importance of having such maps for the designer is indicated by the results that have been obtained in the case of the homogeneous materials [9,28].

**Author Contributions:** Conceptualization and methodology, A.Y.; software, K.T. and P.Z.; validation, P.Z.; formal analysis, K.T. and P.Z.; investigation, A.Y., K.T. and P.Z.; writing—original draft preparation, K.T. and P.Z.; writing–review and editing, A.Y., K.T. and P.Z.; visualization and figures preparation, K.T.; supervision, A.Y.; project administration, K.T.; funding acquisition, K.T. All authors have read and agreed to the published version of the manuscript.

**Funding:** This investigation was performed within the framework of research project No. 2017/27/B/ ST8/01249, funded by the National Science Centre, Poland, and project financing was obtained through the program of the Minister of Education and Science of Poland named "Regional Initiative of Excellence" in 2019–2022, project No. 011/RID/2018/19, amount financing 12,000,000 PLN.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** No new data were created or analyzed in this study. Data sharing is not applicable to this article.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Nomenclature**


### **References**

