Abstract
Additive manufacturing technologies are a key point of the current era of Industry 4.0, promoting the production of mechanical components via the addition of subsequent layers of material. Then, they may be also used to produce surfaces tailored to achieve a desired mechanical contact response. In this work, we develop a method to prototype profiles optimizing a suitable tradeoff between two different target mechanical responses. The mechanical design problem is solved relying on both physical assumptions and optimization methods. An algorithm is proposed, exploiting an analogy between genetics and the multiscale characterization of roughness, where various lengthscales are described in terms of rough profiles, named chromosomes. Finally, the proposed algorithm is tested on a representative example, and the topological and spectral features of roughness of the optimized profiles are discussed.
1. Introduction
The role of the interface between different material constituents/phases is a main research topic in this current era of Industry 4.0, which is radically changing perspectives of industry about the technological manufacturing of components and/or materials [1, 2]. In engineering applications, one often requires the design of surface textures able to satisfy desired target responses and/or to allow rapid manufacturing and morphological changes. This can be useful, e.g., for the inline control of mechanical components [3, 4]. Moreover, the use of additive manufacturing technologies is currently promoting a detailed and unified specific description of roughness, able to take into account several lengthscales in its representation [5, 6].
Nowadays, most approaches to design roughness are based either on mimicking natural surfaces [7, 8], or on performing parametric studies via artificially generated surfaces [9, 10]. In this context, the present work proposes a method to design roughness with the aim of matching two given mechanical contact target responses as close as possible. These two targets are associated with a frictionless elastic normal contact problem. In the article, multiscale roughness is modeled by a finite series of cosinusoids which, using a biological analogy, is called the surface roughness genome. Then, the constituent waves identifying various lengthscales of roughness are properly selected and mixed based on mechanical considerations and suitable optimization tools. The resulting profile achieves a mechanical response that optimizes a suitable multiplicative tradeoff between the two different mechanical targets.
The article is structured as follows. In Section 2, we summarize the contact problem and we introduce the surface roughness genome. In Section 3, we propose an algorithm to generate prototype profiles with the aim of reproducing two given target mechanical contact responses as close as possible. In Section 4, the effectiveness of the proposed algorithm is tested on a representative numerical example. Finally, Section 5 concludes the paper and discusses possible extensions.
2. Contact Problem and Roughness Model
This study refers to the modeling and optimization of the mechanical behavior of two contacting rough surfaces. First, the contact problem and the roughness model are presented. Then, the mechanical interaction of different roughness lengthscales is used to identify two categories of roughness. Finally, an algorithm is proposed to optimize roughness, expressing the objective function of the optimization problem as a suitable tradeoff between two different target responses. To model contact, the frictionless elastic normal case is considered, taking as targets a desired stiffnessload curve, , and a desired contact areaload curve, .
2.1. Contact between Two Rough Surfaces
Following the seminal work of Johnson [11], under the assumption of linear elasticity, the nonconforming contact problem between two rough surfaces is equivalent to a contact problem between a rigid composite rough surface and a halfplane depending on effective elastic parameters; see also Barber [12] for details and proofs. Normal contact is herein controlled by imposing an approaching farfield closing displacement , whose value is computed from the tallest summit of the rough surface, whose height is denoted as . The deformation of asperities produces a normal contact traction distribution in the contact area . The average total contact pressure is computed as the sum of all the contact forces acting in this area divided by the nominal contact area. The derivative of the contact pressure with respect to the imposed displacement defines the contact stiffness .
As an example, in Figure 1, which refers to a twodimensional case, the topography of a numerically generated rough surface is presented (a), together with its deformed configuration during contact with a halfplane (b). In the figure, denotes the value of the dimensionless displacement.
(a)
(b)
In the work, roughness is described in terms of various parameters (see Section 2.2). For each choice of the set of parameters, the contact problem is solved using the Boundary Element Method (BEM); see Bemporad and Paggi [13] for a comparison of various techniques which can be used to impose the satisfaction of unilateral contact constraints.
2.2. Roughness Description over Multiple LengthScales
Referring now for simplicity to onedimensional rough profiles, their topography is herein described in terms of the real part of the Multivariate WeierstrassMandelbrot function (MWM). Ausloos and Berman [14] proposed the MWM function to describe stochastic processes with a larger number of features than in the original 2D formulation early proposed by Mandelbrot [15].
A rough profile modeled via the MWM function has then a topography depending on several parameters:
In the absence of any subscript, the base 10 is used to compute the logarithm. In (1), the parameter is a positive integer, the parameters , , and the parameter . These parameters are called genes, exploiting an analogy with biology. Similarly, the expression surface roughness genome refers to the overall ensemble of genes providing the realization of a profile over several different lengthscales.
The indices and are integers, which determine, respectively, the longest and the shortest wavelengths contributing to the definition of the rough profile . Such a rough profile is realized over an observation length , and sampled in nodes according to a resolution . Therefore, in this model, roughness is realized by varying the index between and , together with the other parameters. Taking into account the reference wavelength , the indices and are expressed as follows:In (2), the notation denotes the largest integer smaller than or equal to its argument. The number of terms in the inner summation in (1) is . The reader is referred to Majumbdar and Bhushan [16]; Wang and Komvopoulos [17]; Wu [18, 19] for further details on the spectral properties relating (2) to (1).
The particular combination of genes associated with each fixed value of the index identifies a rough profile, whose expression readsand which is named chromosome. This profile is associated with the features of roughness defined at the fixed reference lengthscale . In this way, it is possible to describe the features of a single lengthscale of roughness in terms of a chromosome and obtain the complete profile summing up chromosomes:Without loss of generality, the value is assigned to the chromosome with the longest wavelength, which refers to the coarsest realization of the profile.
2.3. Macroscale and Microscale Roughness
The nonlinear interaction of multiple lengthscales of roughness provides the collective mechanical response of a rough profile. For example, the thermal/electric contact conductance is intimately related to the stiffnessload curve , as addressed by Barber [12]. Moreover, as investigated by Paggi and Barber [20], the stiffnessload curve is mainly affected by the long wavelengths contributions. On the other hand, the contact area morphology is determined by the short wavelengths composing the power spectral density (PSD), as confirmed by many authors; see, e.g., Barber [12]; Paggi and Barber [20]; Yastrebov et al. [21].
In the following, the parameters that identify the stiffnessload curve are isolated and emphasized. Here, roughness is characterized as a result of two contributions, according to the role of each chromosome in the definition of the evolution. Considering a given observation lengthscale, the associated rough profile is split into two rough contributions and :
The profile , which represents the macroscale roughness, is obtained by summing up the chromosomes whose evolution correlates with the one of the complete profile with a correlation coefficient .
The profile , which represents the microscale roughness, is associated with the remaining set of chromosomes, having .
This distinction is shown in Figure 2. The macroscale contribution is depicted with a red line and the microscale one with a blue line. Their sum gives the complete rough profile (black line). The stiffnessload curve and the areaload curve of those profiles are shown in Figure 3.
(a)
(b)
Regarding the curves, the ones provided by macroscale roughness and by the complete profile almost overlap. The same observation can be done for the evolution for small pressures. However, the two curves diverge for N/m.
From the operative stand point, chromosomes associated with macroroughness are selected as follows. First, the evolution of a rough profile is computed by solving the frictionless normal elastic contact problem with imposed displacements in its peakvalley amplitude. To solve such a problem, the Nonnegative Least Squares method proposed in Bemporad and Paggi [13] is used. Then, the evolution is computed for each of the chromosomes in the power density spectrum of the profile, and the correlation coefficient between and the curve associated with the complete profile is computed. Chromosomes with , i.e., the ones leading to macroscale roughness, are included in a set of chromosomes . The remaining chromosomes are associated with microscale roughness.
3. Optimal Design of a Rough Profile
In the following, we consider the problem of designing a rough profile in such a way that its stiffnessload curve and contact areaload curve are close, respectively, to a given target stiffnessload curve , and a given target contact areaload curve . To do that, we propose an algorithm, called Mixed Chromosomes CrossOver (MCCO). Loosely speaking, starting from a given genome database, the algorithm mixes chromosomes belonging to genomes in the database leading to macroroughness (which make it possible to achieve the target) with chromosomes belonging to genomes in the database leading to microroughness (which make it possible to achieve also the target). Then, their parameters (genes) are optimized to refine and further improve matching with the target contact responses. In Section 3.1, we describe the MCCO algorithm. Then, in Section 3.2, we report some details on another algorithm, in this case taken from the literature, which is exploited in an inner loop of the proposed MCCO algorithm.
3.1. Proposed Algorithm: Mixed Chromosomes CrossOver (MCCO)
Before going into the details of the proposed MCCO algorithm, we introduce some notations and concepts. First, we define the similarity score which quantifies, for each index , how much the target mechanical response, , and the one associated with the th rough profile (genome) in the database are similar. For , the target and the th response are, respectively, and (where the index refers here to the th genome in the database). For , they are, respectively, and . In (6), denotes the norm, computed on the points used to discretize each mechanical response. The th mechanical response coincides with the target one when holds.
To compute consistently the similarity scores in (6), we have to assume that the two mechanical responses are evaluated on the same range of pressures (from 0 to the same maximum pressure), i.e., that . To do that, for each genome in the database, its amplitude gene is rescaled aswhere the step in (7) is possible because the gap function in the BEM formulation scales linearly with the height field, and so does the pressure; see Johnson [11]. The maximum pressure level is achieved by solving the contact problem with a farfield displacement equal to the target profile amplitude. The same remark holds for the maximum pressure level of the th genome. Moreover, all the genomes in the database refer to the same reference length .
The steps of the MCCO algorithm are reported in Algorithm 1 and are described as follows. Firstly, the similarity score in (6) is computed separately between each target ( or ) and the corresponding mechanical response ( or ) associated with each genome in the database, obtaining the vectors and with components and (Steps )().

