Nondimensional Characterization of the Operational Envelope of a Wet Friction Clutch

Abstract

In recent years, multidisc wet friction clutches are used in demanding powertrains of automatic and dual clutch transmissions targeting high efficiency and smoothness during gearshift. However, the developed flow pattern between the clutch discs is significantly complex and the Computational Fluid Dynamics (CFD) methods employed are quite demanding in terms of computational cost. To deal with this issue semi-analytical solutions were derived, which are limited, however, to specific problems, in order to obtain handy expressions, while also providing insight to the wet clutch physics. Nevertheless, this lack of global validity is counterbalanced by the fact that the governing equations become analytically solvable at specific operational conditions with satisfactory accuracy, provided that the simplifications rendering the truncated terms inactive hold true. In this work, a quantitative way of determining the relative weight of each term of the Navier-Stokes (NS) equations set is presented, based on the post-processing of CFD results using the Buckingham “π-theorem”. The sets of nondimensional numbers created were used to describe and model the physics of the wet clutch.

1. Introduction

Multidisc wet friction clutches have been used widely in automotive transmission systems since the beginning of the 21st century due to their efficiency and flexible operation regarding torque and motion transmission during gearshift. The major difference from their single-disc dry friction counterparts is the use of hydraulic fluid, which is responsible for the hydrodynamic transmission of the clutch. During hydrodynamic transmission, the fluid is subjected to a torsional Couette-type flow, therefore, applying viscous shear stresses on the surface of the driven discs and consequently enabling the gradual increase or decrease of the torque and speed to the values required by the clutch output shaft in acceleration and deceleration situations, respectively. The engagement pressure is applied onto the first from the series of discs (clutch-pack) resulting in squeezing-out action of the fluid between discs. As soon as the fluid has been squeezed out of the gap, the discs come into contact and additional shear stresses are developed due to asperity contact (mixed or dry lubrication). This procedure describes briefly the engagement phase of a wet friction clutch that takes place during gearshift [1,2]. In a similar way, when the clutch is to be disengaged, the compressive force needed for the engagement is lifted and an internal flow is developed, usually from the inner radius of the discs forcing them to assume their normal open position with the help of springs. The main design goal for both transitional states is the minimization of the clutching time (in the order of hundredths of a second) along with achieving smooth and effective torque-speed transmission. An ideal such transition is qualitatively illustrated in Figure 1.

When the clutch is inactive (i.e., no power is transmitted), it is assumed to be in a steady disengaged state. In this position the clutch acts as power loss to the transmission system, due to the drag torque developed on the discs. This matter has been investigated thoroughly [3,4,5] because the efficient calculation and understanding of drag torque can lead to its minimization and, therefore, to further improvement of the wet clutch technology. While the design objectives are different for each state, they do happen to counteract each other sometimes. An approximate expression of drag torque (T) [6] which is valid at low speeds is:

$$\mathrm{T}\propto \mathrm{rA}\left(\mathsf{\mu }\frac{\mathrm{v}}{\mathrm{h}}\right)$$

(1)

where (r) is the reference radius, (A) the surface area, (μ) the dynamic viscosity, (v) the speed, and (h) the film thickness. Since the torque is inversely proportional to (h), its increase could minimize losses, but, on the other hand, clutching time would increase.

It is clear that the efficient modeling, optimization, and design process of a multidisc wet friction clutch implies the solution of the fluid mechanics problem in the gap between two consecutive discs, within the total operational envelope of the clutch. The developed flow pattern is significantly complex, not being able to admit a fully analytical solution of the governing flow equations. Instead, this is achieved using CFD analysis and conclusions are drawn regarding the evaluation and effect of most design parameters on clutch response. The influence of friction disc groove pattern, shape, and size on drag torque as well as on the engagement torque has been investigated in [7,8,9,10], respectively, using CFD models, which were further validated by experimental results. Additionally, many works have been conducted to simulate the two-phase flow of Automatic Transmission (AT) fluid with air in order to predict aeration in disengaged clutches in regimes of high rotational speed where air is inserted between the discs [11,12,13,14]. Especially the flow in the disengaged condition resembles the well-known flow between rotating discs as initially studied by v. Karman and thus a plethora of literature reports are available, which deal with the derivation of velocity and pressure fields, boundary layers, flow instabilities, etc. Launder et al. [15] provided a review article regarding the swirling flows between cylindrical or annular cavities.

