Multi-scale topology optimization of micro-textured gear for interface enriched lubrication and anti-scuffing load-bearing capacity

Introduction

Gears have become critical foundation components in high-end marine engineering equipment due to their outstanding advantages in transmission efficiency, speed ratio range, precision, and fatigue resistance. ¹ Their design and manufacturing in-volve multiple interdisciplinary coupling issues, including precision tooth profile machining, backlash control, time-varying meshing stiffness analysis, installation error compensation, and bearing clearance matching. These aspects reflect the strongly non-linear and extremely complex nature of transmission systems, attracting continuous and in-depth research by scholars worldwide.^2–4 In the field of mechanical engineering, the meshing interface of machined gears appears macroscopically smooth and often anisotropic. However, at the microscopic scale, it is composed of numerous asperities and valleys with varying geometric dimensions, exhibiting a certain direction-al characteristic in interface texture.^5–7 The meshing behavior of gear interfaces can be regarded as the rolling, sliding, and mutual extrusion of these microscopic asperities. Differences in micro-morphology directly reshape the interfacial meshing behavior. On a rough interface, only the protruding asperities actually bear the load. Once the con-tact load exceeds a critical threshold, these asperities collapse, leading to a sharp in-crease in frictional resistance and energy dissipation at the meshing interface. By altering the local contact state, the micro-morphology of the meshing interface plays a crucial role in regulating the load-bearing and anti-scuffing performance of gears. Fractal theory overcomes the size-dependence limitation of conventional approaches by introducing interfacial roughness as a characteristic parameter, enabling a more accurate description of rough interfaces.^8–10 By means of a scale-independent fractal function, the Majumdar–Bhushan (M–B) model is able to reconstruct the morphology of a rough interface composed of micro-asperities.^11–13 Based on the M–B fractal function, a fractal interface contact model is established to correct the backlash error of micro-scale rough interfaces caused by the contact of asperity peaks. The modified backlash is then applied to the gear transmission system to analyze the influence of fractal parameters on the dynamic load-bearing response characteristics of the system.^14,15 A comparison among constant backlash, probabilistically distributed backlash, and fractal backlash demonstrates the distinct advantage of fractal backlash in characterizing the true backlash at meshing interfaces, confirming its higher accuracy and adaptability.^16–18 Based on the fractal theory of textured micro-elements, a meshing torque transmission model for convex–concave asperity rough interfaces is established, enabling the quantitative characterization of how fractal parameters govern the time-varying meshing stiffness (TVMS). The study reveals that TVMS is positively correlated with the fractal dimension but negatively correlated with the characteristic fractal scale coefficient. Among these parameters, the fractal dimension is identified as the dominant factor governing the variation in meshing stiffness.^19–21 A multi-degree-of-freedom dynamics model is established for gear pairs meshing whose contacting tooth flanks possess fractal morphologies that are not simply reducible to two-dimensional profiles. The study reveals that as the fractal dimension increases and the characteristic scale coefficient decreases, the contact stiffness at the interface tends to become more uniform, meshing impacts are suppressed, and the efficiency of vibrational energy dissipation improves, thereby enhancing the smoothness of gear trans-mission.^22,23 Research in this field has largely been based on two-dimensional fractal rough interface theory, which overlooks the morphological characteristics along the tooth width direction, a critical dimension that governs the true contact footprint. Those models fail to capture the actual contact mechanics occurring within real meshing interfaces. ²⁴ The analytical expression for the contact stiffness of a single rough interface is derived from Hertzian theory, while the TVMS between micro-asperity interfaces is obtained through the roughness distribution function. ²⁵ To address this issue, a three-dimensional fractal rough meshing interface model is developed, in which the micro-elements contact model for convex–concave asperities is derived and a distribution function characterizing the rough interface morphology is defined. This approach reveals the influence of fractal parameters on the rough interface gap and the TVMS, thereby providing a theoretical basis for subsequent studies on modifying the mechanical behavior of asperity micro-element contacts.^26,27 Literature review reveals no prior studies on such micro-element fractal models. ²⁸ Existing studies have attempted to extend the three-dimensional M–B equation to anisotropic fractal meshing interfaces by introducing random phases, aiming to characterize the directional features of rough asperity micro-elements. However, the inherently uncontrollable nature of random phases makes it impossible to generate rough profiles with specific distributions of convex–concave asperities on demand, nor can it produce regularly oriented micro-textures.^29,30 This strategy is still insufficient to support a rigorously anisotropic fractal model of micro-scale roughness in gear meshing interfaces.

Current gear meshing optimization methods face fundamental limitations: Taguchi methods struggle with continuous variables and interaction effects; ad-joint-based approaches require objective differentiability, precluding contact nonline-arities; particle swarm algorithms exhibit slow convergence and parameter sensitivity. ³¹ While the solid isotropic material with penalization (SIMP) method has im-proved boundary crispness, it remains confined to single-material distribution. ³² Meanwhile, lubricant distribution is shifting from macro-scale homogeneity to micro-textured multi-scale synergy. Yet existing research predominantly isolates macroscopic and microscopic problems, lacking a unified optimization framework. Fractal roughness theory is inherently descriptive, limited to characterization rather than proactive design. Engineering demands are transitioning from understanding rough interfaces to designing functional interfaces, where micro-textures actively modulate interfacial behavior. This imperative has catalyzed research into micro-texture topology optimization. Fractal theory offers pivotal insights: optimal interfaces require hierarchical multi-scale architectures. Single-scale textures suit only specific operating conditions, whereas biomimetic multi-scale textures accommodate broader operation-al ranges. Micro-texture topology optimization necessitates a multi-level framework spanning macro-layout, meso-geometry, and micro-roughness. Fractal theory elucidates topography-performance relationships, providing both biomimetic inspiration and mathematical formalism. This convergence marks a paradigm shift in tribology from passive adaptation to active control, charting new pathways for high-performance gear transmission systems.

Three mainstream approaches are compared: (a) empirical parametric scanning methods optimizing groove and dimple geometries via single-variable analysis; (b) surrogate model-based optimization (e.g. NSGA-II coupled with CFD) employing single-objective or weighted multi-objective formulations; and (c) the homogenization-based multiscale topology optimization method proposed herein. The comparison criteria include design degrees of freedom, computational efficiency, synergistic optimization capability for micro-texture morphology and distribution, and consideration of lubrication-load carrying coupling mechanisms. The core innovation lies in the proposed approach’s ability to achieve accurate performance prediction without presupposing texture shapes, while obtaining non-intuitive microstructure configurations through topology optimization, thereby overcoming the limitations of decoupled shape-distribution optimization inherent in conventional methods. Table 1 summarizes the distinctions between the proposed homogenization-based multiscale topology optimization framework and existing mainstream approaches. By leveraging homogenized constitutive tensors with spatiotemporal variability, this method enables real-time coupling between micro-texture configurations and macroscopic interface performance, eliminating the sequential workflow of microscopic optimization followed by macroscopic validation required by traditional methods, while simultaneously improving computational efficiency.

OPEN IN VIEWER

Table 1.

Comparative analysis of the proposed homogenization-based multiscale topology optimization method against existing mainstream approaches.

Comparison dimension	Traditional single-scale optimization	Hierarchical optimization	Proposed method
Optimization scale	Single micro- or macro-scale	Macro-micro decoupled	Multi-scale coupled synergy
Computational efficiency	Low (full-scale resolution)	Medium (iterative coupling)	High (homogenization-based dimensionality reduction)
Physical consistency	Cross-scale continuity difficult to ensure	Interface matching challenges	Guaranteed by homogenized tensors
Engineering applicability	Excessive idealized assumptions	Complex parameter transfer	Direct guidance for MET fabrication

