Enhancing Full-Film Lubrication Performance via Arbitrary Surface Texture Design ∗

Minimizing energy loss and improving system load capacity and compactness are important objectives for fluid power systems. Recent studies reveal that micro-textured surfaces can reduce friction in full-film lubrication, and that asymmetric textures can reduce friction and increase normal force simultaneously. As an extension of these previous discoveries, we explore how enhanced texture design can maximize these objectives together. We design surface texture using a set of distinct parameterizations, ranging from simple to complex, to improve performance beyond what is possible for previously investigated texture geometries. Here we consider a rotational tribo-rheometer configuration with a fixed textured bottom disk and a rotating top flat disk with controlled separation gap. To model Newtonian fluid flow, the Reynolds equation is formulated in cylindrical coordinates and solved using a pseudo-spectral method. Model assumptions include incompressibility, steady flow, constant viscosity, and a small gap height to disk radius ratio. Multiobjective optimization problems are solved using the epsilon-constraint method along with an interior-point nonlinear programming algorithm. The trade-off between competing objectives is quantified, revealing mechanisms of performance enhancement. Various geometries are explored and optimized, including symmetric and asymmetric circular dimples, and novel arbitrary continuous texture geometries represented using two-dimensional cubic spline interpolation. Shifting from simple dimpled textures to more general texture geometries resulted in significant simultaneous improvement in both performance metrics for full-film lubrication texture design. An important qualitative result is that textures resembling a spiral blade tend to improve performance for rotating contacts.


