Frictional Heating during Braking of the C/C Composite Disc

An analytical model to determine temperature in a single brake disc of multi-disc system is proposed. The model considers the convective cooling on the lateral surfaces of the disc and structure of composite friction material. Calculations were carried out for a disc made of carbon friction composites material Termar-ADF. The influence of heat transfer with environment, length of bundles with fibers, and concentration of fibers in composite on the temperature of the disc was investigated during single braking with constant deceleration.


Introduction
At the present time, high-speed and heavy-duty vehicles are extensively developed. These include high-speed rail and road transport, machines for military purposes, and, in particular, aircraft. For the successful operation of these categories of vehicles, appropriate braking devices are required that ensure their sufficient maneuverability. To design brakes, it is necessary to fulfill a complex of requirements (ensuring high braking torque and wear resistance, absence of overheating, stability of characteristics, etc.), which is associated with the choice of friction materials, design of the brake unit, bench, and field tests [1].
The braking system of modern passenger as well as military aircraft consists of package of fixed and rotating friction discs located inside the wheel. After switching on the brake, the special pistons compress package, the front surfaces of discs are in contact, and the process of friction braking begins. While the kinetic energy of the aircraft is absorbed, the friction between fixed and rotating discs causes that they are heated up even to 1500 • C, and temperatures of the friction surfaces reach up to 3000 • C. Therefore, the disc's materials must be resistant enough to thermal shock [2].
It is known that the friction polymer materials start to degrade at the temperature of 230 • C, and the degree of degradation increases with temperature within the range 269-400 • C [3,4]. The degradation of the polymer materials causes brake fade phenomena, where the coefficient of friction is reduced with increasing temperature [5]. The high temperature decreases the yield strength and leads to changes in the wear mechanism and the real contact configuration [6]. These phenomena could increase the wear rate of the brake friction material [7]. Therefore, to improve the working life of aircraft brakes, there is an increasing use of carbon friction composites materials (CFCM) that differ in the type of reinforcement, the method of forming the friction layer, the variations of technological heat treatment, etc. For such a material, a certain degree of wear intensity is typical, but, in general, CFCM has significant advantages over friction metal ceramics and plastics. CFCM materials are characterized by less tendency to increase wear with rise in temperature, as well as smaller sensitivity to changes in the specific pressure and sliding speed [8][9][10].
For a multi-disc brake pair made of composite C/C materials, a numerical finite element method (FEM) solution of thermal problem of friction is obtained in [29], and an analytical solution is presented in [30]. That solutions are obtained considering effective properties of disc composite materials. This article provides analytical solution to the thermal problem of friction, considering the structure of the composite material. In detail, the main aim of the article is to obtain an exact solution to the thermal problem of friction, which simulates the process of frictional heating of an individual disc of a multi-disc brake during single braking. It is assumed that the disc is made of composite C/C material containing carbon fibers, collected in bundles, which are located in planes parallel to the surface of friction. At the first stage, an exact solution of the corresponding thermal problem of friction for an orthotropic disc is obtained, considering the effective thermal conductivity coefficients of the material in the axial and radial directions, as well as convective cooling of the lateral surfaces. In the second stage, the procedure for determining the above-mentioned effective properties of the composite at three levels (micro, meso, and macro) is proposed. Based on the obtained solution, the distribution of a transient temperature field in a disc made of Termar-ADF is investigated.