Nonetheless, the numerical treatment of the flow dynamics problem is sometimes quite demanding in terms of computational cost. Especially in simulations including the complete power transmission systems involving wet multidisc clutches, (full disc-pack, discs elasticity problem, transient simulation, etc.) the computational cost becomes excessive since the simulation should include the real-time solution of fluid dynamics problem and calculation of the necessary flow parameters. This matter was reported in [13] and researchers inevitably made use of simplified, case-specific, and usually two-disc models to investigate the problem. On top of that, the need for the extraction of generalized results, able to provide a general insight in the complex dynamics of the clutch along with the excessive theoretical background existing in the analytical and experimental study of flows between rotating discs, has concluded in the extraction of analytical or empirical expressions regarding the calculation of the necessary flow parameters. This procedure is based on the derivation of truncated forms of the Navier–Stokes equations set, which is analytically solvable [3,16].

However, as indicated in this work, those truncated forms are not valid throughout the total operational envelope of the device, as the terms neglected may lead to considerable truncation errors. In this paper, a quantitative way of determining the relative weight of each term of the NS equations set based on nondimensional analysis is presented. Groups of nondimensional parameters constructed using the Buckingham pi theorem and other parameters from the post-processing of CFD results were correlated with each other, leading to the formation of equations describing the clutch operation and, therefore, drawing limits defining the validity of regimes of the available analytical expressions.

2. Fluid Film Model

2.1. Geometry Modifications and Basic Parameters

The simplified model of the fluid film between two consecutive discs of the clutch is presented in Figure 2. The upper disc rotates with angular velocity (ω₁) (representing engine speed) and at the same is moving axially downwards due to the force applied on it by a linear actuator, representing the working principle of a wet multidisc clutch. The lower disc is rotating with angular velocity (ω₂) (representing output shaft speed) and is assumed to be fixed along the axial direction. The discs are considered rigid and flat, the flow outlet pressure is equal to atmospheric, and the inner and outer radius of discs are denoted as (r_in) and (r_out), respectively.

Regarding the flow developed in the gap, its two dominant types are the torsional Couette flow, due to the rotation of the discs, and the outflow towards the openings (P = 0) located at r_in and r_out, due to the squeezing effect occurring as the gap closes. While the first type is very well investigated in the literature as the single most important flow component between rotating discs, the superposition of the squeeze-flow on the previous one is not fully studied yet. Indicatively, in relative studies, the axial (closing) velocity is not considered and the modeling implies fluid supply to the system from an external flow source located at the inner radius [3,4,5,11,13,17]. This model takes into account the gap closing velocity, which however does not lead to an exactly equivalent formulation with respect to the literature due to the fact that the axial velocity of the top disc is a product of an actuation force variable with time and highlights many of the undesirable characteristics of wet clutches (i.e., fluctuations in position due to high fluid inertia as small fluid volume is concentrated and activated via a hydraulic actuator).

Το proceed with the study, the flow type (laminar, turbulent, or transitional) developed in various areas of the gap were examined. Defining the rotational Reynolds number as ${\mathrm{Re}}_{\mathrm{r}}=\mathsf{\omega }{\mathrm{r}}_{\mathrm{out}}^2/\mathsf{\nu }$ and the gap to radius aspect ratio as G = h/(r_out−r_in), the flow type can be identified using the chart provided by Daily and Nece [18] (for rotor-stator configurations). The values of the parameters used in the present study are concentrated in

Table 1. Parameter values.

Symbol	Description	Value	Unit
rout	Outer radius	57.5	mm
rin	Inner radius	20	mm
h	Film thickness	0.1–0.5	mm
$\stackrel{.}{\mathrm{h}}$	Gap closing velocity	0–0.5	m/s
ω1	Disc 1 rot. speed.	1000–5000	rpm
ω2	Disc 2 rot speed.	0	rpm
μ	Dynamic viscosity	0.01	Pa s
ρ	Density	800	kg/m3

.

In the present study, the aeration process was not investigated and, therefore, the gap was considered fully filled with AT fluid, independently of the rotational speed. Regarding the gap closing velocity, its maximum value of 0.5 m/s was included in the simulation range despite its high magnitude, in order to study possible overload conditions of the applied engagement force. Due to the thin film, the developed fluid pressure reached values in the order of tens of bars, and, therefore, the applied force, which had been reported to reach very high peak values in [1], may result in instantaneous dramatic increase of closing velocity. It is noted that closing velocity gets its maximum value shortly after the start of the engagement phase and it gradually decreased asymptotically to zero where discs attained their stacked condition enabling the dry friction transmission process.