INTRODUCTION
Friction is a significant source of energy loss in mechanical components.The influence of surface roughness within lubricated hydrodynamic contacts has been well-studied in the tribology community.Creating micro-dimples at frictional contact interfaces is known to be effective in reducing frictional losses [2][3][4].Furthermore, recent studies proposed a possibility that textured surfaces for sliding hydraulic interfaces can reduce effective friction, improve sealing, and increase the load capacity simultaneously by creating texture in favorable shapes and dimensions [5][6][7].However, the roughened or dimpled surfaces in previous studies were restricted to relatively simple shapes, such as predefined shapes, randomly roughened surfaces, and macroscopic sinusoidal textures [4,[8][9][10][11][12][13][14][15][16].Although most existing studies in the literature include either symmetric texture profiles, or asymmetry only for the top profile of texture rim perimeters, Refs.[13,14,17] showed that generating normal force with an asymmetric depth profile is more effective than with a symmetric depth profile.Due to recent manufacturing advancements, such as additive manufacturing and electric/electrochemical micromachining, textures with more general geometries and scales are realizable for hydrodynamic lubrication surfaces [5,18,19].Here more general texture geometries are investigated.
Film lubrication involves lubricant flow, which is often modeled using the fluid flow governing equations (e.g., Navier-Stokes equations).Solution of general forms of these equations is computationally expensive, limiting the utility of such models for design studies [6,20].In hydrodynamic lubrication, however, the distance between sliding surfaces is much smaller than other length scales in the system; this allows the governing equations for flow to be simplified to the Reynolds equation, which is more favorable for computation [21].In a recent study, the Reynolds equation was solved for the fullfilm lubrication problem with circular dimples using the pseudo-spectral method [22].This numerical method exhibits exponential convergence with increasing order of polynomial basis function, is computationally efficient, and has been validated against experiments [22,23].
Here we are building upon this pseudo-spectral solution method for the Reynolds equation to explore a much wider range of surface texture designs with the primary objective of friction reduction.Arbitrary continuous texture shapes are optimized using effective parameterization techniques.Important existing options for describing curved shape designs include analytical approaches, basis vectors, domain elements, Fourier descriptors, polynomial functions, response surfaces, and splines [24,25].These methods are often used in aeronautical system design (e.g., airfoil shape optimization [24]) and other applications (e.g., photo-voltaic internal reflection texture design [26]).The specific parameterization used in this study involves spline functions defined in two orthogonal directions to represent arbitrary texture shapes in three dimensional space.

PROBLEM DESCRIPTION
The objective of the design problem in this study is to maximize the film lubrication efficiency and effectiveness by designing the shape of the textured surface using systematic methods.The experimental setup in Schuh and Ewoldt [17] revealed that asymmetric-depth-profile dimpled surface textures decrease frictional loss.Figure 1(a) illustrates the problem setup used for the previous experimental research.Based upon this configuration, we perform a more comprehensive computational study here enabled by more flexible texture surface design representation and exploration.A pair of gap-controlled disks are aligned axially, and are separated by a Newtonian fluid.The rotating upper disk is flat, while the stationary bottom disk is textured.All points of the surface in the design domain-illustrated in Fig. 1(b)-are defined using a moderatedimension spline representation.The design domain is a periodic sector of the bottom disk instead of the entire disk.The full disk is divided into N p sectors.The number N p is an arbitrary choice, and N p = 10 is selected for the studies presented here.Future work may involve adjusting N p , or alternate domain representations, but is outside the scope of this article.Simulation results are equivalent to the behavior of a full disk with repeated sectors.This simplification reduces computational expense and is assumed to be reasonable due to the rotational nature of the setup.The sector surface design is assumed to be periodic, i.e., the texture repeats and the boundaries of each sector match to preserve texture continuity.This is consistent with the N p = 10 repeating dimples in the experiments [17].Figure 2 presents a visualization of how a repeated sector represents a periodic texture design for a full disk surface.While the design optimization process and the flow simulation are computed using the single-sector design domain given in Fig. 2(a), the results correspond to physical behavior for a full disk shown in Fig. 2(b) due to periodic boundary conditions.

SURFACE PARAMETERIZATION METHOD
In previous related work, surface textures have been described using very simple parameterizations.For example, one strategy assumes that surface texture features consist only of cylindrical dimples, either with flat or angled lower surfaces, and are parameterized using dimple diameter, depth, and lower surface angle [6,17].Another strategy prescribes a set of allowable surface textures (e.g., circle, ellipse, triangle) and selects among them by comparing performance [8,12,13].Another study designed texture with respect to its top profile view, but maintained the gap height inside the texture at a fixed value [27].This approach provides a general outline shape as a texture boundary, but still does not represent a general height profile for texture design.
A core objective of the present study is to perform design exploration of much more general texture designs as a means to gain greater insight into surface texture design for enhanced efficiency.One strategy is to simply use h i j , the surface height at computational mesh nodes, as the surface design description.While this provides high accuracy and a high-resolution design description, it results in a large-dimension nonlinear optimization problem.Our current implementation of the optimization problem solution requires treating the simulation (i.e., in this study, lubricant flow computation) as a black box, necessitating finite difference calculations.While ongoing work is focused on investigating alternative implementations that leverage problem structure, such a large-dimension design representation is impractical for the present study.In previous work it was determined for this system configuration that the coarsest accurate mesh for simulation involves 26 × 26 nodes (polynomial order N = 25, see Sec. 4.1 for detail) [22].Using h i j directly, accounting for periodic boundary constraints, would therefore require (25 + 1) × 25 = 650 optimization variables.Therefore, a reduced-dimension design representation is needed.As detailed below in Sec.3.4, we use a two-dimensional cubic spline on a coarse mesh (Fig. 3(a)) for design representation, and interpolate to a finer mesh for simulation (Fig. 3(b)).The mesh for a reduced-dimension cubic spline texture shape design representation is shown on the left.This design mesh is much more coarse than the mesh on the right that is required for accurate simulation.The surface geometry defined by the cubic spline is interpolated to determine all the height values at the fine mesh points required for simulation.
We consider several types of textures, of varying design freedom, shown in Fig. 4. The simple angled cylindrical dimple textures from previous work are illustrated in Figs.4(a-b).The aim here is to support exploration of more general texture design shapes, as shown in Figs.4(c-e), with reduced dimension, and to evaluate performance improvements available through these more general texture designs.The simplified cylindrical texture parameterization is reviewed first, followed by a description of a more general two-dimensional spline representation.James T.

Cylindrical textures
The previous experimental setup involved a single cylindrical dimple with fixed location in each sector.The asymmetric texture angle was varied in those previous studies to gain an initial understanding of this behavior [6,17].Here we build upon these initial studies by optimizing this dimpled texture design.The complete set of cylindrical dimple parameters is listed in Table 1.Design variables include radius and depth (x 1 , x 2 ) for symmetric textures, or radius and angle (x 1 , x 3 ) for asymmetric textures.Geometric configurations of these textures are illustrated in Figs.4(a-b) along with the shape parameters.We used a nominal gap h 0 of 0.269 [mm] to compare to experiments; this value of h 0 is fixed and used as lower bound on gap.Also, see Fig. 4(b) for a visualization of how the inclination β in Table 1 is defined in the cylindrical texture with an asymmetric height profile.

Inclined plane spanning full disk sector
We will see that larger radius textures improve performance.With this observation, we can predict that if the texture area spans the full sector area, it may be possible to further improve performance.To demonstrate the effect of expanded texture area, a new geometric parameterization is defined where the entire sector is an inclined plane tilted at angle β, as shown in in Fig. 4(c).Note that this geometry is continuous, but non-smooth.

Other low-order design representations
Polynomial texture: A few additional low-order design representations were investigated in addition to the cylindrical and inclined plane texture parameterizations.First, a uni-modal polynomial function was used to generate the height profile by specifying location of peak function value.Height at both periodic sides and the inner/outer boundaries was fixed to the nominal gap height (h 0 ).Numerical experiments indicated that this class of texture performed poorly, possibly because large local slopes could not be achieved.Detailed results for this case are omitted for brevity.
Radial-basis function texture: Radial-basis functions (RBFs) are used widely in approximating or interpolating functions.RBFs were tested as a texture geometry representation where the height and location of a number of thin-spline RBFs were used as design variables.Several challenges were discovered.Changing texture design in significant ways requires changing the number of RBFs used, which cannot be done during the optimization solution using continuous algorithms.In addition, the number of parameters required is large relative to the range of texture geometries that are accessible.Each RBF requires at least three parameters.RBFs do not provide a low-dimension representation, and do not support efficient design space James T. Allison 4 MD-16-1636 exploration for this problem compared to the spline representation discussed next.Detailed results for this case are also omitted for brevity.

Two-dimensional cubic spline interpolation
Here we introduce a low-dimension texture design representation that supports description of arbitrary continuous geometries.This provides enhanced design flexibility and the possibility of capitalizing on new mechanisms for improving fluid system performance.To describe the full-sector texture design using a limited number of design variables, the height profile is specified at nodes of the coarse mesh shown in Fig. 3(a).As a manufacturability constraint, the inclination angle between coarse mesh nodes is limited to a maximum of 30°.After specifying a low resolution height profile, the two-dimensional spline representation is used to map this profile onto the fine computational mesh shown in Fig. 3(b).This allows the surface texture design to be specified using a limited number of design variables, while still supporting high-resolution simulation.Coarse mesh node heights are the only design variables in this parameterization.Because of the periodic constraint at the left and right sector sides, the design variables include the height profile of only one of two sides.Thus, for (N + 1) × (N + 1) nodes, the number of design variables is (N + 1) × N.Alternatively, a texture symmetry constraint can be imposed.With symmetry, the number of design variables for (N + 1) × (N + 1) nodes is (N + 1) 2 /2 , where • is the ceiling function.Symmetric designs are investigated first, as shown in Fig. 4(d), followed by studies of asymmetric designs, shown in Fig. 4(e), that further improve system performance.

FLOW SIMULATION METHOD 4.1 Reynolds equation
We model the flow following Schuh et al. [23].This model was based on physical experiments with surface textures and Newtonian fluids where cavitation effects were not observed [17].Therefore, the model chosen does not include cavitation, and any use of cavitation is outside the scope of this paper.The Reynolds equation, given in Eqn.(1) below, satisfies both mass and momentum conservation for an incompressible Newtonian fluid under the assumptions that 1) inertia terms are negligible, and 2) the gap is small compared to other length scales in the geometry [28].The pressure field p (r, θ) from Eqn. (1) is used to calculate the velocity field u (r, θ, z) in Eqn.(2).Boundary conditions for velocity are u θ = rΩ and u r = u z = 0 at z = 0, and u θ = u r = u z = 0 at z = h.We use steady state conditions with fixed angular velocity Ω = 10 [rad/s] and constant viscosity η 0 = 1.4 [Pa•s].The Reynolds number, based on the nominal gap, is Re = ρΩR o h 0 /η 0 = 0.033, and this value can be applied to the predefined dimple shape texture cases.For the arbitrary continuous texture cases, the Reynolds number is instead computed based on the range of the gap, and can vary from 0.033 to 0.186.These values correspond to the experimental conditions of Schuh and Ewoldt [17], but the approach is general and results can be made non-dimensional by scaling with respect to viscous effects (p ∼ η 0 RΩ/h 0 ) [5].
The details of simulating the partial differential equation (PDE) in Eqn.(1) are described in Ref. [23].Briefly, the equation is discretized using the pseudo-spectral method, which is a variation of the weighted residual technique (WRT), where the PDE is solved in its variational form [29].The resulting integrals for the variational form of the PDE are solved using Gauss-Lobatto-Legendre (GLL) quadrature, where the function evaluations occur at the GLL points and the quadrature weights are chosen optimally such that integrals of the solution are exact for polynomials of degree 2N − 1, where N is the number of evaluation points (shown as nodes in Fig. 3(b)).For all case studies presented here, a fixed number of computational mesh points, as illustrated in Fig. 3(b), was used for the Reynolds equation solver.This computational mesh density is validated to be sufficient for accurate objective function prediction [22].Since the Reynolds equation model is solved using the same mesh for the entire set of studies, its computational cost and prediction accuracy are independent of design representation fidelity.Full details of this methodology are given in Refs.[22,23,29,30], where the model predictions have been validated against experiments from Refs.[17,22].
A previous study emphasized that the dimensionless parameters are important for interpreting the physics of lubricated sliding contacts [5].However, we do not use dimensionless parameters in flow simulation or in the optimization formulation.Analyses with dimensional variables in this study are still meaningful because similarity is maintained for the flow within the range of Reynolds numbers we are modeling.Results could be made dimensionless by scaling with respect to viscous effects [5].