Statement to the Problem
Multi-disc aircraft braking systems are subjected to heavy mechanical and thermal loads during braking. Typically, such systems are a package of identical discs. A scheme of such multi-disc system and the calculation diagram of the considered problem are presented in Figure 1. A microscopic view of the real structure of such C/C composite can be found in [15,20] and sample photo of single brake disc made of Termar-ADF in [31].
Materials 2020, 13,2691 3 of 20 For a multi-disc brake pair made of composite C/C materials, a numerical finite element method (FEM) solution of thermal problem of friction is obtained in [29], and an analytical solution is presented in [30]. That solutions are obtained considering effective properties of disc composite materials. This article provides analytical solution to the thermal problem of friction, considering the structure of the composite material. In detail, the main aim of the article is to obtain an exact solution to the thermal problem of friction, which simulates the process of frictional heating of an individual disc of a multi-disc brake during single braking. It is assumed that the disc is made of composite C/C material containing carbon fibers, collected in bundles, which are located in planes parallel to the surface of friction. At the first stage, an exact solution of the corresponding thermal problem of friction for an orthotropic disc is obtained, considering the effective thermal conductivity coefficients of the material in the axial and radial directions, as well as convective cooling of the lateral surfaces. In the second stage, the procedure for determining the above-mentioned effective properties of the composite at three levels (micro, meso, and macro) is proposed. Based on the obtained solution, the distribution of a transient temperature field in a disc made of Termar-ADF is investigated.

Statement to the Problem
Multi-disc aircraft braking systems are subjected to heavy mechanical and thermal loads during braking. Typically, such systems are a package of identical discs. A scheme of such multi-disc system and the calculation diagram of the considered problem are presented in Figure 1. A microscopic view of the real structure of such C/C composite can be found in [15,20] and sample photo of single brake disc made of Termar-ADF in [31]. It is assumed that all friction discs have the same internal radius 1 r , external radius 2 r , thickness d 2 , and are made of orthotropic material with different thermal conductivities r K and z K in radial r and axial z directions, respectively. Each of the two outer fixed discs has one friction surface, while the inner rotating and fixed discs have two friction surfaces, and thus they heat up more. Since all the discs are made of the same material, the temperature fields in the internal discs are approximately the same, and their distribution over the thickness of the disc is symmetrical relative to its middle surface. Therefore, the corresponding thermal problem of friction can be formulated not for a whole disc with a thickness of d 2 , but only for its half with a thickness of d , assuming that the friction surface 0 = z is heated by frictional heat flow with the intensity q and the plane of symmetry d z = is adiabatic. Since the thermal conductivity of CFCM in the circumferential direction parallel to the working surfaces of the disc is much higher than in the axial direction z of action of the heat flow, the distribution of the temperature field on each surface d z < < 0 is uniform. Lateral surfaces of the discs  It is assumed that all friction discs have the same internal radius r 1 , external radius r 2 , thickness 2d, and are made of orthotropic material with different thermal conductivities K r and K z in radial r and axial z directions, respectively. Each of the two outer fixed discs has one friction surface, while the inner rotating and fixed discs have two friction surfaces, and thus they heat up more. Since all the discs are made of the same material, the temperature fields in the internal discs are approximately the same, and their distribution over the thickness of the disc is symmetrical relative to its middle surface. Therefore, the corresponding thermal problem of friction can be formulated not for a whole disc with a thickness of 2d, but only for its half with a thickness of d, assuming that the friction surface z = 0 is heated by frictional heat flow with the intensity q and the plane of symmetry z = d is adiabatic. Since the thermal conductivity of CFCM in the circumferential direction parallel to the working surfaces of the disc is much higher than in the axial direction z of action of the heat flow, the distribution of the temperature field on each surface 0 < z < d is uniform. Lateral surfaces of the discs r = r 1 and r = r 2 0 ≤ z ≤ 2d are cooled convectively with constant coefficient of heat transfer h. The mutual overlap of the contacting elements is full, so the nominal contact area A a = π(r 2 2 − r 2 1 ) remains unchanged during braking.
Taking the above-mentioned assumptions into account, distribution of the temperature T across the disc thickness 0 ≤ z ≤ d at a fixed moment of time 0 ≤ t ≤ t s , we may find from solution to the following boundary-value problem of heat conduction [28]: where T a is the initial temperature of the disc; t is the time; t s is the time of braking; f is the coefficient of friction; ρ is the density; c is the specific heat of material; p 0 , ω 0 , and W 0 are respectively initial values of contact pressure, angular velocity, and kinetic energy.
Introducing the dimensionless variables and parameters: We write the problem in Equations (1)-(4) as: where It should be noted that the adoption of the linear model in Equations (1)- (9) to determine the temperature field in the brake disc, imposes some restrictions on the possibility of its use in medium and heavy modes of operation of the brake, when the surface temperature reaches values above 450 and 750 • C, respectively [3]. It is known that at such temperatures the coefficient of friction and the thermo-physical properties of the material may change the initial values. Application of the analytical models in such cases may involve the use of the friction coefficient averaged during braking; this approach has been used in the following numerical analysis. Thermal sensitivity of a material is usually taken into account by converting in its solution the value of the coefficient of thermal conductivity at room temperature to its value at the average volume temperature of the disc during braking. Of course, this approach is approximate, but gives no less sufficiently good compliance of the maximum temperature values with relevant experimental data [21]. To consider the thermal sensitivity of the above-mentioned characteristics, at each time step during braking, requires the development of appropriate nonlinear models, and the usage of numerical methods to solve them [29].
Taking into account the sums of the series [35]: From Equations (16) and (24)- (27), we find the dimensionless temperature of the disc in the form: Bi n cos(λ n ζ), On the heated surface ζ = 0 from the solution in Equation (30), it follows that [27]: where the coefficients Bi n , n = 1, 2, . . . . are determined by Equation (27).
In the limiting case τ s → ∞ , the solutions in Equations (30) and (31) take the form: Next, we verify the achieved solutions. For this purpose, from the solution in Equation (30), we find the derivative: Bi n λ n sin(λ n ζ), On the heated ζ = 0 and adiabatic ζ = 1 surfaces of the disc, from the relation in Equation (34), it follows that: i.e., the solution in Equation (30) satisfies the boundary conditions in Equations (12) and (13). Substituting τ = 0 in Equation (30), we get: Using the sums (28) and (29), we have: Introducing the sum in Equation (38) into Equation (37), we confirms the execution of the initial condition in Equation (14).
Finally, taking into account the notation in Equation (10), we obtain the temperature field in the disc: (30)).
With account of the sums [35] ∞ n=1 cos(λ n ζ) the solution in Equations (40)-(42) takes the form: Substituting ζ = 0 in the solution in Equation (45), we find known formula to calculate the dimensionless temperature of the heated surface of the adiabatic disc [27]:

Effective Thermal Conductivity Coefficients
Let us consider the brake disc made from CFCM Thermal-ADF. It consists of the bundles with fibers incorporated into the matrix. Cylindrical shaped bundles are arranged in layers parallel to the front surfaces of the disc. Effective values of thermal conductivities of such composite K z and K r are determined by means of averaging at the three scale levels: micro, meso, and macro [28].
At the micro-level, a single bundle is considered, which consists of closely spaced fibers embedded in the matrix material ( Figure 2).

Effective Thermal Conductivity Coefficients
Let us consider the brake disc made from CFCM Thermal-ADF. It consists of the bundles with fibers incorporated into the matrix. Cylindrical shaped bundles are arranged in layers parallel to the front surfaces of the disc. Effective values of thermal conductivities of such composite z K and r K are determined by means of averaging at the three scale levels: micro, meso, and macro [28]. At the micro-level, a single bundle is considered, which consists of closely spaced fibers embedded in the matrix material ( Figure 2). The effective coefficients of thermal conductivity of this bundle in the transversal direction (across the fibers) ⊥ K and in longitudinal direction (along the fibers) IÎ K are established using the mixture rule in the following form [36]: where f K and m K are the thermal conductivities of fiber material and matrix material, respectively, and A unit cell of composite at the meso-level has the shape of a rectangular cuboid with length A and squared cross-section with width B . This cell contains a single centered bundle in form of squared cuboid with dimensions a and b , respectively ( Figure 3). It should be noted that the real circular shape of bundle cross-section has been replaced by a square to simplify further calculations, and it does not introduce significant errors in calculations [28]. The effective coefficients of thermal conductivity for selected unit cell of composite in longitudinal II K and transversal ⊥ K directions are determined using the method of heat flux lines linearization [37]. This method consists in the division of the unit cell of composite by adiabatic or isothermal and adiabatic planes and calculating the thermal resistance of each part: The effective coefficients of thermal conductivity of this bundle in the transversal direction (across the fibers)K ⊥ and in longitudinal direction (along the fibers)K II are established using the mixture rule in the following form [36]: where K f and K m are the thermal conductivities of fiber material and matrix material, respectively, and 0 ≤ V f ≤ 1 is the fiber volume ratio. A unit cell of composite at the meso-level has the shape of a rectangular cuboid with length A and squared cross-section with width B. This cell contains a single centered bundle in form of squared cuboid with dimensions a and b, respectively ( Figure 3). It should be noted that the real circular shape of bundle cross-section has been replaced by a square to simplify further calculations, and it does not introduce significant errors in calculations [28].

