Seal surface topography typically consists of global-scale geometric features as well as local-scale roughness details and homogenization-based approaches are, therefore, readily applied. These provide for resolving the global scale (large domain) with a relatively coarse mesh, while resolving the local scale (small domain) in high detail. As the total flow decreases, however, the flow pattern becomes tortuous and this requires a larger local-scale domain to obtain a converged solution. Therefore, a classical homogenization-based approach might not be feasible for simulation of very small flows. In order to study small flows, a model allowing feasibly-sized local domains, for really small flow rates, is developed. Realization was made possible by coupling the two scales with a stochastic element. Results from numerical experiments, show that the present model is in better agreement with the direct deterministic one than the conventional homogenization type of model, both quantitatively in terms of flow rate and qualitatively in reflecting the flow pattern.
Pressure-driven flow through preamble structures such as the percolation of fluids in static seals, see e.g. [1–3] and the flow through fractured porous media, in e.g. [4,5], is a problem that has been addressed in many works before. The numerical solution of this problem is, however, complicated due to its multi-scale nature. Indeed, the flow domain is usually of the order of millimetres or even larger while the level of detail necessary to resolve the narrowest, still influential, constrictions may require nanometre order resolution . As large domains resolved with a high level of detail are often required, the multi-scale nature of the flow problem prevents the usage of deterministic solutions.
Among other solutions, the problem associated with the flow through the gap between two rough surfaces this problem has repeatedly been commonly addressed in the lubrication field by separating the problem into two scales, i.e. a global scale covering the full seal domain and a local scale accounting for the roughness details. This type of scale separation was used by Christensen , who presented a stochastic two-scale approach, Elrod , who applied a perturbation technique and by Patir & Cheng  who presented a statistical approach and derived so-called flow factors, to model the roughness effect on the fluid flow between two rough surfaces and to link the two scales together. This concept has been further developed by others (see [10–13], for instance). Of particular interest is the homogenization technique, which provides a solid mathematical foundation to this approach (e.g. ). An application to the mixed lubrication problem, including a study on leakage, can be found in the two-part paper by Sahlin et al. [15,16]. Other, more complex cases include , where cavitation is included and , where elastohydrodynamic lubrication is considered. In the afore-mentioned mixed lubrication problem, the contact mechanics model has not been treated as a two-scale approach, instead, a periodic roughness description has been adopted . This is a common approach. An exception to that is , where the elastohydrodynamic lubrication problem is formulated by means of separation into two scales. Two-scale approaches have turned out to be applicable in many applications, however, as the gap becomes thinner the flow becomes smaller and the local-scale model requires larger and larger domains to produce a converged value for the flow factors. This makes the two-scale approach lose its effectiveness.
The percolation problem has also been a common issue between those studying the flow through porous media. This problem is essentially equal to the one describing the pressure-driven flow through the gap between two stationary rough surfaces. In this field, it is not uncommon to find stochastic representations of permeability or the porous media itself (see, e.g. [20–24]). Generally, however, the probability distribution of the permeability is an input to the model and is not computed from actual pore-scale measurements. A model for the flow through a porous media that is similar to the model presented here, can be found in . In that work the flow through a fabric pattern is investigated and the permeability distribution is generated stochastically by considering the results of computational fluid dynamics analysis over different fabric configurations.
The stochastic approach applied in porous media flow modelling is also seen as beneficial in the study of small flows between rough surfaces. A reason for this is given in the work by Dapp & Müser , where it is shown that the local-scale pressure drop, at very low flow rates, occurs over one very small constriction only. This implies that the local flow will depend on the geometry of such constrictions and that the total flow can only be described in a statistical manner because of the resolution requirement.
The idea behind this work is to describe the stochastic element by means of a two-scale formulation similar to those presented previously. This is done using the framework of heterogeneous multiscale method (HMM) . The main novelty of this work is that it permits the estimation of the uncertainty of the results due to the random nature of the topography. Another advantage with the present model is that it is possible to restrict the size of the local domain and still obtain a converged solution. Moreover, it is shown that this approach predicts a more realistic flow pattern compared with the flow patterns obtained by using conventional two-scale models for similar local domain size.
The problem to be solved is the one that governs the flow through two surfaces which are compressed against each other. As said previously, the multi-scale nature of the problem prevents the usage of a deterministic solution to such problem. Therefore, a two-scale approach is used.
In order to compute the flow rate and the flow pattern, it is first necessary to compute the deformed shape of the gap between the two surfaces (contact mechanics problem). The deformed gap is due to fluid pressure and asperity contact. In this case, we consider the fluid pressure contribution to be small and it is therefore neglected. This permits separating the problem into two smaller problems that can be solved sequentially: compute the deformed gap and compute the flow rate through that gap.
In the following, the method used in this work is presented. First, an introduction to the HMM framework, used to develop the model, is given. After that the two-scale stochastic model is introduced. This includes the theory behind the two-scale formulation of the flow model and the contact mechanics problem and the introduction of the stochastic element. Also, an overview of the solution procedure is given. Finally, the two-scale flow formulation presented in this work is compared against the well-established homogenization technique.
(a) Heterogeneous multiscale method
The HMM is a general framework used to build two-scale models  which is flexible in the sense that global and local-scale models do not need to be of the same nature. This will allow us to introduce the stochastic element in a more natural way than the more rigid multi-scale homogenization technique, this is, without the imposition of periodicity in roughness. Under the HMM framework, the model is constructed into steps: (i) definition of a global-scale model, which will include some parameter or variable containing the information of the local-scale model and (ii) definition of a local-scale model that permits finding those variables or parameters. The first step is crucial as a wrong definition of the global-scale method can lead to a wrong specification of the local-scale constraints. Although it can be derived from the local-scale data, it is better to use previous (experimental or theoretical) knowledge about the global-scale behaviour. The local-scale model is usually easier to identify, as one can use a more fundamental representation of the reality. However, one must be careful on the selection of boundary conditions for the local-scale problem and the data processing used to transfer information to the global scale. In particular, the boundary conditions in the local-scale model must be specified so that the solution is consistent with the global-scale model. One should understand by consistent that the local scale should be constructed so that the global scale can be seen as a coarse representation of it.
(b) The two-scale stochastic model
In this section, the two-scale stochastic model developed is presented. We start by introducing the two-scale methods for both the flow model and the contact mechanics problem. After this, the stochastic element is introduced in the two-scale formulation. Finally, the solution procedure is outlined.
(i) The two-scale flow model
Before presenting the two-scale formulation, we introduce the deterministic flow model. In this work, we consider the rectangular domain 2.1covering the full seal area. In this domain, the gap between the two rough surfaces, hϵ(x), is resolved at a high level of detail. The subscript ϵ indicates that the gap oscillates rapidly due to the roughness, thus imposing a high resolution in the representation of the domain. In order to solve the fluid pressure distribution through the gap, the lubrication approximation is used. Therefore, the well-known Reynolds equation can be used. For an iso-viscous incompressible fluid and with no relative motion between the surfaces, this equation can be written as 2.2where pϵ is the fluid pressure, which will oscillate rapidly in accordance with hϵ. The problem is completed by posing the boundary conditions. In order to obtain a pressure-driven flow, a pressure drop, Δp, is enforced in the x1-direction by setting Dirichlet boundary conditions. In the x2-direction, periodic boundary conditions are imposed. Once the fluid pressure distribution is known, the total leak rate can be computed by integrating the mass flow over the domain Ωs. In accordance with the Reynolds equation, this is, 2.3were η is the viscosity of the fluid.
As stated previously, it often becomes too computationally expensive to resolve the roughness at the global scale. Provided that the period of the oscillations is sufficiently small, one can separate the problem into two scales: one considering the detail of the oscillations (local scale) and the other accounting for the slow variations in the whole domain (global scale). Following the HMM framework, we start by presenting the global-scale model. Then we will present the local-scale one.
The global scale is solved in a domain Ωc, which is a coarse representation of Ωs. As the Reynolds equation is a mass-conservation law, it is natural to define the global-scale model also as a mass conservation law. Following a finite volume formulation, this can be posed as 2.4where j0 are the global-scale fluxes and (i,j) identify any point of the global scale. For each grid point, the fluxes must, therefore, be estimated from the solutions of the local-scale model. We will show that, in agreement with Darcy’s law, the flux is proportional to the local pressure drop, δp, and we will refer to these constants as permeability, K. Therefore, (2.4) can be rewritten as 2.5where Here the first superscript indicates the direction of the corresponding flux and the second one the direction of the local pressure drop. Note that, in order to express the flux in m3 s−1 and the permeability in m3, (2.5) is scaled by the (constant) viscosity. Together with the boundary conditions, posed equally to those in the deterministic problem, (2.5) permits computing the global-scale fluid pressure. The total leak rate is then 2.6where ∂oΩc is the part of the boundary of Ωc regarded as an outlet.
In the global-scale model, the permeabilities are the parameters containing the local-scale information. We, therefore, define the local scale in order to compute them. Let us start by deriving the fluxes in (2.4). As an example, consider ; the others can be defined analogously. In order to do so, we define a local domain between the neighbour global-scale points, 2.7which has a size Δx1×Δx2. Similarly as in the deterministic formulation, we solve the fluid pressure distribution by means of the Reynolds equation (2.2). We, however, need to pose the boundary conditions in a way that the local-scale model is consistent with the global-scale one. In order to see which conditions must be fulfilled, we assume that the global-scale component of pϵ, p0, varies linearly between the global-scale points, i.e. 2.8where the capital letters refer to the discrete representation in the global scale. In order for it to be consistent, the local-scale solution should be in agreement with the above representation. Therefore, it should satisfy 2.9where and 〈⋅〉 indicates average over the local domain. Other conditions over the mean value of pϵ could be introduced, but the flux would not be affected. Several boundary conditions can be posed so that the constraint (2.9) is satisfied. Among them, the one providing better results is the near periodic condition, i.e. periodicity of pϵ−p0 .
The boundary condition selected, however, introduces an unwanted coupling between the global and the local scale. In order to remove this coupling, we start by defining , with , near periodic in the x1-direction and periodic in the x2-direction, and , defined analogously. Furthermore, we define the non-dimensional variables Therefore, (2.2) can be written in ω as 2.10We can now find a solution by setting the two terms equal to zero independently and thus obtain two problems, both of which are independent from the global scale, i.e. 2.11aand 2.11b
It can be seen that the solution obtained by combining this two problems satisfies the condition (2.9). Once the local problems are solved, we can compute the local flux, i.e. 2.12We can now identify the local- and the global-scale contribution to the flux. By defining the permeability, K as the local-scale contribution, we can write 2.13which leads to the global-scale model (2.5).
(ii) The two-scale contact mechanics model
The fluid flows through the gaps, hϵ, between the contacting surfaces. This gap depends on the applied load. The flow defined previously is solved on a clearance between the surfaces, hϵ, which varies depending on the applied load. Therefore, the deformed shape of the gap must be computed before the flow problem can be assessed. In order to do so, a model based on the one presented in , is used. The equations governing this method can be posed as 2.14a 2.14b 2.14c 2.14dwhere the clearance, hϵ, is defined as the original (no pressure) gap, h0, plus the deformation, u, caused due to the contact pressure, pc, and a rigid body movement, g00. The deformation, u, is computed by the convolution of a coefficients kernel, k, and the contact pressure, pc. This can be seen as adding the deformation caused by the pressure at all points. The equivalent Young modulus, E*, is defined as 2.15where Ei and νi are the Young modulus and Poisson ratio of each surface. Equation (2.14c) establishes a complementarity relation between hϵ and pc, meaning that when there is contact, the clearance is zero and when there is no contact, the pressure is zero. Finally (2.14d) imposes that the average contact pressure is equal to a given total load, W. A perfect plasticity condition is defined by imposing pc≤H, where H is the hardness of the softer material. This problem can be solved by using the algorithm presented in , modified to account for the non-periodicity in the x1-direction (e.g. ).
In the same way as in the fluid computation, the computational time for the contact mechanics problem rises as the domain size grows. Therefore, a two-scale solution is also desired. In order to justify the two-scale formulation of this problem, we need to justify that (i) it is a correct approximation to use a coarser representation of the surface to capture the global-scale trends in the contact mechanics results and (ii) that it is a correct approximation to compute the contact mechanics problem in the local cells by using the nominal pressure pnom=Wlocal as the only information from the global scale.
We start by justifying the point (i). The goal is to show that the utilization of a coarse representation of h0 gives as a result a coarse representation of hϵ and pc. To do so, we apply a Gaussian filter, g, to (2.14a)–(2.14d). We define 2.16as the Fourier transform of a function f(x). The equations then read 2.17a 2.17b 2.17cwhere the capital letters identify the filtered functions. In order to obtain (2.17c), the property of the filter to maintain the mean value has been used. Equation (2.14c) has not been added because it does not hold exactly for the coarsened problem. If one think on the complementarity condition and how the filter works, one can see that it holds approximately provided that the coarsening does not alter the geometry at the global scale. Regarding the plastic limit, the pressure is expected to be lower, and in a more spread area due to the filtering. Therefore, there should be less plastic deformation. This, however, should not significantly affect the result.
Unlike the flow case, where each local cell corresponds to one point in the global scale, this is not the case in this coarse representation. In order to obtain the nominal pressure, pnom, of a local cell, the global-scale result should be averaged over a region of equivalent size.
We now discuss about the correctness of the second statement (ii). It has been shown that the elastic deformation caused by pressure is a long-range effect and it cannot be neglected . However, the requirement here is smaller. The deformation on the local domain needs to be computed only up to a constant (g00). Therefore, by assuming that the elastic deformation caused by pressures far-away from the local domain has a constant value over the considered domain, the local scale contact mechanics problem can be solved independently from the global domain.
In order to assess the validity of this assumption one should note that, according to the St. Vennant’s principle, the contribution of the distant regions is smooth. Then, by assuming that the local cell is reasonably away from the boundaries, similar deformation can be assumed in all directions, leading to a flat overall deformation. Although crude, it is an adequate first approximation. The assumption of flat deformation is clearly not adequate if one refers to the deformation on points near the cell boundaries caused by pressures in points neighbouring the local cell. In such close points, however, it is a reasonable claim that the roughness is similar to that on the cell. Therefore, a periodic roughness assumption is also reasonable.
(iii) The stochastic approach
In the previous sections, a two-scale model has been developed. The stochastic element has not, however, been introduced yet. Indeed, following the presented procedure, four cell problems (two for each direction) must be solved for each grid node in the global scale. Moreover, the results can be sensitive to the measurement used for the computations. Because the surface topography is a random process, two different measurements of the same surface (or two equivalent ones) are expected to give different results. It is, therefore, of greater interest to obtain the results in the form of an expected value plus a confidence interval rather than simply a given value.
If one now takes a look at the permeabilities computed one notices that they can be modelled as a random variable following a log-normal distribution. This distribution is, however, rather broad. Moreover, it is expected to see some spatial correlation between the different permeabilities, caused by fluctuations in the global scale. Therefore, a reference parameter, either average interfacial separation, , or nominal pressure, pnom, is taken from the global-scale contact mechanics computation and the permeabilities are fitted to different log-normal distributions as a function of these parameters. These distributions have a narrower permeability range and, therefore, are more meaningful.
Following this approach, the global scale is therefore not constructed by solving the local problem on each coarse grid point but by randomly generating a permeability value for each of those points. In order to ensure sufficient accuracy, numerous realizations of the global scale are computed following a Monte Carlo approach. The output is, therefore, given as a probability distribution instead of a single value.
It is worth noticing that modelling a gap between two rough surfaces by randomly assign permeabilities from a log-normal distribution (or other similar distributions) is a common practice in the porous media literature (e.g. [20,31]). With the present approach, the parameters for such distribution are no longer obtained experimentally but from modelling the problem at the local scale.
(iv) Solution procedure
A flow chart of the solution procedure for the model is depicted in figure 1. In the following, a more detailed description is given. Both scales are effectively separated during the computations and therefore they are treated separately here also. For clarity, the different domains considered are represented in figure 2
We start by considering the local-scale solution procedure. In order to solve the local problems, the whole measured domain is divided into a set of local cells of size Δx1×Δx2. For the following numerical analysis, these cells are mirrored in order to obtain a periodic roughness. We note that this implies that the permeabilities K12 and K21 will be zero (see  for a justification). However, even without mirroring, these are expected to be small in most cases, also in the surface textures used in this work (figure 3; this should be reconsidered the model is applied to surfaces with different textures). The solution procedure is then as follows:
(i) The deformed gap in the local domain is computed for a range of nominal pressures pnom, as described in §2b(ii). It is advised to compute the deformation for a very small nominal pressure, which shall serve as a zero value during the interpolation.
(ii) The fluid pressure distribution is computed by (2.11) and the permeability is computed according to (2.12) and (2.13). This is done for the range of nominal pressures previously computed and for a specified range of average interfacial separations . A combination of Dirichlet boundary condition in the direction of the pressure drop and homogeneous Neumann boundary conditions in the transverse direction is used instead of the near periodic condition. The reason for this is that this combination is easier to implement and yet equivalent (see  for details).
(iii) The distribution of permeabilities for the given set of cells, and for each reference parameter (pnom or ) are fitted to a log-normal distribution. The set of cells should be large enough in order to obtain a good estimate for the parameters of the log-normal distribution. A random selection of cells is desired to increase robustness against any possible spacial variation of permeability.
Once the permeabilities are obtained from the local-scale, the global-scale model is implemented in the following steps:
(i) A measurement of the topography is used as the input. This measurement domain, Ω, is filtered by using a Gaussian filter and re-sampled to a coarser grid.
(ii) The contact mechanics problem is solved using the coarse representation of the measurement as input.
(iii) The resulting pressure and interfacial separation distribution is further coarsened by averaging over domains of the same size of the local-scale cells. This results in a nominal pressure and an average interfacial separation distribution. It is important to select the same size as of the local cells as the permeability distributions might depend on the size.
(iv) The global domain is covered by repetition of the nominal pressure and average interfacial separation distributions. A global-scale variation might be included in the form of a nominal pressure or average interfacial separation global-scale variation. In order to achieve this, the extended domain is obtained for a range of global nominal pressure or average interfacial separation. The global scale is then obtained by interpolating these extended distributions in accordance with the global-scale trend.
(v) Permeabilities are assigned to each point in the global scale. To do so, the parameters of a log-normal distribution are first assigned to each point by interpolation based on the reference parameter (pnom or ). Then a permeability value is generated using that distribution.
(vii) Points 5 and 6 are repeated until good estimates for the leak rate and its variance are achieved. If the leak rate does not follow a normal distribution, the 95% CI can be computed by performing a large-enough number of realizations so that 95% of them have a leak rate which consistently falls inside a fixed range.
It is important to note that there will be cells with zero local permeability, occurring when there is no available path connecting the inlet and the outlet in the local domain. Those cells cannot be included in the log-normal distribution. Instead, the fraction of cells with zero permeability is stored. Then a number of cells corresponding to this fraction is set to zero in the global scale.
The reference parameter should be chosen between the nominal pressure and the average interfacial separation as the one that better correlates with the permeability. For example, for the turned surface with the topographies depicted in figure 3, permeabilities in the x1-direction correlate well with nominal pressure (except for when contact is lost) while the one in the x2-direction correlate better with the average interfacial separation.
(c) An average permeability approach
Most of the most successful two-scale models for fluid flow are based on the homogenization technique. Therefore, we benchmark our flow model to it. We first note that the proposed boundary conditions make the problem equivalent to that of standard homogenization. Indeed, by making the change of variables and , the local problems read 2.18aand 2.18band the flux reads 2.19which is equivalent to the homogenized results by reinterpreting the flow factors a11 and a21 as the local permeabilities K11 and K21 (e.g. ). We note that this resemblance was already studied in  and that a more rigorous derivation of the consistency between the HMM formulation and the standard homogenization is given in  for the more general time-dependent case.
The main difference lies then in the assumptions regarding the local-scale roughness. In homogenization it is assumed to be periodic, therefore, it is natural to assume that the flow factors (or permeability) are constant or vary smoothly in the global scale or, at least, in sub-domains of the global scale much larger than the global-scale size. One then obtains the values by averaging several cell realizations. Hereafter, we will refer to this procedure as the averaged permeability approach. By contrast, the periodic boundaries in the given approach represent only a consistent way to couple the two scales. This allows accounting naturally for spacial variation of permeability. We have used this feature to introduce a stochastic representation of the surface in §2b(iii).
In the average permeability approach, the pressure-driven flow problem considered will lead to a constant pressure gradient in the x1-direction and a zero pressure gradient in the x2-direction. The total leak rate, QH, is therefore 2.20where is the average permeability in the x1-direction.
The performance of the model is shown in two steps. First, the two-scale model is validated against a deterministic solution. Then, an example of its usage, the model is applied to a test case.
In both cases, two different topographies are used, as shown in figure 3. Topography 1 is highly anisotropic and corresponds to a turned surface. Topography 2 has been sand-blasted to obtain a more isotropic surface (although still with a clear anisotropic component). These topographies have also been chosen in order to obtain two clearly different total flow rates. The leakage studied is the one occurring between these two surfaces and a flat surface. The boundary conditions posed in the global scale (and in the deterministic solution) are periodic in the x2-direction and Dirichlet in the x1-direction, in order to impose a pressure drop Δp. The topography used is mirrored in the x2-direction in order to avoid a jump in topography when applying the periodic boundaries.
The material, assumed to be linear elastic-perfectly plastic, is described by the parameters E1=E2=206 GPa, ν1=ν2=0.3 and H=2.75 GPa. The leak rate is computed by (2.3), (2.6) or (2.20), depending on the methodology. In the presented results, however, it is scaled by the factor η/Δp. The Gaussian filter used for coarsening have a cut-off length 0.01 times the total length of the measurement and the filtered surface is re-sampled to one-sixteenth of the previous number of points. The original lateral resolution for the contact mechanics (both for the deterministic and the local-scale cells) is of 0.896 μm, which corresponds to a grid size of 600×1440 nodes for the topographies presented in figure 3. The tolerance for the load is set to 1×10−3% and the one for the contact plane is set to 10−10 m (see points 4 and 8 of the algorithm presented in  for reference).
(a) Two-scale model validation
In order to do the validation the measured domain Ω is separated into a set of local cells, ω. The size of the local cells in the x1-direction, Δx1 is taken so that each cell corresponds to one wavelength of the main frequency in Topography 1 (figure 3). In order to facilitate comparison, the same sizes have been also used for Topography 2. This means that the measurement domain is divided into six cells, each of which with a length of 0.21 mm. In order to verify convergence with the domain size, three widths in the x2-direction, Δx2, are taken. These are 0.18, 0.09 and 0.045 mm respectively, giving 6, 12 and 24 local-scale cells in that direction. A typical two-scale formulation requires the local domain to be at least one order of magnitude smaller that the global domain , which is not fulfilled in the present validation. This is because of the restriction on the measurement domain size. We note, however, that the error obtained when using a larger domain is expected to be smaller. Therefore, the validation can be directly extended to the more common situation where the two-scale separation is clearer.
The three techniques presented in §§2b and c are compared. In all three cases, the contact mechanics model described in §2 is used to compute the deformed gap. The deterministic approach for contact mechanics is used to compute the deterministic flow, while the two-scale contact mechanics approach is used for the other two techniques. The results for a range of total load W and for the three different widths Δx2 are depicted in figure 4.
The present two-scale model gives a solution that is in good agreement with the deterministic one for both tested topographies, specially for the two wider local cells. The average permeability approach predicts generally a higher leak rate. Although it could be acceptable for Topography 2 when using the higher width, this is not true for Topography 1. Moreover, for a given cell size, the present two-scale model gives significantly better accuracy. The error between the present two-scale formulation and the deterministic solution, defined as 3.1is depicted in figure 5. It can be observed that reasonable values are obtained for the two wider widths (maximum error for Topography 1 using cells of width 0.18 and 0.09 mm is 7 and 11%, respectively, and 5.5 and 19% for Topography 2 and the same widths). In order to identify the main source of error, the error between the deterministic solution and the solution with the present two-scale fluid model using the deterministic computation for the gap as the flow domain is also depicted. For Topography 1, it can be seen that it is much smaller, which identifies the two-scales contact mechanics model as the main source of error for this case. For Topography 2, however, this is not the case. The reason given is that due to the cell selection made for Topography 1, the boundary conditions in the flow model are particularly well suited, as the pressure becomes nearly constant at the edges of the domain due to the larger gap. This is, however, not the case for Topography 2. Thus, it becomes clear that the error and its main source is topography dependent. Despite that, given a careful choice of domain size, a sufficiently good approximation can always be made. The average permeability approach based on the homogenization technique has been proved to produce good results for similar pressure-driven conditions . The discrepancy observed in this work is, therefore, attributed to the small size of the local cells, which makes them not representative of the topography. This is further enhanced when the leak rate is reduced. In order to see why the local-scale cells are not representative when using the average permeability approach but are representative in the present two-scale formulation, one can analyse the flow pattern. The flow pattern obtained using Topography 1, in terms of deterministic flux, is depicted in figure 6a for a high total load W, i.e. for a small leak rate. It can be seen that the fluid follows a particular pattern through the gap. Whenever a path is available, it advances in the x1-direction, which is the direction of the pressure gradient. When this is not possible, or when the path is too small, the flow advances in the x2-direction until a large path is found and flow in the x1-direction is possible again. Furthermore, the advance in the x2-direction is as large as the half of the width of Ω, which is the maximum value due to the surface has been mirrored. In fact, one should expect even larger advances in this direction when computing over a larger domain. This implies that a representative cell for the average permeability approach would need to be at least as large as Ω for it to be possible to capture the flow pattern correctly. In the present two-scale formulation, however, the heterogeneity in the permeability distribution can enforce such flows in the x2-direction even when using small local cells. This can be seen in figure 6b, where it can be observed that the flow pattern is satisfactorily captured by the present two-scales approach. Another fundamental difference of the two approaches lies in how much large values of permeability affect the total leak rate. In the average permeability approach the cells with high permeability will increase the average permeability and therefore significantly affect the total leak rate. In the present two-scale formulation, the leak rate is controlled by the narrowest constriction that must be crossed. This is, in fact, a more correct representation of the problem . Therefore, the influence of the high permeability values is not as large as in the average permeability approach. This is also the explanation for the observed differences in the leak rate prediction shown in figure 4. In figure 7, a similar comparison, using Topography 2, is depicted. The present two-scale approach captures, again the correct flow pattern. In this case, however, the flow in the x2-direction is much less important (although significant enough to affect the results) and the constrictions are more evenly spread. This, in turn, explains why the average permeability approach does not deviate as much as in the previous case.
It is important to emphasize that the average permeability approach only predicts a too high leak rate when it is small, i.e. at high total load. Otherwise, the variation in the permeabilities are not expected to be large and the flow perpendicular to the direction of the pressure drop is expected to be less important. In those cases, the average permeability approach can be considered a good approximation. This can be seen in figure 4, where the prediction from the average permeability approach is close to the deterministic one at the lowest total load and the deviation observed when using Topography 2, which exhibits a higher leakage, is smaller.
(b) Two-scale stochastic model
Once having validated the two-scale part of the model, we show in this section the performance of the model by applying it to a case example. The measurement domain used for this is the same as in the previous section.
Let’s start by considering the local-scale results. Having computed the permeabilities for the set of local cells, these values must be fitted to log-normal distributions. In order to decide which reference parameter to use, we study the correlation between the permeability and the two reference parameters (nominal pressure and average interfacial separation). These relations are shown in figure 8. Focusing first on the correlation between average interfacial separation and permeability, one can observe that the permeabilities in the x2-direction can be approximated by . The same trend is followed by permeabilities in the x1-direction when there is no contact and the separation is large. This approximation will be used for simplicity when obtaining the permeability distributions in the global scale. When the nominal pressure increases, the permeability in the x1-direction is reduced significantly (note here that a large number of cells have zero permeability at high nominal pressures), while the average interfacial separation remains nearly constant. Therefore, we use the nominal pressure instead of average interfacial separation as the reference parameter. As seen in figure 8, there is a large spread in the correlation between permeability and nominal pressure, but one can clearly observe a trend in both mean value (decreasing) and spread (increasing) of permeability. The percentage of cells with zero permeability also increases with increasing nominal pressure. In order to describe the permeability in the x1-direction log-normal distributions are fitted for each value of the nominal pressure.
The log-normal distributions relating permeability and nominal pressure and the approximation are used to build the permeability distributions and the global scale leak rate is computed. An example of a permeability distribution is shown in figure 9 and the resulting flow pattern for a unitary pressure drop in the x1-direction is depicted in figure 10. It can be seen in figure 10a how the flow advances in the x2-direction until an easy path is found in the x1-direction. For the case of Topography 2 (figure 10b) the leakage is much higher, resulting in a less important flow in the x2-direction. One could define this flow as advancing in the x1-direction until a blockage is found, then doing a small shift in the x2-direction to avoid it and continue in the x1-direction.
One can think on these flow patterns, particularly the one corresponding to Topography 1, similarly to the description given by Persson & Yang , where the flow domain is described as a network of paths with critical constrictions randomly located. In the present model, however, no assumption on either the distribution or the size of the paths and the constrictions needs to be imposed.
The leak rate is computed for several permeability distribution realizations. When a sufficient number of realizations have been computed, the leak rate can be described statistically. The probability distribution for the leak rate for three global domain sizes can be seen in figure 11a for Topography 1 and in figure 11c for Topography 2. The results for Topography 1 have been computed using cells of width 0.09 mm, resulting on global domains of 12×48, 12×384 and 96×384 local cells in the x1- and x2-directions, respectively. Cells of width 0.18 mm have been used for Topography 2. When the leak rate is high, or when the global-scale domain is large compared with the local domain size, the distribution tends to be normal. However, as either of the leak rate or the global-scale domain become smaller, the variance of the leak rate distribution increases and it becomes skewed. If a normal distribution can be assumed for the leak rate, it can be characterized by using only the mean value and the standard deviation of that distribution. In the more general case where the distribution is not known, a 95% CI for the leak rate can be computed, as shown in figure 11b and d. Without knowing the distribution, this interval is obtained by performing a large number of computations until 95% of the computed leak rates fall consistently inside a fixed range. Some general trends can be observed with increase of global-scale domain. As it is increased in the x2-direction, a reduction of variability of the leak rate is observed. Both the reduction of variability and the tendency to become normal when increasing the size in the x2-direction are expected results as the effect of extreme values will be smaller and all realizations will be more similar due to the presence of more local cells. It is also notable that the leak rate per unit width is slightly increased in the case of Topography 1. This is because of the larger domain allow areas with high permeability to appear more often in the realizations. This effect is not observed in Topography 2 because its larger leakage allows for a better averaging even at the smallest domain. An increase of the global domain in the x1-direction has the opposite effect and mean leak rate is reduced. The reason for this is that the fluid must cover a longer distance and a narrower constriction is more likely to be encountered. One should expected the leak rate per unit width to converge to a certain value as the size of the domain in the x1 and x2-directions increases. It is noticeable, however, that convergence with size in the x1-direction is far from being reached even at realistic seal dimensions.
4. Concluding remarks
In pressure-driven flows, as the leak rate decreases, the fluid follows a pattern which becomes wide in the direction transverse to the pressure gradient. This obliges utilizing a large local domain in the common two-scale formulations. However, the required size of the local domain can become too large and approach the global domain, which makes it not possible to perform a scale separation of the flow. In order to avoid this problem, a two-scale model based on the HMM framework has been developed. This model does not assume periodic repetition of the topography and, therefore, can capture this wide flow pattern via local variation of permeability. This allows using much smaller local-scale domains. It has been shown that in this way a two-scales model which captures correctly the flow pattern can be developed.
Moreover, the presented two-scale formulation permits the construction of the local scale by considering the local permeabilities a random variable. This feature allows inclusion of the inherent uncertainty coming from the surface topography in the model and estimating the uncertainty on the leak rate due to this cause.
In the case study presented, it has been shown that smaller global-scale domains, as well as smaller leak rates, lead to more uncertainty on the predicted leak rate. It has also been shown that, even for relatively large global domains, the leak rate per unit width is dependent on the domain size. This has been explained as an effect of the random construction of the permeability distributions.
The datasets supporting this article have been uploaded as part of the electronic supplementary material.
F.P.-R. ran the simulations, performed the analysis of the results and drafted the manuscript. S.L. contributed on the designing and analysis of the stochastic part of the model. P.W. contributed in the mathematical formulation of the problem. A.A. and R.L. designed the contact mechanics and flow codes, and participated in the analysis in the results. All of them revised the text and gave final approval for publication.
We have no competing interest.
The authors have been funded by the Swedish Research Council.
- Received January 28, 2016.
- Accepted May 17, 2016.
- © 2016 The Authors.
Published by the Royal Society under the terms of the Creative Commons Attribution License http://creativecommons.org/licenses/by/4.0/, which permits unrestricted use, provided the original author and source are credited.