Boundary conditions
The Reynolds equation model for this study predicts the pressure field with assumptions that recirculation is not present in the velocity field, that the r and θ direction velocity field is a linear combination of simple shear and pressure driven flow, and that the pressure is invariant in the z direction.In the experiments, the outer edge of the textured domain was an oil-air interface with the air at atmospheric pressure.Boundary conditions on the pressure and the velocity at the outer edge can be derived using conservation of normal and tangential stress.It is not possible, however, to implement the boundary conditions on the velocity field directly, because the Reynolds equation only defines and manipulates the pressure field (unlike the Navier-Stokes equations).However, we have explicit equations for velocity fields in terms of pressure as given in Eqn.(2).Therefore, boundary conditions on velocity need to be converted to pressure conditions.
In the design studies presented here it is assumed that texture variations extend all the way to the outermost boundary of the sector.The boundary condition at the outer edge can be a choice between p = 0 (Dirichlet boundary condition) or ∂p/∂r = 0 (Neumann boundary condition).The choice has consequences for the r-velocity u r boundary conditions [23].If p = 0 at the outer edge, then ∂p/∂θ = 0 and this eliminates the pressure driven flow in the θ direction (first term in Eqn.(2b)).Then, if ∂h/∂θ = 0, i.e., texturing, at the outer edge, we have ∂u θ /∂θ = 0, and thus by mass conservation ∂u r /∂r = 0.This means the r velocity component can be non-zero at the outer edge for the condition p = 0.This non-zero r velocity would cause fluid to leave the disk-shaped textured domain (which was not observed experimentally), and would also result in a non-zero shear stress component (τ rz = η (∂u r /∂z)) on the outer free surface.To enforce u r = 0 at the outer edge (no flux), the gradient of pressure in the r direction at the outer edge must be zero (Eqn.(2a)), i.e., ∂p/∂r = 0, the Neumann boundary condition, and ∂p/∂θ is unspecified.The use of periodic and Neumann boundary conditions results in a pressure field distribution quantified in terms of relative pressure with respect to an arbitrary fixed pressure at a certain location.We eliminated this arbitrary shift in the pressure profile by constraining the average value of the pressure to be zero at the outer edge of the texture.This effectively assumes no net pressure drop across the liquid-air interface, e.g., due to surface tension at a curved interface (Laplace pressure drop).This was also eliminated experimentally by calibrating and subtracting this effect [17].

OPTIMIZATION METHOD 5.1 Multiobjective Optimization
Multiobjective optimization involves minimization or maximization of a set of multiple conflicting objective functions.The solution of such a problem is a set of non-dominated (Pareto-optimal) solutions, as opposed to a single optimum point as with single objective optimization.A design point is non-dominated if one objective cannot be improved without degrading at least one other objective [31].
Two primary classes of methods are used to solve multiobjective optimization problems: population-based or scalarizationbased methods.Population-based methods, such as multiobjective genetic algorithms (MOGAs), solve the optimization problem once and generate a set of solutions that form the Pareto set (approximately).These methods often improve the probability of finding global instead of local optima, but may be computationally expensive due to the large number of function evaluations typically required [32].
The second class of methods converts a single multiobjective problem into a set of 'scalarized' single objective problems.These scalarization-based methods are often used in design optimization studies due to their computational efficiency and simplicity [33].The weighted-sum method is the simplest scalarization approach, but cannot identify non-dominated solutions in non-convex regions of a Pareto front.In addition, resulting non-dominated points are often clustered instead of uniformly distributed across the Pareto surface.These limitations can be overcome through other more sophisticated scalarization approaches, such as the ε-constraint method [34], which is illustrated in Fig. 5 and described in Section 5.3.We used the ε-constraint method in this study to enable identification of non-convex Pareto fronts while supporting use of computationally-efficient gradient-based optimization methods.