Effective Thermal Conductivity Coefficients
Let us consider the brake disc made from CFCM Thermal-ADF. It consists of the bundles with fibers incorporated into the matrix. Cylindrical shaped bundles are arranged in layers parallel to the front surfaces of the disc. Effective values of thermal conductivities of such composite z K and r K are determined by means of averaging at the three scale levels: micro, meso, and macro [28].
At the micro-level, a single bundle is considered, which consists of closely spaced fibers embedded in the matrix material ( Figure 2). The effective coefficients of thermal conductivity of this bundle in the transversal direction (across the fibers) ⊥ K and in longitudinal direction (along the fibers) IÎ K are established using the mixture rule in the following form [36]: where f K and m K are the thermal conductivities of fiber material and matrix material, respectively, and A unit cell of composite at the meso-level has the shape of a rectangular cuboid with length A and squared cross-section with width B . This cell contains a single centered bundle in form of squared cuboid with dimensions a and b , respectively ( Figure 3). It should be noted that the real circular shape of bundle cross-section has been replaced by a square to simplify further calculations, and it does not introduce significant errors in calculations [28]. The effective coefficients of thermal conductivity for selected unit cell of composite in longitudinal II K and transversal ⊥ K directions are determined using the method of heat flux lines linearization [37]. This method consists in the division of the unit cell of composite by adiabatic or isothermal and adiabatic planes and calculating the thermal resistance of each part: The effective coefficients of thermal conductivity for selected unit cell of composite in longitudinal K II and transversal K ⊥ directions are determined using the method of heat flux lines linearization [37]. This method consists in the division of the unit cell of composite by adiabatic or isothermal and adiabatic planes and calculating the thermal resistance of each part: where l is the distance traveled by a heat flux Q, S is the area of cross-section of the part perpendicular to the heat flux direction, and K is the thermal conductivity of the material of the same part ( Figure 4). Bearing in mind that the results depend on the manner of the unit cell division, we find effective thermal conductivities of cell using two methods of division by means of adiabatic and isothermal planes. The final results are established as an arithmetic average of two outcomes obtained by different methods. Below, we consider these two methods.
where l is the distance traveled by a heat flux Q , S is the area of cross-section of the part perpendicular to the heat flux direction, and K is the thermal conductivity of the material of the same part ( Figure 4). Bearing in mind that the results depend on the manner of the unit cell division, we find effective thermal conductivities of cell using two methods of division by means of adiabatic and isothermal planes. The final results are established as an arithmetic average of two outcomes obtained by different methods. Below, we consider these two methods.
Division of a unit cell by adiabatic planes 1-1 and 2-2 in the transverse direction is shown in Figure 4a, and in the longitudinal direction in Figure 4b.
where, according to the scheme of combining resistances i R , (Figure 5a), the total thermal resistances ⊥ R and II R of the unit cell are equal to: II  1  ,  II  3  ,  II  3  ,  II  2  ,  II  1  ,  II   3  ,  II  2  ,  II  1  ,  II  II   2   1 1 2  Division of a unit cell by adiabatic planes 1-1 and 2-2 in the transverse direction is shown in Figure 4a, and in the longitudinal direction in Figure 4b.
Based on Equation (48), the effective coefficients of thermal conductivity of selected unit cell of composite in this case are determined in the form: where, according to the scheme of combining resistances R i , i = 1, 2, 3 of individual part (Figure 5a), the total thermal resistances R ⊥ and R II of the unit cell are equal to: and the thermal conductivitiesK II andK ⊥ of the bundle contained in the unit cell are calculated using Equation (47). Substituting the effective thermal resistance ⊥ R (Equations (50) and (51)) and II R (Equations (53) and (540) into Equation (49), we find: Division of the unit cell simultaneous by adiabatic 1-1 and 2-2 and isothermal 3-3 and 4-4 planes, respectively, in transverse and longitudinal directions, are presented in Figure 6. A scheme of combining resistances for this division is shown in Figure 5b. Proceeding as above, we obtain:  Substituting the effective thermal resistance R ⊥ (Equations (50) and (51)) and R II (Equations (53) and (540) into Equation (49), we find: Division of the unit cell simultaneous by adiabatic 1-1 and 2-2 and isothermal 3-3 and 4-4 planes, respectively, in transverse and longitudinal directions, are presented in Figure 6. A scheme of combining resistances for this division is shown in Figure 5b. Substituting the effective thermal resistance ⊥ R (Equations (50) and (51)) and II R (Equations (53) and (540) into Equation (49), we find: (55) Division of the unit cell simultaneous by adiabatic 1-1 and 2-2 and isothermal 3-3 and 4-4 planes, respectively, in transverse and longitudinal directions, are presented in Figure 6. A scheme of combining resistances for this division is shown in Figure 5b. Proceeding as above, we obtain:  Proceeding as above, we obtain: where thermal resistances R ⊥,2 and R II,2 are determined using Equations (51) and (53). Taking into account the effective thermal resistances R ⊥ (Equations (56) and (57)) and R II (Equations (58) and (59)) in Equation (49), we find: The length A and the width B of the unit cell of composite are related with dimensions a and b of the single bundle and the bundle volume ratio 0 ≤ V b ≤ 1 as: Eliminating A in Equation (62), we obtain the following algebraic equation for parameter B: At the macro-level, the fiber bundles are arranged on the planes parallel to the front surfaces of the disc (Figure 7). Materials 2020, 13 , 2691  11 of 20   1   5  ,  II  II,2   4  ,  II  5  ,  II  2  ,  II   5  ,  II  2  ,  II  4  ,  II  II   1  1  2 where thermal resistances (59)) in Equation (49), we find: The length A and the width B of the unit cell of composite are related with dimensions a and b of the single bundle and the bundle volume ratio Eliminating A in Equation (62), we obtain the following algebraic equation for parameter B : At the macro-level, the fiber bundles are arranged on the planes parallel to the front surfaces of the disc (Figure 7). The thermal conductivity of the disc in the axial direction is equal to the effective coefficient of thermal conductivity of the unit composite cell in the transverse direction on the meso-level: where ⊥ K is a mean arithmetic of values calculated from Equations (54) and (60).
When determining the effective coefficient of thermal conductivity of the disc in the radial direction, the orientation of fiber bundles on the planes should be taken into account (Figure 8). The thermal conductivity of the disc in the axial direction is equal to the effective coefficient of thermal conductivity of the unit composite cell in the transverse direction on the meso-level: where K ⊥ is a mean arithmetic of values calculated from Equations (54) and (60). When determining the effective coefficient of thermal conductivity of the disc in the radial direction, the orientation of fiber bundles on the planes should be taken into account (Figure 8). If the bundles are arranged in the radial direction (Figure 8a), then: where II K is a mean arithmetic values obtained from Equations (55) and (61). However, if the bundles are oriented in circumferential direction (Figure 8b), then: where ⊥ K is determined as in Equation (64). In the case of composite discs with chaotically oriented fiber bundles (Figure 8c), the effective coefficient of thermal conductivity in the radial direction is defined as: where II K and ⊥ K are calculated in the same manner as in Equations (65) and (66), respectively.
The empirical Equation (67) assumes that the number of bundles oriented in the radial and circumferential directions are equal. Therefore, the value of r K found using Equation (67) for the chaotic orientation bundles may be established with an error. One of the proposals to modify Equation (67) could be to use statistical methods with the known probability distribution of bundles orientation.

Numerical Analysis
Calculations of the temperature T (Equation (39)) were performed for a disc made of carbon composite material − Termar-ADF with randomly oriented fiber bundles ( ) [27]. The values of others input parameters were equal to: , and C 20°= a T [28].
The change of temperature during braking on the friction surface of the disc is presented in Figure 9. Solid line demonstrates the results obtained using Equations (31) and (39) with effective values of thermal conductivities z K (Equation (64)) and r K (Equation (67)). Dotted line presents appropriate data obtained from measurement of temperature by means of thermocouples [27]. It can be seen that the curves are close to each other. In the initial stage of braking, the temperature achieved in the theoretical solution rapidly increases to maximum value C 466 . The highest value of the experimentally measured temperature is almost the same; however, it is attained a bit later at s 4 max ≈ = t t . In the final stage of braking, the experimental temperature is slightly higher than the value obtained using the theoretical solution. If the bundles are arranged in the radial direction (Figure 8a), then: where K II is a mean arithmetic values obtained from Equations (55) and (61). However, if the bundles are oriented in circumferential direction (Figure 8b), then: where K ⊥ is determined as in Equation (64). In the case of composite discs with chaotically oriented fiber bundles (Figure 8c), the effective coefficient of thermal conductivity in the radial direction is defined as: where K II and K ⊥ are calculated in the same manner as in Equations (65) and (66), respectively. The empirical Equation (67) assumes that the number of bundles oriented in the radial and circumferential directions are equal. Therefore, the value of K r found using Equation (67) for the chaotic orientation bundles may be established with an error. One of the proposals to modify Equation (67) could be to use statistical methods with the known probability distribution of bundles orientation.
The change of temperature during braking on the friction surface of the disc is presented in Figure 9. Solid line demonstrates the results obtained using Equations (31) and (39) with effective values of thermal conductivities K z (Equation (64)) and K r (Equation (67)). Dotted line presents appropriate data obtained from measurement of temperature by means of thermocouples [27]. It can be seen that the curves are close to each other. In the initial stage of braking, the temperature achieved in the theoretical solution rapidly increases to maximum value T max = 466 • C at t = t max = 3.36 s, whereupon it decreases to T = 337 • C at the stop t = t s = 6.8 s. The highest value of the experimentally measured temperature is almost the same; however, it is attained a bit later at t = t max ≈ 4 s. In the final stage of braking, the experimental temperature is slightly higher than the value obtained using the theoretical solution. Changes of temperature during braking on different distances from heated surface of the disc are shown in Figure 10. Calculations were conducted based on the exact solution in Equations (30) and (39). The temperature growth at the very beginning of the process is the fastest, and its maximum value is the highest on the friction surface 0 = z . Increasing the distance from this surface to the inside of the disc, the temperature decreases and the delay effect occurs-elongation of the time max t of maximum temperature achievement. After achieving the highest temperature value at the particular depths in the range mm 6 0 < ≤ z , the temperature drop lasts until standstill, whereas, , the temperature monotonically rises during the whole braking process.  Changes of temperature during braking on different distances from heated surface of the disc are shown in Figure 10. Calculations were conducted based on the exact solution in Equations (30) and (39). The temperature growth at the very beginning of the process is the fastest, and its maximum value is the highest on the friction surface z = 0. Increasing the distance from this surface to the inside of the disc, the temperature decreases and the delay effect occurs-elongation of the time t max of maximum temperature achievement. After achieving the highest temperature value at the particular depths in the range 0 ≤ z < 6 mm, the temperature drop lasts until standstill, whereas, for z ≥ 6 mm, the temperature monotonically rises during the whole braking process. Changes of temperature during braking on different distances from heated surface of the disc are shown in Figure 10. Calculations were conducted based on the exact solution in Equations (30) and (39). The temperature growth at the very beginning of the process is the fastest, and its maximum value is the highest on the friction surface 0 = z . Increasing the distance from this surface to the inside of the disc, the temperature decreases and the delay effect occurs-elongation of the time max t of maximum temperature achievement. After achieving the highest temperature value at the particular depths in the range mm 6 0 < ≤ z , the temperature drop lasts until standstill, whereas, , the temperature monotonically rises during the whole braking process.  The increase of heat exchange with an environment leads to linear decrease of maximum temperature T max of the friction surface of the disc (Figure 11). The highest value T max = 473 • C is achieved at t max = 3.39 s for h = 50 Wm −2 K −1 . Enhancement of convective cooling intensity to h = 250 Wm −2 K −1 results in the drop of the maximum temperature to T max = 460 • C at t max = 3.29 s.
The increase of heat exchange with an environment leads to linear decrease of maximum temperature max T of the friction surface of the disc (Figure 11).  , and it remains steady on this level within mm 30 mm 7 ≤ ≤ a . The composite's ability to conduct heat in the radial direction r K significantly increases from . On the contrary to the described above changes of effective thermal conductivity, the thermal diffusivity is the highest  Changes of effective thermophysical properties of composite material due to increasing the length a of bundles are shown in Figure 12. We see that significant rise of K r and drop of K z and k take place in the range 1 mm ≤ a ≤ 15 mm. Further increase of parameter a causes gentler variations, and for a ≥ 25 mm values of coefficients K r , K z , and k remain almost constant. Since the bundles are located on the planes oriented parallel to the front disc surface z = 0 (Figure 7), their length has slight effect on heat conduction in the axial direction. For the shortest considered length a = 1 mm, the effective coefficient of thermal conductivity in this direction has the greatest value K z = 28.3 W m −1 K −1 . Inconsiderable elongation of bundles results in the slight decrease to K z ≈ 25 W m −1 K −1 , and it remains steady on this level within 7 mm ≤ a ≤ 30 mm. The composite's ability to conduct heat in the radial direction K r significantly increases from K r = 28.7 W m −1 K −1 for a = 1 mm to K r = 63.5 W m −1 K −1 for a = 30 mm. On the contrary to the described above changes of effective thermal conductivity, the thermal diffusivity is the highest k = 1.125 × 10 −5 m 2 s −1 for the composite with the shortest bundles (a = 1 mm). With increase of the bundle length, effective thermal diffusivity drops, attaining minimum value k = 0.985 × 10 −5 m 2 s −1 for a = 30 mm. With increasing bundle length, the maximum temperature of the heated disc surface Figure 13).
As shown in Figure 12, the coefficients of thermal conductivity in radial and axial directions of composite with short fiber bundles take close values and its thermal diffusivity is the greatest. The last leads to a faster heat removal from the surface of friction and, consequently, to relatively low value of maximum temperature. Increase of fiber bundles causes decrease of thermal diffusivity and increase of maximum temperature. Growth of max T is the most rapid over the change range mm 5 1 ≤ < a , when the drop of thermal diffusivity is the largest. With increasing bundle length, the maximum temperature of the heated disc surface z = 0 is rising from the value T max = 439.4 • C for a = 1 mm to T max = 459.2 • C for a = 30 mm( Figure 13). As shown in Figure 12, the coefficients of thermal conductivity in radial and axial directions of composite with short fiber bundles take close values and its thermal diffusivity is the greatest. The last leads to a faster heat removal from the surface of friction and, consequently, to relatively low value of maximum temperature. Increase of fiber bundles causes decrease of thermal diffusivity and increase of maximum temperature. Growth of T max is the most rapid over the change range 1 < a ≤ 5 mm, when the drop of thermal diffusivity is the largest. With increasing bundle length, the maximum temperature of the heated disc surface Figure 13).
As shown in Figure 12, the coefficients of thermal conductivity in radial and axial directions of composite with short fiber bundles take close values and its thermal diffusivity is the greatest. The last leads to a faster heat removal from the surface of friction and, consequently, to relatively low value of maximum temperature. Increase of fiber bundles causes decrease of thermal diffusivity and increase of maximum temperature. Growth of max T is the most rapid over the change range mm 5 1 ≤ < a , when the drop of thermal diffusivity is the largest.  Increase of concentration of bundles in composite results in increase of its thermal diffusivity and conductivity from values k = 0.985 × 10 −5 m 2 s −1 , K r = 63.5 Wm −1 K −1 , and K z = 25 Wm −1 K −1 at V p = 0.5 to k = 3.527 × 10 −5 m 2 s −1 , K r = 155.9 Wm −1 K −1 , and K z = 88.9 Wm −1 K −1 at V b = 0.95 ( Figure 14). Such growth of material ability to dissipate heat from friction surface of the disc causes almost linear drop of maximum temperature value T max and elongation of time t max of its achievement ( Figure 15). For the smallest considered concentration of fiber bundles in composite (V b = 0.5), the maximum temperature T max = 466 • C is the highest and attained at time t max = 3.36 s, while, for the largest value V b = 0.95, we have T max = 289.2 • C and t max = 4.73 s. ( Figure   14). Such growth of material ability to dissipate heat from friction surface of the disc causes almost linear drop of maximum temperature value max T and elongation of time max t of its achievement ( Figure 15). For the smallest considered concentration of fiber bundles in composite (  ( Figure   14). Such growth of material ability to dissipate heat from friction surface of the disc causes almost linear drop of maximum temperature value max T and elongation of time max t of its achievement ( Figure 15). For the smallest considered concentration of fiber bundles in composite (

Conclusions
The calculation scheme is proposed to estimate the maximum temperature of a multi-disc brake system. For this purpose, the boundary-value problem of heat conduction is formulated for a single disc considering the frictional heating on its front surfaces, convective cooling of its lateral surfaces, and the detailed structure of the composite material. This material is made up of layers of fiber bundles arranged in planes parallel to the front surfaces of the disc. Each bundle has a cylindrical shape filled with matrix material reinforced with carbon fibers. The effective thermal conductivities of this composite are calculated at the three scale levels: micro, meso, and macro. The exact solution of the problem is obtained, which allows investigating variations of disc temperature under the effect of input parameters of braking process and friction material parameters (such as the concentration of fibers and their dimensions and orientation). Calculations were carried out for a disc made up of carbon composite material Termar-ADF, widely applied in aircraft brakes. Due to performed numerical analysis, it is established that: • Temperature values calculated on the basis of the obtained analytical solution are consistent with appropriate experimental data presented in [27].

•
During short-lasting (0 < t s ≤ 10 s), heavy loaded (p 0 1 MPa) braking processes, the influence of convective cooling of the lateral disc surfaces on its maximum temperature is inconsiderable. The decrease of the maximum temperature value with increase of coefficient of heat exchange is linear.

•
The increase of the concentration of fiber bundles in the composite results in increase of its thermal conductivity and diffusivity.

•
The high temperature of the disc at fixed intensity of convective cooling can be decreased by using shorter fiber bundles and simultaneously increasing their concentration in the composite material.
The solution was obtained without considering possible variations of friction coefficient and material properties under the influence of temperature. It should be noted that consideration of thermal sensitivity and friction coefficient changes is possible by using the nonlinear models, which require numerical methods, such as FEM, to solve the proper thermal problems of friction. Some steps are taken to do this in [29], although for a different material. Another problem that is difficult to solve, arising in this context, is related to the development of composite structure models with temperature-dependent thermophysical properties of its components. Even in the considered case of constant thermophysical properties, the procedure for determining the effective properties of the composite is quite complicated. We achieved this by using some geometric simplifications, in particular by converting the actual circular cross-sectional shape of the bundle into a square. The development of mathematical models of the composite, taking into account thermal sensitivity of its components and their interaction will be one of the directions of our future research.
Based on these findings, we conclude that the proposed mathematical model allows with sufficient accuracy determining a temperature mode of a multi-disc brake at a single braking. Calculations are executed only for one frictional composite material, but the model proposed can be used for other materials, too. In such a case, it may be necessary to use other methods for obtaining effective thermophysical properties of the composite [38]. Funding: This study was performed within the framework of project financing through the program of the Minister of Science and Higher Education of Poland named "Regional Initiative of Excellence" in 2019-2022, project number 011/RID/2018/19, amount of financing 12,000,000 PLN.

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