Then, the macroscale roughness is firstly identified, finding all the genomes with a similarity score larger than . At this point, these genomes are reduced, limiting to their chromosomes leading to macroroughness. In other words, for each th selected genome, only the chromosomes that affect the target significantly (i.e., the ones for which their stiffnessload response and have correlation coefficients larger than 0.95) are kept. The set of reduced genomes obtained in this way is named (Step ).
Similarly, the genomes with a similarity score larger than are determined. Then, they are reduced, keeping only their chromosomes leading to microroughness. In other words, for each th selected genome, only the chromosomes that affect the target significantly (i.e., the ones for which their contact areaload response and have correlation coefficients larger than 0.95) are kept. The set of reduced genomes so obtained is named (Step )).
The resulting sets and contain information about different sets of chromosomes. All possible new genomes obtained by summing a reduced genome from the set and one from are stored in the set of genomes (Step ). For each such genome , the contact mechanics problem is solved via BEM, determining its and curves (Step (). Moreover, its amplitude gene is rescaled in order to satisfy the maximum pressure requirement, according to (7).
For each genome , two similarity scores are now computed, and . They refer, respectively, to the similarity of each target response ( or ) and the related mechanical response curve ( or ) associated with that genome. A mixed similarity score is also computed, by multiplying these two similarity scores, i.e., (Steps )().
Now, the three genomes in with the largest mixed similarity score are stored in the set (Step (14)). Then, for each genome in , the similarity with the target function (the mixed similarity score) is increased by optimizing the values of the genes belonging to the chromosomes related to macroroughness. As an alternative, the weighted sum of two objectives (i.e., of the squares of the two similarity scores) could be optimized, as done in [22, 23] for different problems. This approach resembles some regularization techniques applied in machine learning problems; see, e.g., Gnecco and Sanguineti [24, 25]; Gnecco et al. [26]. However, for such techniques, the optimal selection of the regularization parameter may not be straightforward. Nevertheless, the weighted product of the two objectives can be a valid alternative. See Gnecco et al. [23]; Gnecco [27] for some numerical examples.
The square of the mixed similarly score is maximized by applying the Globally Convergent Method of Moving Asymptotes (GCMMA) algorithm [28, 29] (Step (16)); see Section 3.2 for further details. Then, starting from the optimized genome, the GCMMA is applied to its remaining chromosomes, which are associated with microroughness (Step (17)), considering the same objective function of Step (11).
Finally, the new genome with the maximum obtained value of the similarity score is identified (Step (19)).
3.2. GCMMA and Its Use inside the MCCO Algorithm
The Globally Convergent Method of Moving Asymptotes (GCMMA) [28, 29] is an algorithm to solve nonlinearly constrained optimization problems. It replaces each such problem with a sequence of approximating optimization subproblems. These are easier to solve, and are still nonlinearly constrained.
The term moving asymptotes refers to the functions which are used in the formulations of the approximating subproblems. Such asymptotes usually change from one iteration of the GCMMA to the successive one. The GCMMA is globally convergent in the sense that, for each possible initialization of the vector of parameters to optimize, the algorithm is proved to be convergent to a stationary point of the initial problem.
In Steps (16) and (17) of Algorithm 1, the GCMMA is applied by taking as objective function not the mixed similarity score , but its square. This is done in order to get a smoother objective function. The GCMMA iterative solution is obtained after a number of steps, starting from an initial choice for the vector of parameters to optimize. For the simulations, has been chosen, in order to keep the computational time per iteration small.
To apply the GCMMA algorithm inside Algorithm 1, a MATLAB implementation has been made, based on the MATLAB functions mmasub.m and subsolv.m, developed by Svanberg himself. Technical details on the implementation of these two MATLAB functions are provided in Svanberg [30]. At each step of the GCMMA, the partial derivatives of the objective function with respect to all the optimization variables are computed numerically, via a small perturbation for each such variable. Moreover, the amplitude parameter is rescaled according to Eq. (7).
The GCMMA algorithm is applied with several different initializations, in order to increase the probability of obtaining a good constrained local maximizer of the original objective function. Finally, the best solution, i.e., the one associated with the greatest (square of the) mixed similarity score obtained during the various optimizations, is produced as output, for both Steps (16) and (17) of Algorithm 1.
To save computational time, the number of variables in each optimization problem is reduced as described in the following. The parameters (genes) , and determine the frequency spectrum, i.e., the interaction among several lengthscales. Hence, they are not modified during the optimization, but they are fixed to their original values. The parameter is rescaled at each optimization step in order to satisfy the pressure requirement, i.e., applying Eq. (7). So, only the phases are considered as parameters to optimize. In the optimization, their values are constrained in the range between of their initial values, in order to preserve the main features of the original chromosomes. Moreover, for both Steps (16) and (17) of Algorithm 1, only the phase genes associated with the subset of chromosomes related, respectively, with macroroughness and microroughness, are optimized.
4. Numerical Experiments
The MCCO algorithm is now tested, designing roughness in a realization length µm. First, an artificial genome database is generated to emulate a priori known database of real rough profiles (see Section 4.1). The two targets and are therefore selected not to match any of the mechanical responses associated with one of the rough profiles in the database. The two selected targets are shown as red lines in Figures 4(a) and 4(b), respectively. Numerical results about the application of the MCCO algorithm to this specific case are reported in Section 4.2.
(a)
(b)
4.1. Genome Database
For all the profiles in the database, the amplitude gene is fixed to . Moreover, the main wavelength is set equal to µm. Twenty different pairs of values for and have been generated according to a Sobol’ sequence [31]. To generate the phase matrix , has been imposed for all the genomes. Three phase vectors of thirty elements each have been generated, with values from 0 to for their components. Such values have been also determined using a Sobol’ sequence. Combining the pairs with the phase vectors, a database made of genomes has been obtained.
Each profile in the database has been discretized using nodes. To formulate and solve each frictionless normal contact problem, the elastic modulus has been set equal to MPa. Then, the mechanical response has been computed in equally displaced farfield displacements starting from the tallest summit of the rough profile, proceeding likewise in Johnson [11]; Bemporad and Paggi [13].
4.2. Numerical Results
The three best genomes obtained from the MCCO algorithm are summarized in Table 1. The table reports, for each such th genome, its mixed similarity score , and the similarity scores and of its mechanical responses with respect to the two targets. Also, it reports the genomes used in the initialization of the GCMMA algorithm in Steps (16) and (17) of Algorithm 1. Such genomes are used to represent, respectively, macro and microroughness. The first index refers to the pair (see Figure 5(a)), whereas the second one to the vector of phases (see Figure 5(b)).
(a)
(b)
The best approximation of the target contact responses is provided by the genome , whose mechanical responses and are shown by blue curves with round markers in Figures 6(a) and 6(b), respectively. Looking at the contact areaload curves reported in Figure 6(b), one can notice that the mechanical responses of all these three solutions diverge consistently from the target curve in the range of contact pressures between N/m and N/m.
(a)
(b)
(a) Macroroughness
(b) Microroughness
(c) Roughness
As far as the stiffnessload curve is concerned (see Figure 6(a)), the and curves show a very good overlap with the target . The second solution is the only one associated with a different genome representing macroroughness, and its trend is not close to in the range of contact pressures between N/m and N/m.
The curve almost overlaps for N/m. Indeed, the genomes and are composed by the same starting sets of chromosomes representing macroroughness. Also the curves and are very similar (see Figure 6(b)), and quite close to the target , apart from an intermediate range of values for the contact pressure.
For both genomes and , the same pair is used to represent microroughness, but a different vector of phases is used in the GCMMA initialization in Step (17) of Algorithm 1. This highlights the dominant role of macroroughness in the frictionless elastic normal contact response.
The three rough profiles obtained from the MCCO algorithm are visualized in Figure 7, which emphasizes the roles of macroroughness (see Figure 7(a)), and microroughness (see Figure 7(b)). Each rough profile in Figure 7(c) is obtained by summing these component profiles (represented by the same colors and line types).
In Figure 7(a), the macroroughness of each solution is depicted. The dashed blue and the dashdotted red profiles, which are associated, respectively, with the genomes and , have a similar macroroughness topography, with small differences, e.g., in the range between N/m and N/m.
As observed before, these solutions are originated from the same set of chromosomes as regards the macroroughness level. Moreover, they have also quite close values of similarity scores (see Table 1).
The microroughness of genomes and is shown in Figure 7(b). The figure shows that they have similar microroughness, as they are associated with the same pair (). The resulting topographies of the genomes and are shown in Figure 7(c) and are also quite similar. Moreover, these two genomes have quite close mechanical responses (see Figure 6).
The topography of genome is quite different from the ones of and , since its macroroughness component differs significantly.
Finally, the spectral decomposition of the three genomes is shown in Figure 8 (in the figure, has the same dimension as ). The difference in amplitude is due to different values for the amplitude gene , which are obtained by applying Eq. (7) with respect to the threshold pressure.
(a)
(b)
(c)
All the results of Fast Fourier Transform (FFT) filtering show a power spectral density of the complete profile containing several peaks. For all the three genomes, the first part of the PSD does not present any evident peak. For the genome , the FFT filtering peaks are close to the discrete spectrum of the genome (see the squares in the central part of the spectrum in Figure 8(b)). Moreover, the squares in Figure 8(b), which represent the microroughness, are six, whereas the peaks observed from the FFT filtering are five.
Looking at the solutions and , which show similar mechanical responses and topographies, one can notice that they present also similar PSDs in the lowfrequency part of the spectrum. However, their highfrequency behaviors are different, because they derive from different genomes.
5. Conclusion
In this paper, the Mixed Chromosomes CrossOver (MCCO) algorithm has been proposed to optimize a roughness topography to match two different target mechanical contact responses, namely the normal contact stiffness and the real contact area ; both are expressed as functions of the average contact pressure . The algorithm has been inspired by an analogy with genetics, according to which the MCCO mixes the macroroughness and the microroughness of two different “genomes”. Then, an application of the algorithm to a representative numerical example has been described, showing very promising results. One can conclude that two particular classes of genomes can be identified and combined to match both the and targets. Particularly, chromosomes associated with macroroughness are fundamental for an optimal design.
Still, the results show that there is room for possible improvement, which could be obtained, e.g., by future variations of the MCCO algorithm. Particularly, the algorithm could be improved in the selection of the chromosomes to combine. Future work will also concern the introduction of additional interactions between the rough surfaces in contact, such as friction, adhesion, wear, and so on. Such mechanical interactions are difficult to predict using BEM, and it is probably necessary to move to the Finite Element Method (FEM) to tackle nonlinear multifield coupled contact and fracture problems with complex and realistic geometries. For example, the case of elastoplastic contact is certainly a fundamental aspect to be investigated, and plastic deformation can be taken into account by using ad hoc computational models.
This work is a first step in the generation of surface topographies achieving desired target mechanical responses. The problem investigated in the article and the MCCO algorithm are potentially useful for the production of surface pressure sensors in sealing applications. In this case, the contact has to be assured for a contact pressure larger than some threshold. At the same time, the electrical conductivity is associated with the curve. Then, the curve might change its trend suddenly above the contact pressure threshold. Finally, the curve may be chosen in such a way to improve sealing.
Another important potential application refers to the problem of surface morphing. In this case, roughness may be modified in time based on external stimuli. Morphing is a quite novel concept, which has been recently applied to smart structures and devices. However, its application to surfaces is a very new research area.
The problem considered in the work is important also in view of the fact that, nowadays, additive manufacturing technologies simplify the realization of topographies, enabling the production of surface prototypes. 3D printing is one of the most fascinating and robust technology, due to its very fast progress and impact on several industrial sectors, which promotes technology transfer and patent applications. Applications of the proposed methodology to surface printing could open new research areas.
Data Availability
The MATLAB data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare that there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
M. Paggi acknowledges the support of MIURDAAD to the Joint Mobility Program 2017 “MultiScale Modeling of Friction for Large Scale Engineering Problems”. G. Gnecco acknowledges the support of FFABR (Fondo per il Finanziamento delle Attività Base di Ricerca) from the Italian Ministry of Education, University and Research (MIUR).