Multiobjective Formulation
The conflicting objectives of the full-film lubrication problem considered here are to 1) minimize the shear load, represented as normalized apparent shear viscosity η a /η 0 , and 2) to maximize the normal force load F N .Friction is a significant source of energy loss for systems involving lubricated hydrodynamic contacts.Friction can be reduced by increasing the gap between the sliding contact surfaces.Increasing this gap, however, degrades load capacity and sealing performance.
Fluid and pressure losses due to poor sealing are very undesirable outcomes, e.g., for hydraulic power systems.Improved sealing requires increased normal force, but increased normal force may increase frictional losses (e.g., fluid squeeze-out, higher contact friction, etc.) [35].These two objective functions therefore conflict, and the solution to this multi-objective design problem will be a set of non-dominated points that quantify trade-off options.Adding certain types of textures to the full-film lubrication problem was found in previous studies to simultaneously improve both of these objective functions.In other words, transitioning to more effective texture classes shifts the Pareto surface toward more desirable objective function values.These previous studies were limited to very simple uniform dimpled textures [4-6, 8, 22, 23].Here we aim to shift the attainable Pareto set even more by considering more general texture designs.We also aim to gain fundamental insights about how best to design surface textures for full-film lubrication applications.
The multiobjective problem formulation used here is: where x is a vector of design variables that represents texture geometry (e.g., spline, etc.) and fluid properties are fixed.In this study, geometric parameterizations are formulated to implicitly satisfy the periodic boundary conditions for the surface shape.For a given value of x, the resulting geometric surface description is then used to determine high-resolution surface height values: h i j (x), ∀{i, j} ∈ D. The height values h i j (•) quantify the surface shape at mesh points needed for the pseudo-spectral method, D is the set of indices for all mesh points in the design domain, and {i, j} are node indices in the two-dimensional mapped mesh space.The objective functions ultimately depend on x.The height values, obtained from x, are used within the simulation to obtain intermediate quantities needed to compute the objective function values.The simulation solves for the pressure distribution p(r, θ) and velocity field u r (r, θ, z), u θ (r, θ, z).
The first objective function, normalized apparent viscosity η a (x)/η 0 , is calculated by a ratio of apparent viscosity (η a (x), depends on design) to actual fluid viscosity (η 0 , fixed, depends on fluid selection) [36].The apparent viscosity is defined from the torque M on the rotating disk (integrated shear stress).For the rotating disk configuration, apparent viscosity can be written as a function of disk torque M. where, R o is the outer radius of textured disk; h 0 is a controlled minimum gap height between the fixed and rotating disks; Ω is the rotating disk angular velocity.Equation (4b) defines how disk torque is calculated, which requires evaluation of the θ-direction (tangential) shear stress (τ θz (r, θ), where z is the vertical coordinate) at the rotating surface across the complete domain.The shear stress calculation requires knowledge of the velocity field, which for the Reynolds equation depends on gradients of the pressure field (∇p(r, θ)) across the entire computational domain of the flow field.For example, the shear stress in the theta direction can be calculated as a function of ∂p (r, θ) /∂θ as given in Eqn.(4c).The second objective function, the negative normal force (−F N ), is calculated by integrating pressure over the domain, and then multiplying this value by the total number of disk sectors (N φ ) as shown in Eqn.(5).

ε-constraint method
The ε-constraint method is a type of scalarization technique for multiobjective optimization.A multiobjective problem is transformed to a single objective problem by retaining just one of the original objective functions, and the remaining objective functions are converted to constraints that bound these other objective function values [34].Figure 5 illustrates how the ε-constraint method minimizes the first objective function ( f 1 (•)), while a constraint prevents the second objective function ( f 2 (•)) from exceeding a bound denoted by red horizontal lines in the objective function space.The objective function space is a multidimensional space f ∈ R m , where R is the real number space and m is the number of objective James T. Allison 7 MD-16-1636 functions.This procedure is repeated, each time with a different bound on f 2 (•).In the example shown in Fig. 5 the optimization problem is solved five times to generate five Pareto-optimal solutions.This strategy supports the use of existing single-objective optimization algorithms in solving multiobjective problems, including the ability to resolve non-convex portions of the Pareto surface, and to generate well-distributed Pareto-optimal points.The multiobjective formulation in Prob.
(3) can be reformulated for ε-constraint solution as follows: where n p is the number of Pareto-optimal points to solve for.One possible strategy, which is used in this study, is to increment (2) (5) (5) The ε-constraint method converts a multiobjective optimization problem to multiple single objective optimization problems.
For each scalarized problem, one objective is minimized while satisfying bound constraints on the other objectives.
6 RESULTS AND DISCUSSION 6.1 Cylindrical textures A full factorial set of cylindrical texture designs, both in symmetric and asymmetric configurations, was generated and evaluated to provide insight into how design variables influence objective functions, and to assess trade-offs.Texture radius (R t ) and depth (h) are the design variables used for the symmetric cylindrical texture study, and texture radius (R t ) and angle (β) are the design variables used for the asymmetric cylindrical texture study.Figure 6 illustrates the contours of the two objective functions for both cases.Figures 6(a-b) correspond to the symmetric texture, and Figs.6(c-d) correspond to the asymmetric texture.
For the symmetric textures, shear load (normalized apparent viscosity) decreases with "big" textures (i.e., larger radius and deeper depth) down to η a /η 0 = 0.714, while no normal force is observed.Thus, the lowest shear load determines the optimal design (circles in Fig. 7).For asymmetric textures, the angle β determines texture depth (β = 0 is a flat surface).Shear load is still minimized for "big" textures (i.e., large radius and large β).Normal force is expected to be maximized at an intermediate angle β [17,23], and here this occurs at β = 3.6°.Normalized shear load ranges from 0.794 to 0.902, and normal force ranges from 0.212 to 0.762 N. Asymmetry is required to generate normal force.
Design variable bounds are chosen such that the Reynolds equation simulation results reasonably agree with results based on the full Navier-Stokes equation.Figure 7 shows both Reynolds equation and Navier-Stokes equation solutions using ANSYS Fluent in terms of normalized apparent viscosity for the symmetric cylindrical texture problem.When the texture depth is 2 mm and 4 mm, the solution of the Reynolds equation deviates 1.65% and 6.93%, respectively, from the solution of the full Navier-Stokes equation.Due to this limitation in the Reynolds equation model, the depth variable should be constrained to preserve accuracy.Furthermore, the full Navier-Stokes predictions indicate that normalized apparent viscosity plateaus beyond a certain depth (note that the plot is log-scale).This means that the frictional performance cannot be enhanced more through increases in depth beyond a certain value.
For both symmetric and asymmetric textures, the first objective function (minimizing normalized apparent viscosity) has monotonic dependence on texture depth and radius as shown in Figs.6(a, c) (the objective decreases with increasing depth and radius).In the asymmetric configuration study the angle was varied between 0-20°.Within this range, increasing the angle corresponds to increased average texture depth, and increased lubricant volume.This observation, along with the results presented here, indicate that the total volume removed for a texture is important for reducing friction.
The behavior of the second objective function (maximizing normal force) is different for the symmetric and asymmetric configurations.Symmetric cylindrical design variables do not influence normal force because geometric symmetry in the sliding direction results in negative and positive pressure distributions that countervail each other through the expanding-contracting channel gap height.This matches the results reported in previous experimental and theoretical studies that measurable normal forces were not detected in symmetric cylindrical surface textures [17,21].
The second objective function varies widely when changing the asymmetric configuration design variables.Figure 6(d) shows that the normal force depends on both the texture angle β and the radius R t .In particular, the best normal force value occurred at a specific texture angle value, β = 3.6°, with the texture radius R t at the upper bound, 4.0 mm.Moving away from this point improves normalized apparent viscosity, but degrades normal force.The result of this optimization study, therefore, is a set of non-dominated designs that express the trade-off between these two objective functions.We can improve normal force through asymmetric textures, but at the cost of degraded apparent viscosity.Detailed insight of the flow physics (shear stress and pressure fields) are shown in Fig. 8 and Fig. 9. Figure 8 compares symmetric and asymmetric optimal designs (minimum η a /η 0 for symmetric cylinder and maximum F N for asymmetric cylinder).The pressure field results in Fig. 8(c) show that the pressure is counterbalanced across the vertical line of symmetry.Thus, the overall pressure acting in the direction normal to sliding is canceled, resulting in no normal force.This numerical result explains why design variables for symmetrical textures have no effect on the normal force, as shown in Fig. 6(b).