For the above problem, it was calculated that G = 0.013 and 3.1 $\times$ 10⁴ ≤ ${\mathrm{Re}}_{\mathrm{r}}$ ≤ 1.3 $\times$ 10⁵ (for h = 0.5 mm). Therefore, it can be assumed that the flow was laminar, which agreed with the remarks made in [15,19].

2.2. Flow Equations

The flow was considered steady, incompressible, and isothermal, and the fluid was Newtonian, so (ρ) and (μ) are constant. The velocity vector was ${\mathrm{u}=[\mathrm{u}}_{\mathrm{r}}{\mathrm{u}}_{\mathsf{\vartheta }}{\mathrm{u}}_{\mathrm{z}}{]}^{\mathrm{T}}$ and the continuity and Navier–Stokes equations in the stationary reference frame are written in polar coordinates r, $\mathsf{\vartheta },$ and z, taking into account that ∂⁄∂ϑ = 0 due to the axisymmetry of the problem [20].

Continuity equation:

$$\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }\left({\mathrm{ru}}_{\mathrm{r}}\right)}{\mathsf{\partial }\mathrm{r}}+\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }\mathrm{z}}=0.$$

(2)

Momentum equation in radial direction:

$$\mathsf{\rho }\left({\mathrm{u}}_{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }\mathrm{r}}-\frac{{\mathrm{u}}_{\mathsf{\vartheta }}^2}{\mathrm{r}}{+\mathrm{u}}_{\mathrm{z}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }\mathrm{z}}\right)=-\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{r}}+\mathsf{\mu }\left[\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }{\mathrm{r}}^2}+\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }\mathrm{r}}+\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }{\mathrm{z}}^2}-\frac{{\mathrm{u}}_{\mathrm{r}}}{{\mathrm{r}}^2}\right].$$

(3)

Momentum equation in angular direction:

$$\mathsf{\rho }\left({\mathrm{u}}_{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }\mathrm{r}}+\frac{{\mathrm{u}}_{\mathrm{r}}{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathrm{r}}{+\mathrm{u}}_{\mathrm{z}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }\mathrm{z}}\right)=\mathsf{\mu }\left[\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }{\mathrm{r}}^2}+\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }\mathrm{r}}-\frac{{\mathrm{u}}_{\mathsf{\vartheta }}}{{\mathrm{r}}^2}+\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }{\mathrm{z}}^2}\right].$$

(4)

Momentum equation in axial direction:

$$\mathsf{\rho }\left({\mathrm{u}}_{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }\mathrm{r}}{+\mathrm{u}}_{\mathrm{z}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }\mathrm{z}}\right)=-\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{z}}+\mathsf{\mu }\left[\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }{\mathrm{r}}^2}+\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }\mathrm{r}}+\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{z}}}{\mathsf{\partial }{\mathrm{z}}^2}\right].$$

(5)

The boundary conditions for the Equations (2)–(5) are:

$$\mathrm{P}\left({\mathrm{r}}_{\mathrm{in}}\right)=\mathrm{P}\left({\mathrm{r}}_{\mathrm{out}}\right)=0$$

(6)

$${\mathrm{u}}_{\mathrm{r}}\left(\mathrm{r},\mathrm{h}\right){=\stackrel{.}{\mathrm{h}}\mathrm{and}\mathrm{u}}_{\mathrm{r}}(\mathrm{r},0)=0$$

(7)

$${\mathrm{u}}_{\mathrm{r}}\left(\mathrm{r},\mathrm{h}\right){=0\mathrm{and}\mathrm{u}}_{\mathrm{r}}(\mathrm{r},0)=0$$

(8)

$${\mathrm{u}}_{\mathsf{\vartheta }}\left(\mathrm{r},\mathrm{h}\right){=\mathsf{\omega }}_1{\mathrm{r}\mathrm{and}\mathrm{u}}_{\mathsf{\vartheta }}\left(\mathrm{r},\mathrm{h}\right){=\mathsf{\omega }}_2\mathrm{r}.$$

(9)

The nondimensionalization of the flow equations took place, introducing the following nondimensional terms.