This study proposes a homogenization-based multi-scale collaborative optimization methodology for gear MTMIs. An optimization objective function targeting macroscopic interfacial performance is formulated within the multi-scale theoretical framework, and a topology optimization strategy is employed to achieve optimal MET configuration design. The homogenization analysis characterizes the modulatory effects of MET on local flow and stress fields, while the topology optimization breaks through conventional parametric geometric constraints, enabling the generation of irregular MET topological patterns adaptive to operating conditions, thereby ensuring optimal macroscopic ASLBC.

Multiscale modeling of micro-textured meshing interfaces with micro-element texture configurations

The line contact thermo-elastohydrodynamic lubrication (TEHL) model proposed is illustrated in Figure 1. It incorporates rolling–sliding motion, micro-texture, and roughness characteristics to simulate textured interface lubrication meshing and load-bearing behavior that more closely approximates real engineering conditions. The anisotropic 3D fractal rough meshing interface is characterized by equivalent ellipsoidal METs with directional distributions. A 3D anisotropic fractal interface model representing multi-DOF nonlinear dynamics is constructed, establishing mathematical mapping between microscopic topographical parameters (fractal dimension, anisotropic stretching, asperity aspect ratios) and macroscopic dynamic response. The meshing interface possesses multi-scale micro-texture morphology characteristics. The rectangular micro-texture configuration is selected as the fundamental geometric primitive in this study, considering its (i) decoupled dimensional controllability for parametric sensitivity analysis, (ii) compatibility with precision laser manufacturing processes, and (iii) establishment of benchmark data for subsequent comparative evaluation with optimized bio-inspired textures. Rectangular equal-section grooves are adopted as the basic micro-texture units, with a geometric width of W , a depth of H , and a spacing T between adjacent grooves. The micro-groove units are arranged perpendicular to the direction of the interface meshing line.

OPEN IN VIEWER

Figure 1.

Line contact elastohydrodynamic lubrication model.

Based on the geometric characteristics of micro-textures, a priori design strategy is proposed to achieve synergistic optimization of feature scale and spatial distribution. Under the practical constraints of limited resolution and configuration freedom in existing micro-texturing technologies for meshing interfaces, this strategy offers engineering feasibility. Advanced manufacturing techniques enable the in-situ construction of high-precision, periodically reproducible complex microstructures on tooth interfaces, providing new degrees of freedom for lubrication performance regulation. By extending the hydrodynamic lubrication design of meshing interfaces to a harmonizable scale, the influence mechanism of micro-texture configuration patterns on the lubrication load-bearing coupling response at the interface is systematically characterized using homogenization theory. Furthermore, topology optimization is introduced to inversely design the micro-unit configuration pattern itself. This constitutes the core breakthrough direction of the present study, aiming to establish a cross-scale mapping relationship between micro-texture configuration and lubrication performance, thereby achieving high-precision customization of complex texture patterns and synergistic enhancement of lubrication and load-bearing performance. This perpendicular distribution configuration has been proven to significantly enhance the load-bearing capacity of the meshing interface. At the meshing point, the transient line contact of a gear pair can be equivalent to a pair of cylindrical rollers with curvature radii R 1 and R 2 , elastic moduli E 1 and E 2 , shear moduli G 1 and G 2 , and circumferential velocities ω 1 and ω 2 , as shown in Figure 2.

OPEN IN VIEWER

Figure 2.

Equivalent cylindrical roller model for the transient line contact of gear pairs.

The size of micro-texture elements typically ranges between 1.0 and 100.0 μm. In contrast, the MET on the meshing interface generally falls within the range of 0.1–1.0 μm, as illustrated in Figure 3. The macroscopic geometric dimensions of gear pairs, such as the radius of curvature at the meshing point, are measured in millimeters. As a gear contact interface problem, the macroscopic lubrication load-bearing response can be synergistically regulated through optimization of inter-face micro-textures. When appropriate micro-textured units are applied to the geo-metric features of the contact interface, the dynamic lubrication behavior, within a multi-scale mechanical framework, becomes dominated by the global microscale thin film. This thin film evolves into the classical Reynolds equation at the characteristic scale of the micro-textures.

OPEN IN VIEWER

Figure 3.

Micro-textured gear pair rough line contact elastohydrodynamic lubrication model: (a) Micro-textured gear meshing pair and (b) Micro-element convex-concave asperity roughness model.

Based on the assumption of an incompressible Newtonian fluid at the interface and ignoring unilateral constraints from cavitation contact, the thin film model discards all physical nonlinearities, providing a clean starting point for constructing a multi-scale optimization framework. This formulation simultaneously outputs both dissipative and non-dissipative information at the interface, the non-dissipative component manifests as the self-generated fluid pressure, the solution to the thin film equation, which directly contributes to the net load-bearing capacity. A multiphysics model correlating MET distribution parameters with the ASLBC of IEL is established by comprehensively considering the coupled characteristics of stress fields, film thickness distributions, and temperature fields. This study reveals two non-standard key characteristics of multiscale hydrodynamic lubrication, challenging the conventional understanding of scale-optimization theory. A tensor-based spatiotemporal evolution model is developed to characterize the homogenized response of micro-textures, elucidating the adaptive lubrication and load-bearing mechanisms induced by varying interface conditions. The work provides theoretical support for the design of high-load-capacity lubrication systems.

Multiscale collaborative optimization of micro-textured configurations for meshing interfaces

Considering the design issues of MET configurations at the macroscopic scale, a universal analytical framework is developed. This study considers the characteristics of stress, film thickness, and temperature fields, among other factors, to establish a model that correlates MET distribution parameters with the load-bearing performance of the IEL. The evolution of IEL load-bearing capacity under the influence of multi-scale MET parameters is optimized, and the coupling relationships between macroscopic parameters, such as load distribution, speed, and structural configurations, and microscopic parameters, including micro-asperity peak heights, texture orientation, and autocorrelation lengths, are explored, as shown in Figure 4. The former tracks the continuous decay of the instantaneous film thickness over time, while the latter describes the directional gradient of the film thickness h 0 along spatial coordinates. Both types of problems satisfy homogeneous Dirichlet conditions for pressure P 0 on the global boundaries. The upper interface is smooth and remains in a constant slip motion state, while the lower interface is embedded with MET and remains stationary. For medium layer film lubrication, the two interfaces are parallel to highlight the pure transient squeeze effect; for ASLBC, the interfaces may be inclined along both directions, with an appropriate vertical scaling to clearly reveal the coupling characteristics between the macroscopic geometric configuration and the MET.

OPEN IN VIEWER

Figure 4.

Interfacial hydrodynamic contact exhibits two typical macroscopic regimes: medium layer film lubrication and the ASLBC problem, both reflecting the spatiotemporal evolution of film thickness within the meshing zone.

Regarding the lubrication enrichment issue at the meshing interface, the upper interface moves only tangentially with a velocity U , and the variation in lubricant film thickness primarily stems from spatial distribution effects. Due to the inclination angle between the two meshing interfaces, the average film thickness h 0 of the lubricant medium exhibits a spatial gradient along the direction of hydrodynamic pressure in the wedge-shaped gap. For the squeeze film problem of MET lubricant medium layers at meshing interfaces, the two interfaces remain parallel and in contact. The variation in film thickness is entirely time-dependent, with the upper interface undergoing purely vertical motion at a prescribed velocity U ¯ , thereby causing the average film thickness h 0 to change over time. These two model problems correspond, respectively, to the two fundamental forms of film thickness variation by which the lubricant medium at meshing interfaces generates contact pressure and thus provides load-bearing capacity. Based on the optimization scheme to be presented in this paper, the MET design for any macroscopic problem that simultaneously involves these two basic variations can be promptly developed on this foundation.