Symmetric cylindrical texture design exploration study
For the asymmetric configuration, the pressure distribution is not counterbalanced and can produce a net normal force.Here the positive overall pressure generates a positive normal force.Figure 9 shows how the pressure distribution changes with angle β.At the angle of maximum normal force (β = 3.6°), the range of pressures is at its maximum.The maximum normal force increases with β up to β = 3.6°, but then decreases with angle beyond β = 3.6°.The interface between positive and negative pressure consistently moves from left to right with increasing angle.This is because when the angle of asymmetry is greater than 45°, the resulting surface texture is geometrically similar to a texture with an angle of (90°− β) with the flow in the opposite direction, and it has been shown previously that when the direction of motion changes for asymmetric surface textures, the sign of the normal force also changes [17,22].

Inclined plane spanning full sector area
Figure 10 shows a result of an inclined plane texture design that involves a slow expansion followed by rapid contraction.In the pressure distribution contour plot (Fig. 10(c)), a thick contour line (zero-pressure line) highlights the boundary between positive and negative normal force regions.The positive normal force region has a larger area than the negative normal force region.Integrating this pressure distribution over the sector area results in a positive net normal force.This is the general observation for cases where slow expansion is followed by rapid contraction.Figure 11(a) shows the Pareto fronts for the inclined plane and cylindrical textures.It compares design solutions in the objective function space, including both symmetric and asymmetric configurations for cylindrical textures.By expanding the inclined region of the texture, both objectives can be improved simultaneously.It is clear that normalized apparent viscosity is reduced significantly when shifting to the full inclined plane geometry.Since this design representation has only a single design variable, (β), it is possible to visualize the response of both objectives with respect to β as shown in Fig. 11(b).This plot shows a clear conflict between the two objective functions.Negative inclined angle β creates flow in rapid expansion followed by slow contraction, while positive inclined angle β results in flow in slow expansion followed by fast contraction.The sign of the inclined plane angle β determines the sign of normal force F N , but the magnitude of the normal force is symmetric between positive and negative inclined angle β.Furthermore, apparent viscosity η a /η 0 has even symmetry with respect to β.The region corresponding to Pareto-optimal designs is identified by the white area in Fig. 11(b).The gray area in this plot indicates the set of suboptimal designs that are dominated by the Pareto-optimal designs.(c-d).In (b), the two objectives (blue: η a /η 0 , green: F N ) are plotted together as a function of β for the inclined plane design.The gray area indicates the region of suboptimal designs.