$$\tilde{\mathrm{z}}=\frac{\mathrm{z}}{\mathrm{h}},\tilde{\mathrm{r}}=\frac{\mathrm{r}}{{\mathrm{r}}_{\mathrm{out}}},{\tilde{\mathrm{u}}}_{\mathrm{r}}=\frac{{\mathrm{u}}_{\mathrm{r}}}{{\mathrm{V}}_{\mathrm{r}}},{\tilde{\mathrm{u}}}_{\mathsf{\vartheta }}=\frac{{\mathrm{u}}_{\mathsf{\vartheta }}}{{\mathrm{V}}_{\mathsf{\vartheta }}},{\tilde{\mathrm{u}}}_{\mathrm{z}}=\frac{{\mathrm{u}}_{\mathrm{z}}}{{\mathrm{V}}_{\mathrm{z}}},\tilde{\mathrm{P}}=\frac{\mathrm{P}}{\mathsf{\rho }{\mathrm{V}}_{\mathrm{r}}^2}$$

(10)

where

$${\mathrm{V}}_{\mathrm{r}}=\frac{{\mathrm{r}}_{\mathrm{out}}\stackrel{.}{\mathrm{h}}}{\mathrm{h}}{,\mathrm{V}}_{\vartheta }{=\mathrm{r}}_{\mathrm{out}}{\mathsf{\omega }}_1{,\mathrm{V}}_{\mathrm{z}}{=\mathrm{G}\mathrm{V}}_{\mathrm{r}}.$$

As previously mentioned, there do exist analytical expressions for the pressure field and velocity distributions within the fluid film, and they were derived from truncated forms of the Navier–Stokes equations. To further facilitate nondimensionalization, only one velocity component, namely V_ϑ, was used for all velocity components. However, this over-simplification led to non-acceptable expressions, such us in the case of (u_ϑ), where, if all terms except (∂²u_ϑ⁄∂z²) were neglected, it followed a linear distribution with respect to (z) [3], which in the general case is erroneous, as illustrated in Figure 3. Referring to Figure 3, the order of magnitude of the momentum equations terms was presented for two working points of the clutch, corresponding to open and almost closed gap. It was evident that the variation of their values was significant and, thus, no global conclusions may be extracted concerning their truncation, except for the case of (∂P⁄∂z) in Equation (5) and the second derivatives with respect to (z) in the second part of Equations (2) and (3), respectively, which are dominant.

3. Analytical Solutions and CFD Results

3.1. Analytical Solutions

The majority of the available analytical expressions for pressure and velocity distributions within the film were derived mainly for flows between rotating discs and they do not take into account the squeezing effect that exists in a wet clutch. For the needs of this stud, they were modified accordingly to include the axial closing velocity of the gap $\stackrel{.}{\mathrm{h}}$ instead of an inflow to the system. The expression of radial velocity, u_r(r,z), which contained the pressure gradient (dP⁄dr), was used in Equation (2), which was integrated with respect to (z) under the boundary conditions of Equation (7). Therefore, the second order differential equation for pressure was obtained, which was integrated with respect to (r) under the boundary conditions of Equation (6) and the final expression for pressure was derived. The two cases considered corresponded to the following system of truncated Navier–Stokes equations [3,16].

$$-\frac{{\mathrm{u}}_{\mathsf{\vartheta }}^2}{\mathrm{r}}=-\frac{1}{\mathsf{\rho }}\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{r}}+\frac{\mathsf{\mu }}{\mathsf{\rho }}\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }{\mathrm{z}}^2},0=\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }{\mathrm{z}}^2},0=\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{z}}$$

(11)