The average film thickness h 0 of the dielectric layer film exhibits periodic variation between its maximum h 0 max and minimum h 0 min values. In the lubricating film model of the dielectric layer, the two interfaces remain parallel, and changes in film thickness manifest as purely time-dependent, resulting from the prescribed vertical motion of the upper interface at a given velocity ∂ h 0 ∂ t . In the ASLBC model, the upper interface also undergoes motion, but this motion is purely tangential sliding (set the sliding velocity to U ¯ ). Consequently, variations in the average film thickness h 0 exhibit spatial characteristics, governed by the inclination angle formed between the two interfaces. These two typical problems correspond to two distinct mechanisms for generating interfacial fluid pressure and the types of film thickness variations that affect ASLBC. Within a multi-scale framework, the design of METs based on combinations of these two fundamental variants can be efficiently achieved through the optimization solution proposed in this study. The micro-texture height variation, characterized by a zero-mean distribution designated as H ¯ (where H ¯ represents the texture-induced topographical evolution height), exhibits periodic variations within individual texture micro-elements. All METs across the meshing interface conform to this fundamental pattern. This assumption essentially endows the entire meshing interface with a unified and reproducible micro-textural morphology, rather than permitting geometry that varies randomly with macroscopic contact position. ³³ The local film thickness within METs governs the asperity response characterization at the micro-scale level, distinguishing between convex–concave asperities. This quantitative relationship can be expressed as:

h ( x , y , t ) = h 0 ( x , t ) − H ¯ ( y )

(1)

where, x represents the two-dimensional macroscopic position of the meshing interface, and y refers to a specified point in MET. It should be noted that variations in the MET of the meshing interface imply that the film thickness h becomes position-dependent, that is, h ( x , y ) . In equation (1), the x-dependence is neglected when addressing the lubricant-film problem, whereas the t-dependence is neglected when addressing the meshing-load problem. This assumption enables a classical homogenization analysis based on asymptotic expansion, through which the macroscopic influence of the MET on the meshing response can be captured. For hydrodynamic lubrication, the microscopic physics is accurately described by the Reynolds equation, under the general condition that both contacting interfaces of the meshing interface are rough and undergo sliding motion, a rigorous contact analysis can be performed. The problem discussed herein is formulated exclusively in a manner appropriate to the current category, culminating in a definitive statement of the analytical results. If Q ff denotes the macroscopic fluid flux per unit area at the meshing interface, it can be expressed as:

Q ff = − A ∇ P 0 + C U ¯

(2)

where, P 0 denotes the contact pressure generated in the fluid between the meshing interfaces, is the macroscopic solution to the lubrication load-bearing problem, which can be described as:

∇ · Q ff = − ∂ h 0 ∂ t ( in Ω )

(3)

Based on the macroscopic two-dimensional interface Ω model, homogeneous Dirichlet boundary conditions ∂ Ω are applied uniformly. ³⁴ These conditions are interpreted as the pressure relative to the ambient pressure (any exceptions will be stated explicitly). ³⁵ The macroscopic constitutive tensors { A , C } characterize the homogenized response of the micro-texture, their magnitudes depend on how the lubricant film thickness h varies with the spatial variable x and the temporal variable t , and are therefore dictated by whether the problem at hand micro-textured interface lubrication or ASLBC. In either scenario, the tensors are non-constant, which increases the computational cost of MET. For lubricated interface problems with micro-textures, a single value A applies across the entire field (in this case, U ¯ = 0 , so C is irrelevant). Since the reference film thickness h varies with time, the tensor must be recalculated at each macroscopic time integration point by solving the MET problem. In ASLBC problems, the reference film thickness h exhibits spatial variations across the interface, leading to a pointwise constant distribution of { A , C } (due to ∂ h 0 / ∂ t = 0 ). During the evaluation of the finite element weak form, it is still necessary to solve two MET problems (each defined as an independent problem corresponding to one tensor) at each macroscopic spatial integration point for recalculation. To provide a unified description of these characteristics, define μ as the fluid viscosity (set to 0 . 15 Pa · s or all cases), and the following coefficients are specified (under isoviscous Newtonian conditions, temperature is treated as a passive scalar for validation purposes only, excluding thermo-viscous coupling effects ³⁶ ):