Arbitrary continuous texture designs
The study is now extended to arbitrary continuous surface elevation changes parameterized using spline interpolation.Objectives of this investigation include understanding how to improve performance further through more sophisticated texture designs, and what physical mechanisms make any improvements possible.
Using the coarse mesh described in Section 3 and Fig. 3(a), designs are represented using a (N + 1) × N matrix of surface elevation values, accounting for periodic boundary constraints.For the reference case this matrix has dimension 6 × 5, with a resulting design representation dimension of 30 for asymmetric arbitrary continuous surface texturing.We first consider symmetric designs (Fig. 4(d)) which have a lower design representation dimension due to the symmetry constraint ( (5 + 1) 2 /2 = 18).Then, we allow asymmetric designs by removing the symmetry constraint (Fig. 4(e)).We also vary the design space resolution from N = 3 to N = 7 to examine the dependence of optimal performance on resolution.
James T. Allison 11 MD-16-1636 Figure 12 shows the optimal design with symmetry constraint.It is nearly a flat plate which maximizes the allowable gap height.As with the symmetric cylindrical texture result, the pressure distribution is balanced and the normal force is zero.As before, the problem reduces to a single objective function problem due to symmetry and insensitivity of normal force to texture design.This is a significant result.Even with arbitrary surface topography, symmetry must be broken to generate a normal force from textured surfaces.Gap reduction has no benefit in terms of normal force, so the optimization algorithm is free to minimize shear by setting the surface elevation to the maximum gap at all locations, resulting in the flat surface shown in Fig. 12. Near r = 0 the surface elevation does increase slightly.This is because the sensitivity of the normalized apparent viscosity with respect to gap height is negligible in this region.We found that a completely flat surface enhances the objective value by only 0.016% compared to the optimization solution.(c) p FIGURE 12: Optimal texturing of an arbitrary but symmetric surface approaches a flat plate.Symmetric texturing is unable to produce normal force.The optimal solution maximizes the gap everywhere to decrease the viscous shear load.
When the symmetry constraint is removed, the sliding surface textures can generate non-zero normal forces as with the asymmetric cylindrical texture designs.Figure 13 shows asymmetric texture geometries and the corresponding pressure fields for six different Pareto-optimal designs.Designs are displayed in order of increasing normalized apparent viscosity value (and also increasing normal force).The designs shown in Fig. 13 are sampled uniformly from the Pareto set obtained by solving Eqn.(6a) with splines N = 5.The corresponding Pareto set is shown as blue squares (ref.case) in Fig. 14, which is observed to be non-convex in the objective-function space.Using an approach such as the ε-constraint method was essential for resolving non-convex portions of the Pareto set.
In some portions of the design domain, the texture slope limitation has a significant impact on results.It may be helpful (e.g., for manufacturability) to add constraints to limit texture slope in further design studies.This could be implemented approximately using linear constraints on design variables.With the 30-variable spline representation, this slope constraint strategy requires 2 × (6 − 1) 2 = 50 linear constraints.All spline interpolated studies presented here use these linear constraints to limit the inclined angle of the textured profile.
The geometric result shown in Fig. 13(a-f) converged to a shape that is similar to a spiral blade.The spiral blade-like texture profile has a portion colored in yellow (higher heights, Figs.13(a-f)) that acts as a converging channel directing flow radially inward.This increases pressure near the disk center and generates a positive net normal force.This mechanism helps explain how the asymmetric spiral blade-like surface design can help increase normal force in this rotational configuration.The directed flow toward the disk center may also help reduce leakage, an important practical consideration for fluid power systems.This increased pressure at the center helps support axial loads and load applied on the sealing components, while the low hydraulic pressure near the outer rim of the disks helps to contain the fluid within the gap between sliding disks.These design results match the known previous experimental observation that a spiral groove enhances load-capacity for a parallel flat surface bearing configuration [9].James T. Allison 12 MD-16-1636 Surface texture geometry, asymmetric Computed pressure fields for above surface textures (g) η a /η 0 = 0.299

Comparison of optimal designs
Pareto-optimal designs shown in Fig. 14 show how the objective functions can be improved significantly by transitioning from simple dimpled textures to arbitrary continuous texture geometries.Note that the vertical axis is F N as opposed to −F N , to be more intuitive, so points closer to the upper left are more desirable.Only normalized apparent viscosity (η a /η 0 ) changes in the symmetric cases; the optimal solution maximizes the operating gap to the greatest extent possible.For the asymmetric case, switching from dimpled cylindrical to arbitrary spline-based textures shifts the Pareto sets significantly toward the upper left, indicating roughly an order of magnitude performance improvement while enforcing the manufacturability constraint.Enhanced normal force generation is observed as more spline control points are added as design variables (increasing texture design resolution).In the region of nearly zero normal force, however, increasing the number of spline control points, (N + 1) × N, does not provide meaningful improvement in reducing normalized apparent viscosity, because the optimal surfaces are flat with maximum gap.James T. Allison 13 MD-16-1636 Figure 15 compares texture profile level sets for full disk geometries at F N = 3 N.For each disk level set, the corresponding value of normalized apparent viscosity is given.As Fig. 14 illustrated previously, increased design mesh resolution results in optimal designs with lower normalized apparent viscosity values.While the same level of normal force is generated (F N = 3 N) by these five designs, the blade shape becomes sharper and the landscape profile becomes more detailed as the design mesh becomes finer.Figure 16 illustrates comparisons of optimal cross-sectional depth profiles for spline-based designs across a range of texture resolutions (N = 3, 5, 7) for four different radial positions r.Among solutions in the Pareto set, the optimal design point generating the normal force of F N = 3 N is chosen for this comparison.At all three resolutions, the same trend is observed.Namely, texture profiles consist of asymmetric expansion-contraction channels.In addition, these periodic asymmetric profiles shift to the right with decreasing radial position, resulting in spiral texture geometries.While low-dimension spline representations (N = 3) produce texture geometries that can capitalize on the physical effects discussed above, increasing N results in complex waviness in designs and improved performance.Higher-resolution representations, however, increase computational expense (discussed next).

Computational expense comparison
All problems were solved using dual Intel ® Xeon ® X5650 CPUs with a total of 12 physical cores.Gradient evaluation within the interior-point algorithm [37] is computed in parallel to utilize all CPU cores.All cases converged to solutions with a function value tolerance of 10 −8 .A summary of average computation time for single scalarized problems is given in James T. Allison 14 MD-16-1636 Table 3. Computation time for the complete multi-objective problem is approximately the number of Pareto-optimal points for a given curve in Fig. 14, multiplied by the average computation time for each interior-point (IP) solution reported in Table 3.
To compare the computational cost across the range of spline design resolutions, Table 3 lists the cost per single IP run (i.e., the time to find one point in the Pareto set).The computation time required to calculate both objective functions for a single design was typically less than 0.1 seconds.This is independent of N, because the computational domain (Fig. 3(b)) is the same for all design resolutions N. The solution time for a single IP run, however, increases significantly with N. For higher-resolution cases (e.g., N = 5, 6, 7), functional performance increases as shown in Fig. 14, but this must be weighed against the increase in computational expense quantified in Table 3.An appropriate choice may depend on a number of factors, including application and whether the studies are supporting early-or late-stage design decisions.