$$-\frac{{\mathrm{u}}_{\mathsf{\vartheta }}^2}{\mathrm{r}}=-\frac{1}{\mathsf{\rho }}\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{r}}+\frac{\mathsf{\mu }}{\mathsf{\rho }}\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathrm{r}}}{\mathsf{\partial }{\mathrm{z}}^2}{,\mathrm{u}}_{\mathrm{r}}\frac{\mathsf{\partial }{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }\mathrm{r}}+\frac{{\mathrm{u}}_{\mathrm{r}}{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathrm{r}}=\frac{\mathsf{\mu }}{\mathsf{\rho }}\frac{{\mathsf{\partial }}^2{\mathrm{u}}_{\mathsf{\vartheta }}}{\mathsf{\partial }{\mathrm{z}}^2},0=\frac{\mathsf{\partial }\mathrm{P}}{\mathsf{\partial }\mathrm{z}}$$

(12)

The major differences of Equations (11) and (12) are the terms in the first part of the momentum equation in the angular direction, which differentiate the variation of u_ϑ(r,z) with respect to (z), being a first and fifth order polynomial, respectively, according to [3,16]. Since the researchers in [3,16] assumed that the raised pressure within the fluid film was associated with the radial velocity profile or equivalently with the inflow to the system (Q), they both yield the same Ordinary Differential Equation (ODE) for the pressure:

$$\frac{\mathrm{dP}}{\mathrm{dr}}=\frac{-6\mathrm{Q}\mathsf{\mu }}{{\mathsf{\pi }\mathrm{rh}}^3}+\mathsf{\rho }\mathrm{r}\left({\mathsf{\omega }}_2^2{+\mathsf{\omega }}_2\mathsf{\Delta }\mathsf{\omega }+\frac{3}{10}{\mathsf{\Delta }\mathsf{\omega }}^2\right).$$

(13)

In this work, instead of the inflow (Q), the closing velocity ($\stackrel{.}{\mathrm{h}}$) was considered and, therefore, the pressured ODE had the form:

$$\frac{{\mathrm{d}}^2\mathrm{P}}{{\mathrm{dr}}^2}+\frac{1}{\mathrm{r}}\frac{\mathrm{dP}}{\mathrm{dr}}=\frac{12\mathsf{\mu }\stackrel{.}{\mathrm{h}}}{{\mathrm{h}}^3}+2\mathsf{\rho }\left({\mathsf{\omega }}_2^2{+\mathsf{\omega }}_2\mathsf{\Delta }\mathsf{\omega }+\frac{3}{10}{\mathsf{\Delta }\mathsf{\omega }}^2\right).$$

(14)

where Δω = ω₁−ω₂. Equation (14) was solved under the boundary conditions of Equation (6), yielding the analytical pressure distribution:

$$\mathrm{P}\left(\mathrm{r}\right)=\frac{1}{4}\left[12\frac{\mathsf{\mu }\stackrel{.}{\mathrm{h}}}{{\mathrm{h}}^3}+2\mathsf{\rho }\left({\mathsf{\omega }}_2^2{+\mathsf{\Delta }\mathsf{\omega }\mathsf{\omega }}_2+\frac{3}{10}{\mathsf{\Delta }\mathsf{\omega }}^2\right)\right]\left[{(\mathrm{r}}^2{-\mathrm{r}}_{\mathrm{out}}^2)+\left({\mathrm{r}}_{\mathrm{out}}^2{-\mathrm{r}}_{\mathrm{in}}^2\right)\frac{\ln \left(\frac{{\mathrm{r}}_{\mathrm{out}}}{\mathrm{r}}\right)}{\ln \left(\frac{{\mathrm{r}}_{\mathrm{out}}}{{\mathrm{r}}_{\mathrm{in}}}\right)}\right].$$

(15)

3.2. Computational Model

The validity of the pressure distribution corresponding to Equation (15) was investigated through a comparative study with the results from a CFD model where the full set of equations was used for solving the problem. The geometry of the model along with the boundary conditions (Equations (6)–(9)) are presented in Figure 4. The CFD environment used was ANSYS Fluent and the model was treated as 3D. However, in order to save computational time, the problem was modeled as 2D axisymmetric including swirl. For both approaches a mesh independency study was conducted (Figure 5) and the results were compared and found to be in good agreement.

The tangential momentum equation for 2D swirling flows may be written as [21]:

$$\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }}{\mathsf{\partial }\mathrm{z}}\left({\mathrm{r}\mathsf{\rho }\mathrm{u}}_{\mathrm{z}}\mathrm{w}\right)+\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }}{\mathsf{\partial }\mathrm{r}}\left({\mathrm{r}\mathsf{\rho }\mathrm{u}}_{\mathrm{r}}\mathrm{w}\right)=\frac{1}{\mathrm{r}}\frac{\mathsf{\partial }}{\mathsf{\partial }\mathrm{z}}\left[\mathrm{r}\mathsf{\mu }\frac{\mathsf{\partial }\mathrm{w}}{\mathsf{\partial }\mathrm{x}}\right]+\frac{1}{{\mathrm{r}}^2}\frac{\mathsf{\partial }}{\mathsf{\partial }\mathrm{r}}\left[{\mathrm{r}}^3\mathsf{\mu }\frac{\mathsf{\partial }}{\mathsf{\partial }\mathrm{r}}\left(\frac{\mathrm{w}}{\mathrm{r}}\right)\right]-\mathsf{\rho }\frac{{\mathrm{u}}_{\mathrm{r}}\mathrm{w}}{\mathrm{r}}$$