{ a = h 3 12 μ b = 0 . 5 h

(4)

Herein, with the introduction of the gradient term:

{ Λ i j = ∂ Ω j ∂ y i λ i j = ∂ ω j ∂ y i

(5)

where, the two unknown vectors, ω and Ω , characterize the unit-mesh problem within the MET that needs to be solved, subject to periodic boundary conditions.

{ ∂ ∂ y i ( a Λ i j + b δ ij ) = 0 ∂ ∂ y i ( a λ i j + a δ i j ) = 0

(6)

where, δ ij denotes the components of the identity tensor I T . Defining the spatial average of this quantity over the MET domain as {·}, the constitutive tensor is then determined by the following relation:

{ C ij = { a Λ i j + b δ ij } A ij = { a λ i j + a δ i j }

(7)

An equivalent description is also acceptable as follows:

{ C ij = { a ( λ k i + δ jk ) ( λ k j + δ jk ) } A ij = { a ( λ k i Λ k j + Λ i j ) + b ( λ j i + δ ij ) }

(8)

Considering the above derivation, an explicit mapping between MET configuration parameters, macroscopic lubrication performance, and ASLBC can be systematically established. This provides an interpretable and quantifiable theoretical foundation for subsequent sensitivity analysis in MET design. It follows that A is a symmetric positive definite tensor, characterizing the inherent stiffness and dissipation properties of the MET. In contrast, C generally loses symmetry, reflecting the non-conservative effects arising from the coupling between slip and local curvature at the meshing interface. For an ideal homogeneous meshing interface (where the h ¯ = 0 condition holds), under given Reynolds function and MET characteristic scale ratios, rapid integrated prediction of MET geometry, lubrication performance, and ASLBC can be achieved. In this case, the macroscopic constitutive tensor reduces to an explicit Reynolds-type formulation:

{ A hs = h 0 3 I 12 μ C hs = 0 . 5 h 0 I

(9)

The coupling mechanism between the lubricant medium layer film and meshing load-bearing capacity, along with the corresponding MET parameter characterization, is expounded. The numerical solution for the homogeneous case ( h ¯ = 0 ) of the MTMI is obtained. The pressure distribution patterns under this condition are analyzed and visualized using both two-dimensional and three-dimensional interface plots. (Left image) Homogeneous interface pressure P hs for lubricant medium layer film (homogeneous interface meshing load-bearing capacity L hs = 100 N ) and (right image) homogeneous interface pressure P hs for meshing load-bearing capacity ( L hs = 500 N ).

Under the homogenized condition, the analytical results for both problems illustrated in Figure 4 are compiled in Figure 5 to serve as a reference for subsequent investigations. The homogenization approach is usually adopted to characterize the MTMIs under TEHL, with the Reynolds equation serving as the governing equation. When multi-scale MET effects are considered, the lubrication behavior of meshing interfaces becomes highly sensitive to the representative parameters.^37–39 This modeling procedure normally requires the lubricant film thickness h 0 to be vanishingly small compared with the macroscopic characteristic length L m ( h 0 / L m → 0 holds). In deriving the Reynolds equation from the three-dimensional Navier–Stokes equations, the following simplifying assumptions are conventionally introduced: the pressure variation across the film thickness is neglected, micro-scale inertial effects and other nonlinear MET terms are omitted, and the medium-layer film thickness is taken to evolve gradually.

OPEN IN VIEWER

Figure 5.

Leveraging the coupling mechanism between the lubricant medium layer film and meshing ASLBC.

Characterization of design variables in MET geometry-performance relationships

Since the non-uniform film-thickness h of the lubricant medium layer within the MET pattern triggers two coupled problems of type { A , C } whose design-variable sensitivities have already been systematically addressed in the literature, only the essentials are recapitulated here. Based on the lubrication medium layer film thickness formula (1) to achieve a unified description of the MET height, all MET units associated with macroscopic properties periodically recur on the textured interface. A normalized design variable s ∈ [ 0 , 1 ] is introduced on the meshing interface to describe the MET distribution, which is discretized for further processing. Let s be defined as the MET nodes or degrees of freedom s i on each unit. By discretizing s within a specified radius R , the geometric dimensions (length, depth, etc.) of the MET can be conveniently characterized, yielding the morphological variable ρ ∈ [ 0 , 1 ] with the degree of freedom of the MET being constant ρ I . Relying on standard topology optimization theory, the MET height H ¯ ∈ [ h ¯ min , h ¯ max ] can be expressed as:

H ¯ = h ¯ min − ρ λ ( h ¯ min − h ¯ max )

(10)

Additionally, it is stipulated that the average value across the entire domain of the textured interface micro-element within the meshing region shall be zero.

{ H ¯ } = 0

(11)

With specific design parameters set, when λ = 1 , the value of ρ will transition smoothly between 0 and 1, thereby forming a MTMI with a continuously varying height. This configuration can be referred to as a [ 0 − TO − 1 ] gradient design. On the other hand, when λ = 3 , the value of ρ exhibits a sharp jump between 0 and 1, resulting in a pronounced contrast in the height distribution of the meshing interface. This type of design may accordingly be termed a [ 0 − OR − 1 ] step design.

The morphological filtering operator F maps the design variable field s to the density field ρ = F ( s ) , characterizing the distribution morphology of the MET adapted to a given topology design. Specifically, the model employs a linear filtering in [ 0 − TO − 1 ] mode and an exponential erosion filtering in [ 0 − OR − 1 ] mode, with both filter radii measured in terms of the number of METs. For any given distribution of s , ρ = F ( s ) can be precomputed, and the local field h = h 0 − H ¯ (where H ¯ is the texture-induced evolution height) can be constructed. This allows solving the single MET problem to evaluate the macroscopic elasticity tensors { A , C } and their sensitivities with respect to the design variables { ∂ A ∂ s I , ∂ C ∂ s I } . According to the variational analysis of the unit cell problem as expressed in equation (6), these sensitivities can be directly and analytically computed through area integrals over the single MET region, eliminating the need for additional solutions of adjoint sensitivity equations. This approach enhances computational efficiency during multiscale optimization iterations. The derivation is inherently consistent with the homogenized response function defined in equation (8).

{ ∂ A ij ∂ s I = [ ∂ a ∂ s I ( λ k i + δ ik ) ( λ k j + δ jk ) ] ∂ C ij ∂ s I = [ ∂ a ∂ s I ( λ k i Λ k j + Λ i j ) + ∂ b ∂ s I ( λ j i + δ ij ) ]

(12)

To compute the required local sensitivities ∂ a ∂ s I and ∂ b ∂ s I , direct derivation and determination can be performed based on equations (4) and (10), utilizing the implicit construction form of a specific filter.

Sensitivity analysis of numerical discretization for multiscale characterization of micro-textured interfaces

In both the IEL and ASLBC macro-scale problems, the meshing interface Ω is defined as a square domain ( L m × L m ) with side length L m = 200 μ m . The solution of the macro-scale problem on this interface is implemented using a finite element mesh of N M × N M , employing exclusively bilinear quadrilateral elements with Gauss–Legendre quadrature of order N G × N G performed on each element (where N G = 2 ). In the IEL problem, micro-scale analysis must be conducted at ( N M N G ) 2 integration points. The total number of nodal degrees of freedom for this mesh will be denoted by the order of magnitude N M 2 in subsequent discussions. Conversely, in the ASLBC problem, the pressure P 0 varies smoothly along the meshing interface; adopting N M = 20 ensures adequate solution accuracy while effectively controlling the computational cost of the two-scale analysis. To solve the micro-scale problems for computing the effective coefficients { A , C } , the MET is discretized using a mesh of N M × N M bilinear elements (with N M = 40 ). The actual side length L m of an individual MET does not influence the homogenized response and may be selected arbitrarily, based on the assumption that the Reynolds equation remains valid at the micro-scale. Should this fundamental assumption fail, the three-dimensional Stokes equations are adopted at the micro-scale instead.

In the ASLBC problem, the upper interface undergoes periodic reciprocating motion, thus requiring computations to be performed over the entire meshing cycle. ⁴⁰ This cycle is discretized into N T = 10 time steps; should time integration be required, an explicit Euler scheme is employed, with all time-dependent quantities evaluated N T times accordingly. All macroscopic objective functions being formulated on the basis of the numerical solution P 0 to the macroscopic problem necessitates the computation of the sensitivity ∂ P 0 / ∂ s I of with respect to the MET description parameters s I (denoted as ∂ P 0 / ∂ s I , where s I corresponds to N m 2 degrees of freedom). This is achieved by first solving equation (3) for P 0 given a prescribed distribution of s , and differentiating P 0 with respect to s I to obtain the corresponding sensitivity problem.

( ∂ A ∂ s I ∇ P 0 + A ∇ ∂ P 0 ∂ s I − ∂ C ∂ s I U ¯ ) · ∇ = 0

(13)

On portions of the boundary ∂ Ω , the sensitivity ∂ P 0 / ∂ s I problem is governed by homogeneous Dirichlet or Neumann boundary conditions, as dictated by the corresponding conditions for equation (3). Following the convergence of the ∂ P 0 / ∂ s I iteration, direct evaluation of the macroscopic objective function sensitivity, formulated in terms of the initial parameters P 0 , becomes feasible. The computational burden inherent to such sensitivity analyses can be substantially mitigated through adjoint-based approaches. Of particular significance in steady-state IEL scenarios (i.e. time-independent formulations) is that the ASLBC at the meshing interface functions as a key macroscopic performance indicator L BC , thereby accentuating the importance of MET sensitivity data for guiding downstream texture configuration optimization.

∂ L BC ∂ s I = ∫ Ω ∂ P 0 ∂ s I d Ω

(14)

Let P 0 denote the array comprising all degrees of freedom of P 0 , K is the matrix corresponding to the A ∇ P 0 term, and K consist of all known quantities associated with the C U ¯ term. The boundary integral term is expressed as L BC = ω · P 0 , where w represents the array of weight functions associated with the numerical integration of the finite element discretization of P 0 . Denoting the degrees of freedom array of ∂ P 0 / ∂ s I by ∂ P 0 / ∂ s I , the discrete form of equation (13) is written as:

K ∂ P 0 ∂ s I + ∂ K ∂ s I P 0 = ∂ F ∂ s I

(15)

As such, equation (14) may be expressed alternatively as:

∂ L BC ∂ s I = ω · ∂ P 0 ∂ s I = ( ω K ) · [ ∂ F ∂ s I − ∂ K ∂ s I P 0 ]

(16)

The reduction in computational cost is directly manifested in the requirement to solve only a single linear system (i.e. the system regarding the adjoint variable ω / K ). However, under the current formulation of the multiscale MET optimization problem, the computational efficiency gains achieved through adjoint sensitivity analysis are not significant. Taking the ASLBC problem as an example, each MET problem of size ( N M N G ) 2 requires independent processing on a ( N m × N m ) micro-texture element mesh to construct the derivatives ∂ F / ∂ s I and ∂ K / ∂ s I . Regardless of whether the direct method or the adjoint sensitivity method is employed, this process must be repeated N m 2 times to determine all ∂ P 0 / ∂ s I , a step that dominates the overall computational cost. Consequently, the approach of directly solving equation (13) to obtain ∂ P 0 / ∂ s I is currently preferred. The advantage of this method lies in its ability to directly compute sensitivities for arbitrary macroscopic objective functions without constructing specific adjoint problems tailored to particular target functions.

A progressive meshing strategy is adopted, with a baseline element size of 0.1 mm specified for the bulk fluid domain, while local refinement is applied to critical regions such as gear meshing contact zones, reducing the element size to 0.02 mm. Hexahedral structured grids are employed to better capture boundary layer flow characteristics. Five prismatic layers are arranged adjacent to the gear meshing interfaces, with a first-layer thickness of 0.001 mm and a growth rate of 1.2, satisfying wall-resolved requirements. This approach ensures an optimal balance between computational accuracy and efficiency.

Homogenization-based multiscale topology optimization for micro-textured interfaces

Accounting for global operating condition constraints, the objective function φ is defined as a functional of the design variable s with respect to the MET distribution over the MTMI.^41,42 The optimization aims to minimize this functional subject to the constraint condition χ ( s ) = { H ¯ ( s ) } = 0 . Accordingly, the standard mathematical formulation of the MET configuration design optimization problem is expressed as:

Minimize the objective function φ ( s ) , subject to the equality constraint χ ( s ) = 0 , and the bound constraint s ∈ [ 0 , 1 ]

(17)

The objective function value ϕ and its sensitivity with respect to the design variables s , denoted as ∂ φ / ∂ s I , are employed to iteratively update the initial distribution of s , and the iterative process proceeds until the objective function value stabilizes to a practically constant level. In most cases, convergence is achieved within 1000 iterations. This study aims to construct an optimal MET configuration that yields target output values { A # , C # } closely approximating the prescribed performance metrics. To this end, the following objective function is formulated:

φ ( s ) = 1 2 ( [ | A ( s ) − A # | ] 2 [ | A # | ] 2 + [ | C ( s ) − C # | ] 2 [ | C # | ] 2 )

(18)

Figure 6 illustrates the considered MET configuration design optimization problem and its associated objective functions. If the target values are determined based on the known homogenized MET configuration of the MTMI, implying physical realizability, the process essentially equates to MET reconstruction. Conversely, should the target values be prescribed independently of the MET configuration inherent to the meshing interface, the pursuit of approximating the equivalent constitutive tensor C # with pronounced asymmetric characteristics transforms the analysis of MET parametric properties into an evaluation methodology for designing micro-textured meshing interfaces with unconventional macroscopic properties. This methodology aims to establish a synergistic correlation mechanism between IEL and ASLBC. To address the multi-physics coupling design challenges of micro-textured meshing interfaces, a multi-scale topology optimization model for MET is established based on homogenization theory. By employing macroscopic–microscopic decoupling strategies, equivalent performance parameters of micro-textures are embedded into the macroscopic optimization framework, constructing a cross-scale collaborative optimization method that maximizes ASLBC subject to IEL constraints. This approach achieves concurrent optimization of MET configurations and distribution parameters, providing theoretical foundations for homogenized micro-texture design of meshing interfaces.

OPEN IN VIEWER

Figure 6.

The MET configuration optimization problem and its corresponding objective functions.

A dual coupled convergence criterion based on relative residuals and variable increments is developed for the numerical solution of nonlinear equation systems. Considering the mesh discretization scales and physical field magnitudes, the convergence thresholds are set to 10 − 6 for relative residuals and 10 − 8 for variable increments. These tolerance settings are intentionally specified as slightly stricter than the inherent discretization error levels of the computational domain, achieving a rational balance between numerical accuracy and computational efficiency.

The characterization and analysis of the MET are predicated upon the instantaneous configuration of the IEL medium film.^43–45 For MET characterization, the local film thickness variation is quantified by geometric parameters { h ¯ min , h ¯ max } . A benchmark solution is first established for the MET configuration under homogeneous interface conditions, yielding the maximum homogeneous pressure P hom max as the reference baseline. The optimal MET design is then pursued to achieve the desired macroscopic pressure distribution P 0 # ( x ) within the computational domain Ω , enabling precise matching of complex pressure distributions. This is accomplished by deploying spatially varying MET configurations h ¯ ( x , y ) across the meshing interface to generate physically realizable inputs for forward analysis. A target binary (0–1 type) MET configuration is prescribed on the meshing interface, and direct numerical simulation is performed to compute the corresponding macroscopic pressure distribution P 0 # ( x ) . The pronounced discrepancy between P 0 # , max and P hom max substantiates the enhancement effect of the MTMI on the IEL contact ASLBC. Subsequently, an IEL-correlated ASLBC objective function for the meshing interface is constructed, establishing explicit functional relationships between IEL and the critical parameters of ASLBC. This achieves the mathematical characterization of multi-scale lubrication and load-bearing problems for textured interfaces, yielding:

φ ( s ) = ∫ Ω ( P 0 ( s ) − P 0 # ) 2 d Ω ∫ Ω ( P 0 # ) 2 d Ω

(19)

The objective function is driven to zero to obtain a micro-texture design that yields an identical pressure distribution to the target micro-texture.^46–48 This MET problem formulation also enables visual assessment of MET performance within a multiscale analysis framework, specifically regarding its capability to reconstruct a given micro-texture. It is noted that the maximum pressure P 0 max closely approximates P 0 # , max , with the characteristic dimensions and directionality of the micro-texture design bearing strong resemblance to those of the target micro-texture. Nevertheless, perfect agreement between the two is not achieved, attributable to the inherent non-uniqueness in microstructural design. This non-uniqueness becomes more pronounced in the 0–TO–1 configuration, where the expanded degrees of freedom in interface topology lead to more significant deviations between the target and designed micro-textures, despite the exact matching of P 0 # , max and P 0 max . The influence of the 0–TO–1 micro-texture, quantified by the discrepancy between P hom max and P 0 max , is less pronounced than that of the 0–OR–1 micro-texture, owing to the smoother topological variations of the former. Consequently, it can be anticipated that if the optimization objective is to maximize the load-bearing capacity of the interface, the 0–OR–1 configuration represents the superior choice under the specified design constraints (currently defined by { h ¯ min , h ¯ max } and the constraint { H ¯ } = 0 ). An alternative approach to circumvent the challenge of generating realizable target pressure distributions P 0 for optimization purposes is to adopt an indirect representation of the target distribution. In the MET problem, the same temporal instant is considered. However, rather than matching a local quantity of interest such as pressure, the objective shifts to matching a global performance metric, the load-bearing capacity. The target load capacity L 0 is prescribed as a scaling factor F L of the load capacity associated with the homogeneous interface L hom . In this case, L 0 # = L hom F L holds, and the associated objective function, which is to be driven to zero, is described as:

φ ( s ) = ( L 0 ( s ) − L 0 # ) 2 ( L 0 # ) 2

(20)

To intuitively characterize the multi-scale MET configuration variations, representative profiles are extracted along selected cross-sections for comparative analysis and visualization. The geometric parameters of the MET are specified as { h ¯ min , h ¯ max } = { − 0 . 5 , 0 . 5 } (μm), with the results presented in Figure 7. For 0–OR–1 (binary/stepwise) designs, the design objective is to achieve a sharp transition of the MET height H ¯ between its maximum and minimum bounds. Under lower load conditions with a load factor of F L = 1 . 2 , the actual MET height H ¯ exhibits a gradual transition characteristic. This phenomenon arises from the necessity to simultaneously satisfy the constraint condition { H ¯ } = 0 ; maintaining a sharp transition would inevitably result in a significant increase in the load factor F L . When F L is elevated to twice the load-bearing capacity of the homogeneous interface, the transition feature of MET height H ¯ becomes pronounced, thereby validating the stepwise characteristic of this MET configuration design under higher loads.

OPEN IN VIEWER

Figure 7.

Comparative analysis of multi-scale MET configuration design: (a) Horizontal direction position 0-OR-1 MET configuration design and (b) Horizontal direction position 0-TO-1 MET configuration design.

Meanwhile, the design objective of the 0–TO–1 (gradient/transition) is to generate smooth texture micro-elements with continuously varying MET height H ¯ . However, as anticipated, substantial load-bearing capacity can only be achieved when the transitional characteristics of MET height H ¯ become sufficiently sharp. With increasing F L , the interfacial topologies of both MET configurations exhibit convergent evolution: isolated METs progressively deepen with corresponding thickening of their sidewalls. How to achieve MET configuration design with target load-bearing capacity represents a critical practical challenge in the field of heavy-duty gear transmission systems. In numerous mechanical engineering applications, the pursuit of enhanced IEL performance and improved meshing ASLBC constitutes the core objective of MET topology design. It is particularly noteworthy that macro-scale MET dimensional parameters typically exhibit pronounced time-varying characteristics; therefore, integrating these two critical dimensions (IEL and ASLBC) into a unified framework warrants systematic elaboration.

Spatiotemporal coupling optimization of MET for enhanced meshing load-bearing capacity

This study investigates the issue of maximizing load-bearing capacity under the spatiotemporal variation of MET parameters. Two distinct operating conditions are imposed on the interface rolling–sliding velocity ∂ h 0 / ∂ t ^49–51:

For both operating conditions, the rolling–sliding motion at the meshing interface initiates from the maximum film thickness h o max , with the film thickness varying within the range { h o min , h o max } = { 1 . 0 , 1 . 25 } μm. To facilitate comparative analysis, the amplitude of film thickness fluctuation is set to { h o min , h o max } = { − 0 . 5 , 0 . 5 } ( in μ m ) . Numerical calculations are performed only for half a motion cycle, as illustrated in Figure 8. A total of 100 time steps are employed to generate the curves, ensuring sufficient resolution for capturing the meshing load-bearing variations. It should be noted that the maximum load for the continuous (sinusoidal) curve occurs prior to T / 2 , as the velocity decreases concurrently with the gap reduction. In contrast, the velocity remains constant for the discrete (square-wave) curve. For both cases, the lines serve as cross-validation: in the discrete-continuous (square-wave and sinusoidal) scenarios, sinusoidal (or square-wave) indicates that the METs originally obtained under the sinusoidal (or square-wave) curves have been applied to calculate the meshing load-bearing variation. An integral optimization strategy spanning the entire meshing cycle is employed to achieve globally optimal control of MET temporal characteristics, maximizing the load-bearing capacity of micro-textured interfaces.

OPEN IN VIEWER

Figure 8.

Evolution of load-bearing capacity in MTMIs throughout the IEL cycle: multi-scale MET configuration design and optimization for dual-velocity regimes: (a) Discrete stepwise velocity curves and (b) Continuous transient velocity curves.

The Reynolds cavitation boundary condition is enforced, wherein the cavitation pressure is prescribed as either ambient or saturated vapor pressure. The numerical implementation employs a formulation based on the penalty function method.^18,52,53

φ ( s ) = ∫ T ∫ Ω ( P 0 ( s ) ) 2 d Ω dt ∫ T ∫ Ω ( P Re f ) 2 d Ω dt

(21)

To transform the optimization problem into a tractable minimization form and accelerate convergence, a negative sign is introduced. For normalization and scaling purposes, the reference pressure is set as P Re f = P hom . The constraint χ ( s ) = { H ¯ } = 0 holds universally and remains valid throughout the entire meshing time domain Ω . Consequently, both φ and χ are functions solely of the time-independent design variable distribution. Employing an explicit Euler scheme for temporal integration, the optimization algorithm requires, at each macro-iteration, the sequential solution of the macroscopic problem and sensitivity analysis across N T time steps for a given design variable distribution of s .

The optimization results are summarized in Figure 7. Interestingly, although the macro-scale velocity field distributions differ significantly, the obtained MET configuration designs are nearly identical. This indicates that for this problem, the optimal MET is relatively insensitive to macroscopic operating conditions (velocity fields). The cross-validation results in Figure 8 further support this conclusion: the MET configuration optimized under one velocity field distribution still exhibits nearly overlapping load–time curves under another velocity condition. Even under the 0–TO–1 design constraint, the optimization results show distinct geometric sharpening features, which are crucial for maximizing load-bearing capacity. Based on the aforementioned MET configurations, the time-varying load history over a half-period is calculated using a multi-scale analysis method. The varying curves in Figure 8 demonstrate that micro-textured interfaces significantly enhance the meshing load-bearing performance. It should be noted that although the 0–TO–1 design exhibits sharp features, the introduction of a linear filter limits the degree of sharpening, resulting in its instantaneous load-bearing capacity at any meshing moment being consistently lower than that of the 0–OR–1 design. The degree of optimization approximation for each MET design can be evaluated by calculating the load factor. Figure 7 presents the load coefficients (load factors) corresponding to each MET configuration design. These values reveal the physically achievable upper limit for the target value L 0 that can be set in the optimization problem. By varying solely, the velocity parameter U ¯ under default operating conditions, the maximization of meshing ASLBC is achieved for both 0–OR–1 and 0–TO–1 MET configuration designs, with comparisons made against the homogeneous interface baseline. The maximum meshing ASLBC is realized through interface micro-dimples (rather than micro-protrusions), as micro-dimples exhibit superior modulation capability over macroscopic responses under identical constraints. For the 0–OR–1 configuration design in this case study, the continuation optimization method is exclusively employed, the penalization factor λ in equation (10) is progressively increased from 1 to 3 to eliminate intermediate density regions of the design variable s within the [0, 1] interval. For the 0–TO–1 configuration design, the MET orientation exhibits a slight yet continuous adjustment trend to ensure the attainment of maximum meshing ASLBC. Increasing the number of time steps N T , beyond the default value of 10 would incur higher computational costs. However, numerical experiments demonstrate that this default setting already yields sufficiently accurate optimization results. As illustrated in Figure 9 (where, for simplicity, only the 0–OR–1 type design is presented), even N T = 5 proves adequate for achieving satisfactory precision. This observation is attributed to the fact that micro-textures optimized for a single time instant (i.e. N T = 1 ) exhibit characteristics comparable to the summarized in Figure 7.

OPEN IN VIEWER

Figure 9.

The temporal domain is partitioned into 100 time steps to accurately resolve the transient characteristics of meshing load-bearing distributions: (a) Discrete stepwise velocity curves (0-OR-1) and (b) Continuous transient velocity curves (0-OT-1).

Nevertheless, retaining the time integration in the objective function with a moderate N T remains essential, as it eliminates the arbitrariness associated with selecting a specific time point and mitigates potential bias from an individual temporal snapshot. Given that micro-texture performance demonstrates low sensitivity to specific velocity profiles, the default value of N T = 10 is expected to remain applicable across diverse velocity conditions, thereby providing a reliable cost estimate for solving such transient problems. Even at the currently employed macroscopic spatial resolution, the computational cost associated with optimizing MET geometry configurations exceeds that of the present framework by more than two orders of magnitude. This indicates that while temporal and spatial variations in the homogenized response can be addressed within a conceptually similar algorithmic framework, the latter presents significantly greater numerical challenges, a direction that warrants focused investigation in future work. Furthermore, the geometric parameters governing MET characteristics remain unexplored in the current study and will be systematically examined in subsequent research.

Multi-scale geometric configuration and spatial distribution of MET in modulating the ASLBC

For the IEL problem with MET geometry and the ASLBC in meshing interfaces, the spatial variations of homogenized coefficients { A , C } must be considered as functions of the average medium film thickness h 0 under prevailing operating conditions. Under the default homogeneous boundary conditions (sliding velocity U ¯ = 1 . 0 m/s, nominal medium film thickness range { h 0 min , h 0 max } = { 0 . 2 , 1 . 0 } μm), the multiscale design yields the maximum homogenized pressure P hom max and the corresponding load-bearing capacity L hom . The influence of MET geometric parameters and distribution characteristics on the load-bearing performance of textured meshing interfaces can be quantitatively evaluated through two dimensionless indicators: the pressure factor F P and the load factor F L ⁵⁴ :

{ F P = P 0 max P hom max F L = L 0 L hom

(22)

The MET geometric parameters are specified as { h ¯ min , h ¯ max } = { − 0 . 2 , 0 . 1 } (μm). For subsequent multi-scale optimization analysis of MET parameters, the maximization of meshing ASLBC is designated as the objective.

φ ( s ) = − ∫ Ω ( P 0 ( s ) ) 2 d Ω ∫ Ω ( P Re f ) 2 d Ω

(23)

where, P Re f denotes the pressure distribution obtained at the film thickness of h = h 0 − h ¯ min , satisfying the condition P Re f > P 0 . The filter radius is set to R = 8 for the 0–OR–1 scheme, whereas it is maintained at the default value of R = 4 for the 0–TO–1 scheme.

Analyzing the temporal and spatial variations of the initial film thickness h 0 at the meshing interface, along with the corresponding changes in macroscopic interface variables, incurs substantial computational costs. The incorporation of spatial variability analysis necessitates solving numerous microscopic problems and performing sensitivity calculations, making it difficult to rapidly assess the correlation between various problem parameters and optimal MET configurations. Once the MET configuration is specified, the film thickness H ¯ in equation (1) becomes a fixed value, leaving h 0 as the sole remaining parameter governing the macroscopic interfacial response { A , C } . This parameter may vary across complex global patterns in either the temporal or spatial domain, that is, throughout the entire meshing cycle or across the whole meshing interface. To circumvent this difficulty, a Taylor expansion approach will be employed. Global variation patterns are irrelevant, as the determination of the macroscopic response tensors { A , C } for a given MET configuration requires only the local value of the initial film thickness h 0 , specifically, at a particular meshing instant or meshing position. Consequently, although performing a Taylor expansion of a general macroscopic function G with respect to time or space would be challenging and highly problem-dependent, expanding with respect to h 0 about a reference value h 0 Re f is straightforward and, in principle, applicable across diverse problem settings. The n th-order Taylor expansion is then expressed as ⁵⁵ :

G ( h 0 ) = G ( h 0 Re f ) + ∑ i = 1 n ( h 0 − h Re f ) i i ! ∂ i G ∂ h 0 i | h 0 Re f

(24)

Furthermore, to maintain compact notation, let γ TI = ∂ G ∂ s I ; it then follows that the sensitivity of the Taylor expansion is identically given by expanding the sensitivity:

γ TI ( h 0 ) = γ TI ( h 0 Re f ) + ∑ i = 1 n ( h 0 − h Re f ) i i ! ∂ i γ TI ∂ h 0 i | h 0 Re f

(25)

Numerical simulation analysis identifies n = 4 as the optimal choice that balances computational accuracy and efficiency, thereby ensuring the fidelity of the homogenized meshing interface mathematical model with a minimum requirement of n ≥ 3 . Since analytical expressions for the Taylor expansion derivatives of order are intractable, a fourth-order accurate finite difference scheme is adopted for numerical approximation. This scheme necessitates function evaluations at the reference film thickness h 0 Re f and six symmetrically distributed offset points with displacements of { ± δ h 0 , ± 2 δ h 0 , ± 3 δ h 0 } , where δ h 0 = 0 . 01 μm, resulting in a total of D = 8 computations. These calculations involve the evaluation of { A , C } and their partial derivatives with respect to s I , namely { ∂ A ∂ s I , ∂ C ∂ s I } . For time-dependent problems, the computational cost reduction factor for ASLBC solution and MET sensitivity analysis is N T / D . The cost of macroscopic sensitivity analysis remains invariant under both time-dependent and spatially-varying conditions. For spatially-varying problems, the cost reduction factor for ASLBC solution with MET-engaged meshing increases to ( G N M ) 2 / D , thereby enhancing the efficiency of parameter space exploration.

The meshing ASLBC problem with MET configuration is revisited to determine the expansion point, employing the 0–OR–1 microtexture as a representative case where the tensors { A , C } exhibit significant variation across the interface. The Taylor expansion approximations of these tensors are subsequently compared against the exact homogenized response. The results in Figure 10 indicate that, for this MET configuration, h 0 Re f should be selected far from either extreme to obtain a P 0 distribution that closely approximates the exact homogenized pressure distribution. The 0–TO–1 design results are almost identical to those obtained previously, whereas the 0–OR–1 interface topologies now lie closer to the 0–TO–1 counterparts. The good qualitative and quantitative agreement demonstrates that this approximation method serves as an effective tool for rapidly exploring the design space. Further MET parametric analyses (not shown), encompassing various microtexture parameters and IEL problems, confirming that this choice remains valid across most cases. Wherein expansion points of h 0 Re f ∈ [ h 0 min , h 0 max ] = [ 0 . 2 , 1 . 0 ] μm are employed in conjunction with exact homogenized meshing interface responses, with particular attention to the temporal evolution of contact pressure characteristics along the interface midline.

OPEN IN VIEWER

Figure 10.

Comparative analysis is performed between Taylor expansion-based predictions of the initial parameters P 0 .

The MET configuration follows the 0–OR–1 optimization outcome. To achieve higher accuracy in MET design evaluation, an exact homogenization error-free treatment is employed under specified operating conditions. Taylor expansion approximation methods are utilized to quantify the influence of MET parameters previously treated as constants, aiming to maximize the meshing ASLBC of the micro-textured interface for each case. Investigations reveal that the performance superiority of micro-textured meshing interfaces is largely independent of the specific film thickness distribution. Across various film thickness regimes (results not shown), the near-circular cavity configuration consistently emerges as the optimal MET topology throughout the meshing time domain. This finding offers significant computational advantages, as it eliminates the need for complex topology optimization procedures and enables the adoption of parametric MET configuration schemes with reduced degrees of freedom.

Although previous studies indicate that most macroscopic problem parameters exert negligible influence on optimal MET configuration design, the imposition of different boundary conditions in solving equation (3) yields markedly distinct outcomes. The homogeneous Dirichlet condition of zero pressure at domain boundaries has been adopted as the default in all preceding cases. When this is replaced by setting the normal component of the interfacial fluid flux (equation (2)) to zero at the boundary edges, the resulting MET configurations exhibit pronouncedly elongated characteristics. When employing zero-flux boundary conditions, the optimal texture configuration exhibits an asymmetric distribution with a shift of approximately 12%–15% toward the outlet edge. This offset effectively suppresses cavitation zone expansion through enhanced thermal effects. Under zero-Dirichlet boundary conditions, the texture distribution appears relatively symmetric, it may induce non-physical negative pressure phenomena in the outlet region. ⁵⁶ The influence of the MET morphology parameter set { h ¯ min , h ¯ max } is examined, wherein both minimum and maximum height values increase simultaneously. Constrained by the condition { H ¯ } = 0 across all designs, this translation necessitates cavity widening. The consequent substantial enhancement of pressure coefficients and load-bearing coefficients further underscores the impact of optimal MET configurations on the macroscopic meshing ASLBC response of lubricated interfaces. A discrete-continuum multiscale approach is developed to bridge micro-scale MET topology (discrete 0–1 design) with macro-scale IEL performance (continuous field analysis). Dynamic IEL load-bearing models are established under discrete stepwise velocity profiles (characterizing start-stop or gear-shift events) and continuous transient velocity curves (representing smooth acceleration/deceleration), respectively, to quantify the differential impact of velocity jumps versus gradual ramps on the meshing ASLBC.

Optimized multiple linear regression algorithm with parameter correlation

The D–W statistics for Data Model 1 (with meshing load-bearing capacity as the dependent variable) and Data Model 2 (with interface friction coefficient as the dependent variable) are 1.902 and 1.243, respectively, indicating the absence of autocorrelation in both linear regression models. The results presented in Table 2 are derived from the Durbin–Watson (D–W) test analysis. Further analysis reveals that MET parameters exhibit a relatively weak influence on the meshing load-bearing capacity, whereas their impact on the interface friction coefficient is considerably more pronounced. The coefficients of determination (R²) for Data Model 1 and Data Model 2 are approximately 0.069 and 0.531, respectively, with approximate standard deviations of 3372.61 and 0.0054.

OPEN IN VIEWER

Table 2.

Results of Durbin–Watson autocorrelation test.

Data model	R 2	R	Adjusted (R2)	Deviation	D–W test
Model 1	0.069	0.252a	0.067	3372.61	1.902
Model 2	0.531	0.658a	0.529	0.0054	1.243

The multi-dimensional dataset encompassing MET parameter combinations, meshing load-bearing capacity, and interface friction coefficients is imported into the Scientific Platform Serving for Statistics Professional (SPSSPRO) analytics platform, wherein the Technique for Order Preference by Similarity to Ideal Solution (TOPSIS) method is employed for multi-objective decision analysis. The meshing load-bearing capacity is designated as a positive indicator (benefit-type criterion), while the interface friction coefficient is designated as a negative indicator (cost-type criterion). The weights for each evaluation index are objectively determined through the entropy weight method. Based on the aforementioned analytical outcomes, the regression coefficients for the multiple linear regression (MLR) equations are obtained, as presented in Table 3. The mathematical expressions for Data Model 1 and Data Model 2 is formulated as follows:

{ X 1 = 1 . 138 × 10 5 − 30 . 40 Y 1 + 3 . 586 Y 2 − 2 . 274 × 10 3 Y 3 + 0 . 357 × 10 3 Y 4 X 2 = 0 . 103 − 3 . 302 × 10 − 4 Y 1 − 1 . 458 × 10 − 5 Y 2 − 1 . 210 × 10 − 2 Y 3 − 1 . 031 × 10 − 3 Y 4

(26)

where X 1 and X 2 represent the gear meshing load-bearing capacity and interface friction coefficient characteristics, respectively, whereas Y 1 , Y 2 , Y 3 , and Y 4 correspond sequentially to the MET dimension, area ratio, depth, and geometric configuration.

OPEN IN VIEWER

Table 3.

Parameter estimation in multiple linear regression model.

		USC		SC			CS
Data model		Coefficient (B)	SD	Delta (δ)	t-Statistic	Significance (S)	Tolerances (T)	VIF
Model-1	Constant term	1.138 × 105	352.106	–	289.739	0.000	–	1.005
	MET depth	-30.40	32.438	-0.027	-0.898	0.367	0.993	1.000
	MET size	3.568	0.573	0.197	6.985	0.000	1.000	1.000
	Area ratio	2.274 × 103	647.658	-0.092	-3.416	0.001	1.000	1.005
	MET geometry	0.357 × 103	76.752	0.096	3.548	0.001	0.993
Model-2	Constant term	0.103	0.000	–	229.483	0.000		1.005
	MET depth	-3.302 × 10−4	0.000	-0.193	-7.918	0.000	0.993	1.000
	MET size	-1.458 × 10−5	0.000	-0.624	-25.289	0.000	1.000	1.000
	Area ratio	-1.210 × 10−2	0.001	-0.298	-13.891	0.000	1.000	1.005
	MET geometry	-1.031 × 10−3	0.000	-0.161	-6.864	0.000	0.993	1.005

CS: covariance statistic; SC: standardized coefficient; SD: standard deviation; USC: unstandardized coefficient; VIF: variance inflation factor.

The significance value of MET depth in Data Model 1 is 0.367 (greater than 0.05), indicating that MET depth is not a statistically significant factor, as presented in Table 3. Considering the parameter importance across Data Model 1 and Data Model 2, the sensitivity ranking of meshing load-bearing capacity to MET parameters follows the order: dimension, area ratio, geometric configuration, and depth. Similarly, the sensitivity ranking of interface friction coefficient to MET parameters is dimension, area ratio, depth, and geometric configuration. The parameter sensitivity of interface friction coefficient obtained in this study aligns with the findings, demonstrating that MET depth exerts a less pronounced influence on interface friction coefficient compared to area ratio. This ranking differs from the results derived through the sequential factor analysis (SFA) standard deviation method, which may be attributed to the limited dataset employed in SFA. Following TOPSIS optimization, the optimal meshing load-bearing capacity and interface friction performance are achieved under the following conditions: characteristic dimension of 500 µm, area ratio of 40%, MET depth of 5.0 µm, and geometric configuration of transverse slit MET (TSl-MET). The optimized performance metrics are meshing load-bearing capacity of 1.106 MPa and interface friction coefficient of 0.0679. Compared with the untextured interface, the friction coefficient is reduced by 13.7%. Furthermore, relative to the results obtained via the SFA method, the friction coefficient exhibited an additional reduction of 4.58%.

Discussion of the results

The synergistic macroscopic response of IEL and ASLBC can be directionally modulated through micro-texturing. Optimizing the MET topology with the maximization of ASLBC at the meshing interface as the objective parameter holds significant engineering value. Existing research predominantly focuses on meso-scale MET morphological design, necessitating the development of full-scale numerical analysis methods for contact pressure distribution tailored to macroscopic MTMIs. However, as the characteristic scale of MET configurations diminishes, this analytical approach faces severe computational efficiency challenges. Current studies predominantly conduct topology design based on preset interface MET configurations, with geometric features largely confined to simple morphologies such as circular/square dimples and cylindrical protrusions, which inadequately exploit the potential of micro-textured IEL in modulating meshing ASLBC.

Qualitative comparisons with published experiments validate the multiscale numerical model. (1) Nonlinear saturation of MET depth effects: The predicted ASLBC saturation trend agrees with MTMI data in Wang et al. ²⁰ Their tribological tests showed friction reduction plateaus when dimple depth exceeds ∼3–5 μm, aligning with the model-predicted diminishing returns in load-bearing capacity. Excessive depth increases flow resistance and attenuates hydrodynamic effects, yielding marginal ASLBC gains. (2) Optimal area density: The model predicts optimal MET area density of 15%–20%, consistent with the ∼18% reported for gear meshing load-bearing interface in Wang et al. ³¹ Despite differing metrics (ASLBC maximization vs friction minimization), both studies converge to similar ranges, suggesting a universal optimum for textured interface lubrication. Trend consistency across key parameters validates the model’s physical rationality and reliability. This agreement demonstrates effective capture of core lubrication mechanisms, providing credible guidance for experimental design and engineering applications.

This paper presents an anisotropic three-dimensional fractal rough contact model of discrete micro-elements that more closely approximates real gear meshing interfaces, and unveils the underlying mechanisms governing their dynamic load-bearing characteristics. The study further uncovers the periodic evolution law of gear transmission errors induced by the fractal characteristics of interface micro-elements. By integrating textured micro-elements into the dynamic system of the meshing interface, the accuracy of parameter identification is improved, and the influence mechanism of anisotropic interface roughness on the dynamic behavior of textured micro-elements is clarified. This work provides a quantifiable theoretical basis and a feasible technical pathway for optimizing precision machining parameters of meshing interfaces.

Although the gradient-based optimization algorithm employed cannot strictly guarantee local optimality, large-scale numerical simulations and cross-validation studies demonstrate that this method can stably converge to MET configurations with significant performance improvements, and sensitivity analyses of key design parameters further confirm the engineering applicability of the optimization results. Future research could be extended toward the coupling optimization of texture orientation and operating conditions, with the objective of establishing more generalized and universally applicable design guidelines.

Conclusions

The proposed MET design framework opens new avenues for computational MTMI engineering, offering both numerical efficiency and topological flexibility. However, further investigation is required to fully evaluate and expand its application potential.

Despite the significant improvements in macroscopic ASLBC achieved by the topologically optimized non-regular MET patterns, their practical manufacturability merits consideration. Although the generated configurations depart from conventional parametric geometries, the resulting topological contours are predominantly characterized by smooth, gradually varying profiles rather than high-frequency, discontinuous features. This inherent characteristic renders the optimized designs highly compatible with advanced ultra-precision surface texturing techniques, particularly laser micromachining with focal spot modulation or mask-projection micro-ablation. Furthermore, homogenization analysis reveals that the optimized stress fields alleviate localized stress concentrations at texture boundaries, thereby reducing the risk of fatigue spalling during service. Consequently, the proposed methodology not only delivers superior tribological performance in the numerical domain but also yields a tangible, processable geometry readily translatable to physical gear surface engineering applications.