Limitations
The design studies here are performed under certain assumptions and simplifications described in Sections 2, 3, and 4, which place inherent limitations on the results of this design study.These limitations are due to simplification of governing equations, properties of certain texture design representations, and limited design representation resolution.
First, the Reynolds equation model is valid only when the aspect ratio of the confined flow is sufficiently wide (Fig. 7), but is favorable for design optimization studies due to the associated computational efficiency.Thus, design study results are limited by the approximate simulation model, which may dominate true physical limitations.In addition, cavitation is not considered in our Reynolds equation model.Although some previous studies accounted for cavitation effects [8,9,12,27,38], it was not observed in our experiments [17].Cavitation effects could be included in future work using the same optimization framework presented here with an appropriate system model.
Second, our study did not look at the coefficient of friction µ * , which is one popular way to identify an optimal design.The coefficient of friction is defined as the shear stress divided by the normal force.This measure is equivalent to the value of our first objective η a /η 0 divided by our second objective F N .The resulting optimal design using µ * as the objective function would produce a single point on the Pareto front generated here.In more general design problems, design utility may not be represented accurately by the coefficient of friction.For example, if small F N is needed but reducing η a /η 0 is more important, this is not captured by minimizing the coefficient of friction.
Third, design studies are limited to a sector-one tenth of a full disk-with a rotational periodic boundary repeated in every φ = 2π/10 radians.Although periodic boundary conditions approximate accurately full disk behavior, the assumption of periodicity implicitly constrains the designs explored.In addition, different sector sizes (e.g., {2π/6, 2π/8, 2π/12}, etc.) may produce different results and new insights.Future studies should investigate different sector sizes, as well as elimination of the periodic sector assumption and design the full disk texture directly.Furthermore, different kinematics (e.g., linear sliding motion) other than the rotating disk configuration will have different optimal designs, since pressure cannot built up near a center-of-rotation region.Future studies should address different practical configurations.
Fourth, while performance was improved by moving to a free-form texture design representation, spline parameterization still requires at least tens of design variables.Optimization studies were performed using efficient gradient-based methods with limited multi-start to balance computational expense and the improving probability of finding global optima.More effective global methods may be impractical due to problem dimension, and because the simulation is utilized here as a black-box model.
Increasing design resolution using the current optimization strategy may be impractical, but increased resolution is required to investigate true performance limits accessible via texture design.Several limitations were observed when using the James T. Allison 15 MD-16-1636 current low-dimension texture design representation, including generation of very sharp edges or sudden non-smooth slope changes.Alternative optimization strategies that support higher-resolution design representation should be investigated.One possible strategy is to use design abstraction techniques, such as generative design algorithms [39], image warping methods [40], or other surface representation techniques [41].A direct optimization using the simulation mesh values h i j may be another promising strategy, but practical solution will require novel large-scale optimization formulations that capitalize on problem structure to manage the significantly increased problem dimension (e.g., thousands of optimization variables).
Ongoing work is exploring the possibility of using convex programming techniques to support significant increases in design resolution while maintaining computational tractability (e.g., [42]).

CONCLUSION
Full-film lubrication problems have been studied extensively in the tribology community, often with the objective of improving normal load capacity and decreasing shear load resistance.Using textured contact surfaces is known to be effective for reducing frictional energy loss.More recently, an experimental study has revealed that asymmetric depth profile texture patterns can improve performance further by generating net normal forces to carry higher loads [17].Based on these previous discoveries pertaining to full-film lubrication, this article presented an investigation of how the transition to more sophisticated texture designs can improve performance in exceptional ways, and has offered insights into the underlying mechanisms that enable these improvements.The studies presented here focused on simultaneous improvement of apparent viscosity and normal force.
Predefined texture shapes, such as cylindrical dimples, limit potential performance improvements that is achievable with more general surface texture design.A general design pattern was clearly present among the Pareto-optimal results for the asymmetric spline-based design study.The resulting designs resemble a spiral blade, and observations of the flow and pressure fields provide insight into the mechanisms leveraged by these designs to improve performance.
The ultimate objective of these and other related efforts is to realize new levels of performance and efficiency for fluid power systems, and to reveal physical and design principles that support these improvements.Ongoing efforts complement the results presented here, including development of new models and solution techniques that will enable computationallyefficient design of very high resolution texture geometries by capitalizing on problem structure.Next steps include investigation of full-film lubrication problems in configurations that are more representative of practical engineering systems, investigating texture design representations based on generative algorithms, incorporating manufacturability and cost considerations, and transitioning from Newtonian to non-Newtonian fluid system models.For example, a linear sliding motion can be another configuration to represent practical reciprocating contact surfaces [43].Different targets for normal load or sliding speed may affect optimal solutions [44].Optimization results have provided a rich set of design data from which qualitative design insights can be extracted.Future work will iteratively deepen, validate, and enhance this qualitative understanding, leading to valuable design knowledge that may impact technologies, such as bearing systems, reciprocating mechanisms (e.g., automotive engines), energy systems, and fluid power systems.The mesh for a reduced-dimension cubic spline texture shape design representation is shown on the left.This design mesh is much more coarse than the mesh on the right that is required for accurate simulation.The surface geometry defined by the cubic spline is interpolated to determine all the height values at the fine mesh points required for simulation.(c-d).In (b), the two objectives (blue: η a /η 0 , green: F N ) are plotted together as a function of β for the inclined plane design.The gray area indicates the region of suboptimal designs.FIGURE 12 Optimal texturing of an arbitrary but symmetric surface approaches a flat plate.Symmetric texturing is unable to produce normal force.The optimal solution maximizes the gap everywhere to decrease the viscous shear load.

FIGURE 1 :FIGURE 2 :
FIGURE 1: Textured surface design problem setup adapted from experiments conducted by Schuh and Ewoldt [17].The periodic asymmetric dimpled textures used in the previous study are replaced here with arbitrary continuous texture shapes.(a) Front view of experimental setup, (b) top view of one textured surface periodic sector.

FIGURE 3 :
FIGURE 3: Mesh points for (a) design representation and (b) flow simulation.The mesh for a reduced-dimension cubic spline texture shape design representation is shown on the left.This design mesh is much more coarse than the mesh on the right that is required for accurate simulation.The surface geometry defined by the cubic spline is interpolated to determine all the height values at the fine mesh points required for simulation.