(16)

where (w) is the swirl velocity. The flow was considered steady for each gap closing velocity, $\stackrel{.}{\mathrm{h}}$, laminar, and the sum of residuals was set below 10⁻⁶ for all equations. For the needs of this study, Disc 2 was assumed fixed. SIMPLE algorithm was used for the couple of pressure and velocity, and second order upwind spatial discretization scheme was used for the momentum conservation and pressure equations. The results of the pressure distribution versus radius are presented in Figure 6 for six operational points of the clutch. Examining Figure 6, the deviation of the pressure analytical expressions to the CFD results is observed, especially in the case where the film thickness had the value of the normal open gap (0.5 mm). This happens because the terms neglected in Equations (11) and (12) hold a considerable influence in these operating conditions. As the film thickness decreases (0.1 mm), the agreement of the analytical results to the CFD-obtained ones improved, as expected because terms including first- and second-order derivatives with respect to (z) became dominant.

In Figure 7 and Figure 8, pressure contours and streamlines together with velocity vector fields are presented for four of the above six cases. The pressure distribution was parabolic, which is also depicted in Figure 6, and its maximum value was in a zone between r_in and r_out, where flow lines were less dense, since the radial component of velocity was zero and, therefore, a zero-flow zone was developed at this region. Regarding the influence of rotational speed, it was concluded that when squeezing motion was significant, which was characterized by the term ($\stackrel{.}{\mathrm{h}}$/h³) [3], then variation of rotational speed did not affect significantly developed pressure. In contrast, as film thickness increased, the correlation between Couette and squeeze mechanisms changed with increasing of rotational speed while fluid pressure tended to decrease.

4. Nondimensional Analysis and Post-Process of CFD Results

4.1. Creation of Nondimensional Numbers

The purpose of this study was the characterization of the behavior observed in Figure 6 regarding the pressure variation with respect to variables h, $\stackrel{.}{\mathrm{h}}$, and ω₁. This was achieved through the creation of nondimensional parameters using the Vaschy–Buckingham “π-theorem“ [22]. According to that theorem, if a physical process is described by a group of N natural quantities, it is possible to reduce the number of independent parameters to N-M nondimensional quantities, where M is the total number of the basic dimensions (length (L), time (T), mass (M), absolute temperature (Θ), etc.) that appear in the process. Regarding the clutch operational envelope, the most important parameters describing the pressure variations are eight: Density ρ (M)(L)⁻³, dynamic viscosity μ (M)(L)⁻¹(T)¹, film thickness h (L), closing velocity $\stackrel{.}{\mathrm{h}}$ (L)(T)⁻¹, top disc annular velocity ω₁ (T)⁻¹, inner and outer radius r_in, r_out (L), and the integral of pressure on the surface of the discs F (M) (L) (T)⁻². If (ρ), (ω₁), and (h) are selected as the repeating variables, then according to the “π-theorem”, five nondimensional parameters may be created:

$$\begin{gathered} {\mathsf{\Pi }}_1{=\mathsf{\mu }\mathrm{h}}^{{\mathrm{a}}_1}{\mathsf{\omega }}^{{\mathrm{b}}_1}{\mathsf{\rho }}^{{\mathrm{c}}_1} \\ {\mathsf{\Pi }}_2{=\stackrel{.}{\mathrm{h}}\mathrm{h}}^{{\mathrm{a}}_2}{\mathsf{\omega }}^{{\mathrm{b}}_2}{\mathsf{\rho }}^{{\mathrm{c}}_2} \\ {\mathsf{\Pi }}_3{=\mathrm{r}}_{\mathrm{in}}{\mathrm{h}}^{{\mathrm{a}}_3}{\mathsf{\omega }}^{{\mathrm{b}}_3}{\mathsf{\rho }}^{{\mathrm{c}}_3} \\ {\mathsf{\Pi }}_4{=\mathrm{r}}_{\mathrm{out}}{\mathrm{h}}^4{\mathsf{\omega }}^{{\mathrm{b}}_4}{\mathsf{\rho }}^{{\mathrm{c}}_4} \\ {\mathsf{\Pi }}_5{=\mathrm{F}\mathrm{h}}^5{\mathsf{\omega }}^{{\mathrm{b}}_5}{\mathsf{\rho }}^{{\mathrm{c}}_5}. \end{gathered}$$