FIGURE 4 :
FIGURE 4: Alternative texture shape design representations.A cross-sectional view of each representation is shown in the top row, while a top view of each textured sector is shown in the bottom row.The gray area denotes the textured region (non-gray regions indicate unchanged flat surfaces).The top flat plates are rotating in the direction of the thick blue arrows, while the bottom textured surfaces are fixed.(a) Symmetric cylindrical texture, (b) asymmetric cylindrical texture, (c) asymmetric planar texture spanning full sector area, (d) arbitrary continuous texture with symmetry constraint, (e) arbitrary continuous texture.

FIGURE 6 :
FIGURE 6: Objective function contours based on comprehensive sampling of the cylindrical texture design (symmetric and asymmetric textures).Trade-offs are apparent, and the Pareto sets are illustrated in the figures.Asymmetry improves performance capability, especially in terms of normal force.

FIGURE 7 :
FIGURE 7: Comparison of the Reynolds equation and the Navier-Stokes equation solutions in terms of the normalized apparent viscosity for the symmetric cylindrical textures as a function of depth.

FIGURE 10 :FIGURE 11 :
FIGURE 10: Geometry, shear stress, and pressure field of the surface texture with inclined plane spanning full sector area.

FIGURE 13 :
FIGURE 13: Geometry and pressure field of optimized surface textures with spline design representation (spline resolution N = 5).Asymmetry is permitted.Six designs are uniformly sampled from 27 designs in the Pareto set to illustrate trends.(a)-(f) are geometric texture profile level sets, and (g)-(l) are pressure field level sets.Note that these Pareto-optimal designs are located in Fig. 14 (blue squares).

FIGURE 14 :
FIGURE 14: Comparison of optimal designs in the objective function space for all design studies presented in this article.Performance indices are improved significantly by increasing design flexibility via the texture surface representation.

FIGURE 15 :
FIGURE 15:  Comparison of geometric texture profile level sets of optimal designs for N =3, 4, 5, 6, 7 at F N = 3 N.As design resolution increases, the thickness of the blade-like feature becomes sharper, and a lower normalized apparent viscosity value can be achieved.

FIGURE 16 :
FIGURE 16: Comparison of depth profiles: cross-sectional view of optimal designs, spline design representation cases with N =3, 5, 7, and F N = 3 N.The horizontal axis indicates angular (θ) location within a sector, and vertical axes correspond to texture elevation.Each subfigure corresponds to a different spline N. Positive thrust occurs when the opposing upper flat surface moves in the positive θ-direction.

FIGURE 1
FIGURE 1  Textured surface design problem setup adapted from experiments conducted by Schuh and Ewoldt[17].The periodic asymmetric dimpled textures used in the previous study are replaced here with arbitrary continuous texture shapes.(a) Front view of experimental setup, (b) top view of one textured surface periodic sector.

FIGURE 2 FIGURE 3
FIGURE 2 An example of texture height profile contours shown for a design domain sector (a) and for a full disk (b).Darker color corresponds to a lower surface level (larger gap height)

FIGURE 4 FIGURE 5
FIGURE 4 Alternative texture shape design representations.A cross-sectional view of each representation is shown in the top row, while a top view of each textured sector is shown in the bottom row.The gray area denotes the textured region (nongray regions indicate unchanged flat surfaces).The top flat plates are rotating in the direction of the thick blue arrows, while the bottom textured surfaces are fixed.(a) Symmetric cylindrical texture, (b) asymmetric cylindrical texture, (c) asymmetric planar texture spanning full sector area, (d) arbitrary continuous texture with symmetry constraint, (e) arbitrary continuous texture.

FIGURE 6
FIGURE 6  Objective function contours based on comprehensive sampling of the cylindrical texture design (symmetric and asymmetric textures).Trade-offs are apparent, and the Pareto sets are illustrated in the figures.Asymmetry improves performance capability, especially in terms of normal force.

FIGURE 7
FIGURE 7 Comparison of the Reynolds equation and the Navier-Stokes equation solutions in terms of the normalized apparent viscosity for the symmetric cylindrical textures as a function of depth.

FIGURE 8
FIGURE 8  Geometry, shear stress, and pressure fields of optimized cylindrical texture designs in symmetric and asymmetric configurations.Symmetric texture R t = 4 mm, h = 2.5 mm.Asymmetric texture R t = 4 mm, β = 3.6°.

FIGURE 9
FIGURE 9  Geometry and pressure field for asymmetric cylindrical textures with various angles (β).Six β angles are sampled within the range of inclination limits.(a)-(f) are geometric texture profile level sets, and (g)-(l) are pressure field level sets of corresponding β angles.

FIGURE 11
FIGURE 11  Optimal designs for the cylindrical design exploration and the inclined plane spanning full sector area studies.The solid red stars in (a) is the Pareto set for the cylindrical texture identified from the contours of Figs.6(c-d).In (b), the two objectives (blue: η a /η 0 , green: F N ) are plotted together as a function of β for the inclined plane design.The gray area indicates the region of suboptimal designs.

FIGURE 13
FIGURE 13  Geometry and pressure field of optimized surface textures with spline design representation (spline resolution N = 5).Asymmetry is permitted.Six designs are uniformly sampled from 27 designs in the Pareto set to illustrate trends.(a)-(f) are geometric texture profile level sets, and (g)-(l) are pressure field level sets.Note that these Pareto-optimal designs are located in Fig.14(blue squares).

FIGURE 14
FIGURE 14  Comparison of optimal designs in the objective function space for all design studies presented in this article.Performance indices are improved significantly by increasing design flexibility via the texture surface representation.

FIGURE 15 FIGURE 16
FIGURE 15  Comparison of geometric texture profile level sets of optimal designs for N =3, 4, 5, 6, 7 at F N = 3 N.As design

Table 1 :
Cylindrical texture surface parameters

Table 2 :
Spline interpolated textured surface parameters

Table 3 :
Comparison of average computation cost per single scalarized optimization problem and system performance obtained from ε-constraint method and interior-point (IP) algorithm

TABLE 1
Cylindrical texture surface parameters TABLE 2 Spline interpolated textured surface parameters TABLE 3 Comparison of average computation cost per single scalarized optimization problem and system performance obtained from ε-constraint method and interior-point (IP) algorithm