(17)

Finally, the following five nondimensional numbers describing the operation of the clutch are obtained:

$$\begin{gathered} {\mathsf{\Pi }}_1=\frac{\mathsf{\mu }}{{\mathrm{h}}^2{\mathsf{\omega }}_1^2\mathsf{\rho }} \\ {\mathsf{\Pi }}_2=\frac{\stackrel{.}{\mathrm{h}}}{\mathrm{h}{\mathsf{\omega }}_1} \\ {\mathsf{\Pi }}_3=\frac{{\mathrm{r}}_{\mathrm{in}}}{\mathrm{h}}{,\mathsf{\Pi }}_4=\frac{{\mathrm{r}}_{\mathrm{out}}}{\mathrm{h}} \\ {\mathsf{\Pi }}_5=\frac{\mathrm{F}}{{\mathrm{h}}^4{\mathsf{\omega }}_1^2\mathsf{\rho }}. \end{gathered}$$

(18)

4.2. Post-Processing of CFD Results

In order to describe the variation of pressure with respect to h, $\stackrel{.}{\mathrm{h}}$, and ω₁, it is necessary to draw a relation between Π₂ and Π₅ in the form of Π₅ = f (Π₂). To achieve that, a multiparameter CFD analysis of the variation of the pressure integral with respect to above variables was conducted for a limited number of cases (five different values for each variable leading to the study of 125 separate cases in total, as shown in

Table 2. Values of h, $\stackrel{.}{\mathrm{h}}$, ω₁ used for calculation of F.

h [mm]	$\stackrel{.}{\mathrm{h}} [\mathbf{m}/\mathbf{s}]$	ω1 [rpm]
0.5	0.5	1000
0.4	0.375	2000
0.3	0.25	3000
0.2	0.125	4000
0.1	0	5000

.

For a better understanding of the range of (Π₂), Figure 9 illustrates its magnitude, calculated for the parameter values given in Table 2.

In Figure 10 and Figure 11, the comparison of the value of (Π₅) obtained from the CFD results and the calculated (Π₅) from the analytical pressure distribution (Equation (15)) is presented with respect to (Π₂). For the needs of illustration, the results are given for a specific value of (ω = 3000 rpm), while (h) and ($\stackrel{.}{\mathrm{h}}$) were used as independent varying parameters in the plots. Selecting a specific value of (ω) was justified from the small influence it had on the developed pressure compared to (h) and ($\stackrel{.}{\mathrm{h}}$), as depicted in Figure 7 and Figure 8. Due to increased magnitude of (Π₅), Log₁₀ scale was used for its axis and, therefore, negative values were omitted. Furthermore, two sets of plots are presented in Figure 10 and Figure 11, where one of (h), ($\stackrel{.}{\mathrm{h}}$) was varying and the other one was assumed fixed, having the values of Table 2. The case where (h) was varied and ($\stackrel{.}{\mathrm{h}}$) was assumed constant (Figure 10), was associated with the utilization of a constant closing velocity which is maintained by the clutch actuator during the engagement process. Instead, the case where ($\stackrel{.}{\mathrm{h}}$) was varied and (h) was assumed constant (Figure 11) corresponded to a constant engagement force.

Investigating Figure 11, it can be concluded that there was a continuous deviation between the analytical expression and CFD results for the developed pressure within the fluid film over the total regime of (h) and, therefore, it is proposed that CFD analysis should be used in the case where a constant force model is simulated. On the other hand, in the case depicted in Figure 10, a constant limit of (Π₂) under which the CFD and analytical results agreed may be extracted at Π₂ ≤ 4, combined with closing velocity less than $\stackrel{.}{\mathrm{h}}$ ≤ 0.125 m/s (Figure 10-Bottom right). This trend can be used by researchers to optimize the use of pressure analytical expressions in order to calculate the necessary actuator force in wet clutch models and save computational cost.

5. Conclusions

In this study a quantitative analysis was performed in order to evaluate the validity of the currently available analytical solutions of the governing flow equations (obtained from truncated sets of Navier–Stokes equations) to prescribe accurately the physics of a wet friction clutch throughout most of its operational envelope. The investigation was focused on the deviation of the pressure distribution within the fluid film and a comparative study was performed between the results of the analytical expressions and the CFD analysis. Finally, via the creation of newly introduced nondimensional numbers, the characterization of the clutch operation was achieved.