## Abstract

In thin-gauge cold rolling of metal sheet, the surface roughness of work rolls (WRs) is known to affect the rolled sheet surface morphology, the required rolling load, and the roll wear. While modeling of rough surfaces using statistical asperity theory has been widely applied to problems involving semi-infinite solids, the application of asperity distributions and their elastic-plastic behavior has not been considered in roll-stack models for cold sheet rolling. In this work, a simplified-mixed finite element method (SM-FEM) is combined with statistical elastic-plastic asperity theory to study contact interference and coupling effects between a rough work roll (WR) surface and the roll-stack mechanics in cold sheet rolling. By mixing equivalent rough surface contact foundations, Hertz foundations, and Timoshenko beam stiffness, an approach is created to efficiently model interactions between the micro-scale asperities and the macro-scale roll-stack deformation. Nonlinearities from elastic-plastic material behavior of the asperities and the sheet, as well as changing contact conditions along the roll length, are also accommodated. Performance of the multi-scale SM-FEM approach is made by comparison with a continuum finite element virtual material model. 3D studies for a 4-high mill reveal new multi-scale coupling behaviors, including nonuniform roughness transfer, and perturbations to the sheet thickness “crown” and contact force profiles. The described multi-scale SM-FEM approach is general and applies to rough surface contact problems involving plates and shear-deformable beams having multiple contact interfaces and arbitrary surface profiles.

## 1 Introduction

Because work roll surface topography plays a major role in both roll-bite tribology and product surface quality in the cold rolling of metal sheet or strip (depicted in Fig. 1), researchers have attempted to build models to understand the micro-scale contact interactions between the two relevant bodies of interest, namely, the work roll and the sheet (strip), wherein at least one body, usually the work roll due to its greater hardness, is assumed to have rough surface topography.

In 2014, Wu et al. [1] summarized two feasible methods to incorporate the contact mechanics of rough surfaces into rolling models—the finite element method/analysis (FEM/FEA), and the method based on inhomogeneous deformed plane-strain theory. Within the rough surface roll-bite region, both of these methods can predict normal and shear stresses, lubrication behavior and friction coefficient, as well as transfer characteristics of the roll roughness to the sheet. Hence, both methods provide some understanding of the resulting sheet surface quality, as reported in Refs. [2–6]. The main challenge in using continuum FEM, however, stems from the extremely high computational cost when attempting to simulate a true micro-scale surface topography. On the other hand, the plane-strain theory method is limited to 2D representations of the rolling process, and it thus neglects any possible out-of-plane coupling effects arising from 3D bulk deformations of the rolls such as bending, shearing, and variations in Hertzian flattening along the roll axis direction. Even for smooth contact conditions, the bulk-body deformations are challenging to predict due to the complexity of both roll-stack contact configurations, e.g., 20-high cluster mills [7], and the complexity of sheet thickness profile control mechanisms, e.g., continuously variable crown (CVC) roll shifting mechanisms [8]. In addition to these respective shortcomings of the continuum FEM and inhomogeneous deformed plane-strain theory methods, no published research dealing with surface roughness in the cold rolling of metals has yet considered the elastic-plastic constitutive behavior (which comprise “stiffness” characteristics of the asperities) as being distinct from the bulk elastic stiffness of the work rolls. Nonetheless, when rolling thin-gauge sheet below ∼0.3 mm, researchers have reported that mill operators can “sense” roughness effects on the sheet, possibly because the relative roll surface deviations may be large enough such that the asperities absorb sufficient contact energy to noticeably change the relationship between the contact force [9] and the sheet thickness reduction or elongation [4], thereby possibly even affecting the sheet flatness. Asperity stiffness therefore potentially represents an important factor to better understand the overall 3D roll-strip contact mechanics relationship, which, would thus require consideration of the bulk elastic behavior of the rolls, the elastic-plastic deformation of the sheet, and the elastic-to-highly plastic behavior of the roll asperities, in accordance with magnitude of the applied rolling load.

### 1.1 Statistical Asperity Contact Modeling.

Statistical asperity contact modeling, which is discussed below and incorporated into the simplified-mixed finite element method (SM-FEM) in this study, provides for the most popular technique to predict asperity interfacial normal and shear stiffness during rough surface contact. The theory allows for statistical inference of contact behaviors using probability density functions to characterize the geometric nature of asperities (i.e., rough surfaces). As described in Sec. 2, Statistical asperity theory is applied in this work to create an equivalent rough surface foundation stiffness contribution based on mathematical integration of the area-normalized asperity stiffness over the “roll-bite” contact region between the upper work roll and the sheet in a 4-high cold rolling mill. The resulting equivalent, rough surface foundation is adapted into the SM-FEM approach using an efficient series combination that accounts for both elastic and plastic asperity stiffness behaviors. As detailed later, the stiffness for nonlinear elastic and plastic behaviors is obtained from derivatives of the corresponding load versus displacement conditions.

The first statistical asperity model was introduced by Greenwood and Williamson [10]. Referring to Fig. 2(a), their approach applied Hertz elastic contact theory [11] to predict the contact behavior of a single asperity while simultaneously adopting a statistical modeling technique [12] to accommodate distributions of asperities. The asperities were assumed to be hemispherical with heights following a Gaussian distribution. Because a fully elastic approach was taken, however, these early models only provided adequate prediction under small relative displacement (or small interference). Material constitutive models for elastic-plastic asperity deformations subsequently became an important field of study and continue to be researched today. Kogut and Etsion [13,14] developed empirical constitutive relations based on observations from continuum finite element analysis. Proposed in their method were four deformation regions, including fully elastic, contact subsurface elastic-plastic, contact surface elastic-plastic, and fully plastic. Available experimental data illustrated that Kogut and Etsion’s model exhibited strong agreement with the contact stress and true contact area of the rough surface [15]. Xiao and Sun [16] recently introduced a “virtual material” approach that involved including a thin, relatively softer, nonphysical material layer to represent the rough surface contact interface (see Fig. 2(b)). A virtual elastic modulus, shear modulus, and Poisson ratio were derived using Hertz’s analytical normal contact relation [17] in combination with classic statistical asperity contact theory. Xiao and Sun’s results showed good agreement with Kogut and Etsion’s statistical asperity theory method in which the material constitutive behavior regions were classified empirically with the aid of finite element analysis [14].

Aside from normal contact flattening, asperities develop shear stresses during sliding between the adjacent bodies, with or without lubrication [18]. In cold rolling, the asperity/lubrication mixed regime governs friction [19], and both “dry” and lubricant-based frictional conditions show strong relationships to roughness, rolling speed, reduction, and forward slip [20]. The friction coefficient is also a function of contact force [21]. Note that because friction effects and asperity shear interactions are indirectly represented in the rolling force versus thickness reduction (from industry data in this work), the model predictions avoid the requirement to separately address lubrication and shear effects.

### 1.2 Existing Surface Roughness Representations in Rolling Models.

In many practical rough surface contact problems, the macro-scale or bulk 3D deformations have coupling influences on the micro-scale asperity deformations, thus affecting the true contact area, maximum contact force and surface stress for prescribed relative displacements or for free vibration [22]. Such coupling effects, which semi-infinite solid or “slab” analyses *cannot* reveal, may also occur in the cold rolling process wherein the macro-scale, bulk-body deformations in the roll-stack are characterized by bending, shearing, and variations in Hertzian flattening along the roll lengths. Still, rolling mill studies to date have not adequately examined such possible coupling effects. For instance, the experimental investigation of a two-high “skin pass” mill parametrically summarized only the in-plane coupling between work roll surface roughness and macro-scale work roll flattening (at a section of the roll and sheet) based on rolling force, work roll diameter, sheet reduction, and roll velocity [23,24]. The results thus revealed 2D relationships between rolling parameters and roughness transfer, but neglected consideration of any 3D coupling effects, such as bending of the rolls, which might lead to variation in roughness transfer along the roll axis direction. Indeed, even though it is well-known in the cold rolling process that the bulk-body bending, shearing, and nonuniform flattening of rolls significantly affects contact load distributions between the rolls and sheet along the roll axes, the computational expense of 3D continuum finite element analysis has led researchers to compromise with more feasible 2D (plane-strain) models of the asperity contact mechanics within the roll-bite and of the subsequent roll roughness transfers to the sheet. For example, Kijima [3] created a 2D finite element model that included an elastic-plastic strip (sheet) and an elastic outer layer surrounding a rigid work roll. The study noted that for the same unit rolling force and roughness topography a relatively larger roll transfers more roughness to the sheet than does a smaller roll. Since the roll roughness was a modeled as fully elastic in Kijima’s simulations, however, the results showed greater rolling force prediction compared with the experiments (due to stiffer asperities), and the difference was accentuated when the roll surface was rougher. This overprediction in force for a prescribed relative roll displacement is similar to that exhibited by Greenwood and Williamson’s fully elastic statistical asperity model [10]. Similarly, Wu et al. [6] used a plane-strain continuum finite element model to study the non-steady-state texture transfer during and shortly after initial “threading” of the sheet into the roll-bite, as well as texture transfer from the roll to the strip during steady-state rolling. The roll was modeled as rigid and included a random surface profile deviation along the contact arc with the sheet, for which plastic deformation of the latter was governed by the Tresca yielding criterion. As expected, and similar to Kijima’s results, Wu’s rigid roll assumption overpredicted contact stress, force, and roughness transfer for a prescribed relative displacement. Since elastic yielding behavior of asperities differs from the roll body, it requires separate consideration in cold rolling because the asperities deform significantly more than the elastic roll body.

### 1.3 Surface Roughness Roll-Stack Model Need in the Cold Rolling Process.

Considering related studies to date, an efficient approach does not currently exist to combine micro-scale roll roughness with 3D macro-scale roll displacements to achieve a coupled, multi-scale roll-stack deformation model. Recently, however, the highly efficient SM-FEM has been applied with industrial measurements of work roll diameter profiles to help understand the effects of micro-scale roll grinding errors on the sheet thickness profile and contact force distributions in a 4-high cold rolling mill. The results showed strong agreement to those from large-scale continuum FEM, and were obtained at significantly lower computational cost (a few seconds versus several hours with identical computing hardware) [25]. Although micro-scale roll diameter profile deviations were considered, the model was unable to elucidate any coupling effects between roll surface roughness and the bulk roll-stack deformations for the following three reasons: (1) spatial resolution of roll diameter deviations along the roll axis was about 10 mm (per industry data), which represents too large a scale than required to model surface roughness; (2) the SM-FEM formulation only accounted for geometric changes in the work roll diameter, and thus it did not consider the stiffness of the roll surface profile deviations as being distinct from the bulk stiffness of the rolls; (3) constitutive behavior of roll surface was assumed to be purely elastic. Limitations 2 and 3 represent similar simplifications to those in both Kijima’s [3] and Wu et al.’s [6] models, and if used in this study would thus lead to overprediction of the rough-interface contact force for prescribed displacement boundary conditions, since plastic asperity deformation would not be considered. In this study, therefore, an improved SM-FEM model is formulated and proposed that incorporates statistical asperity theory as well as *distinct* elastic-plastic constitutive behavior for the asperities. As is shown, the improved method provides the capability to efficiently simulate the coupled roll-stack deformations and the surface roughness contact mechanics for general cold rolling mill configurations, as well as any rough surface contact problem involving plates and/or shear-deformable beams having multiple contact interfaces and arbitrary surface profile geometries.

The remainder of this paper is organized as follows: Sec. 2 details the mathematical implementation of the asperity stiffness in the SM-FEM formulation. Section 3 gives a comparison of results for the new roughness-capable SM-FEM, the original SM-FEM, and large-scale continuum finite element simulation for a 2D case (since modeling of a 3D multi-scale surface roughness/roll-stack mechanics problem is still computationally infeasible using continuum FEM and prohibits comparison). In Sec. 4, 3D case studies for a 4-high cold rolling mill are demonstrated using the new SM-FEM, wherein the work roll includes different surface roughness magnitudes; the simulated contact force distributions and sheet thickness profiles reveal coupling effects between the micro-scale asperity deformations and the macro-scale roll-stack deformations. Remarks on the improved, roughness-capable SM-FEM and its implications for practical applications are in Sec. 5.

## 2 Multi-Scale Model Formulation

Described next is the procedure for assembling and solving the nonlinear global system of equations using the SM-FEM for static, elastic-plastic asperity problems involving mutually contacting cylinders and plates. In contrast to this proposed improved SM-FEM, with the original SM-FEM approach, slab methods, and continuum FEM, surface profile deviations or roughness were modeled in terms of geometric surface deviations during contact interference calculations. With the SM-FEM formulation proposed here, the surface roughness layer is modeled analogously to (but different from) a virtual material. Statistical asperity stiffness is used here to represent equivalent elastic-plastic constitutive behavior, as depicted in Fig. 3. The actual rough surface geometric deviations are therefore not explicitly required with the new formulation, which enables efficient solutions that are also capable of modeling 3D rough surface cold rolling problems.

Mathematically, the equivalent elastic-plastic asperity stiffness is obtained by integrating individual quasi-parallel asperity stiffnesses along the roll-bite contact arc region and then creating an equivalent work roll/sheet foundation stiffness through series combination of the asperity stiffness integral with the bulk Hertzian work roll elastic stiffness and the equivalent elastic-plastic stiffness of the sheet. Figure 3 depicts components of the series stiffness combination. Integrating the asperity stiffness along the small contact line (*z*-direction) effectively aggregates the influence of roll-bite asperities on the vertical deformation (and sheet thickness) at each transverse location *x* along the sheet width; this integration is analogous to the forward slip^{3} that occurs between the work roll and the sheet during rolling. In fact, the contact arc-based integration of the roughness is well suited to the aim of this work, which is to better understand the impact that roughness has on the roll-stack behavior and subsequent exit thickness profile of the sheet (which are a function of all asperities in the roll bite), rather than to understand the details of roughness effects *within* the roll bite itself. In addition, note that series combinations in both discrete spring/beam models and continuous foundation/beam models have been widely used in both static roll-stack modeling [7,26,27] as well as dynamic modeling [28]. Because of nonlinear, Hertz-type bulk elastic deformations of the rolls, elastic-plastic deformation of the sheet, and highly plastic deformations of the asperities under typical operational loads, the proposed SM-FEM equations are solved using a modified Newton–Raphson iterative procedure that updates the discretized equivalent series foundation stiffness according to the modified constitutive behaviors of the asperities, sheet, and rolls during progressive computations of the unit contact force. The resulting nonlinear, multi-scale contact model of roll surface roughness and roll-stack deformation incorporates coupled Winkler foundations, 3D Timoshenko beam elements, and statistical asperity distributions into a global stiffness system requiring iterative solution. Details on the new formulation appear below, where the original SM-FEM equations [7,29,30] are first briefly introduced.

*K*_{G}, represents direct superposition of all continuous Winkler foundation stiffness contributions,

*K*_{F}, and all bending and shearing stiffnesses,

*K*_{T}, as indicated in Eq. (1). Detailed definitions of both

*K*_{F}and

*K*_{T}are given in the Appendix

*K*_{G}, and the global load vector,

**, are explicitly stated as functions of the global nodal displacement vector,**

*f***, which contains all six translational and rotational degrees of freedom at each node**

*u*### 2.1 Equivalent Foundation Stiffness.

*K*_{F}, stems from strain energy, $U\u03f5$, for two adjacently contacting beams/plates based on a contact interference,

*δ*, between points on their respective axes at axial contact coordinate,

*x*:

*k*

_{feq}, defined by Eq. (4), is a series combination of the bulk contact stiffness of body 1, namely, $kf1$, the bulk contact stiffness of body 2, namely $kf2$, and the asperity stiffness, $kfa$, integrated along the arc of contact between bodies 1 and 2 (see Fig. 3, where body 1 is the sheet and body 2 is the work roll):

*δ*

_{2}, is defined by Eq. (6)

*b*is the contact width that varies with unit contact force,

*f*, according to Eq. (7), and where

*E*is the elastic modulus,

*v*

_{p}is the Poisson ratio, and

*D*is the diameter (assumed infinity for the sheet)

Note that in this paper the contact width is taken as the roll-bite contact arc due to the relatively small thickness reduction; for larger reductions, a more appropriate contact width dimension should be used, the simplest (which assumes a arc circular) being the classic Hitchcock deformed roll diameter relation [33,34].

The equivalent foundation modulus for the rolled sheet, $kf1$, can be identified from the secant or tangent in the functional relationship between unit rolling force and sheet reduction at the anticipated operating point [35]. Use of the secant instead of the tangent may reduce accuracy but can promote constancy in $kf1$, and thus aid computational stability in the presence of significant hardening behavior in the rolling force versus reduction relation during iterative solution [36]. Here, $kf1$ is computed from the secant relationship between unit contact force and thickness reduction (engineering strain) from industrial measurement (e.g., Fig. 4), although this relation may be obtained from any suitable roll-bite force model.

*k*

_{fa}, is based on integration of the asperity stiffness per unit contact area,

*κ*

_{a}, along the contact arc (see Fig. 3)

*κ*

_{a}is found by differentiating the contact pressure,

*p*, with respect to interference,

*δ*:

The contact pressure *p* is obtained using either the Greenwood and Williamson (GW) [10] or the Kogut and Etsion (KE) [14] statistical asperity method formulation (see Eqs. (10) and (15)). Although, as mentioned earlier, many papers have shown that the GW method only applies for small contact interferences because it assumes fully elastic asperity behavior, surface hardening (localized strain-strengthening) can cause contact interfaces to behave akin to larger elastic stiffness during the load cycle increment [37]. Hence, since combined wear and plastic deformation of asperities occur in cold rolling, where the latter may manifest to some degree in micro-surface hardening, both the GW and KE models are incorporated into and evaluated with the improved SM-FEM.

*p*

_{GW}, and asperity stiffness per unit area, $\kappa aGW$, are expressed by Eqs. (10) and (11):

In Eqs. (10) and (11), $d\xaf$ is the asperity displacement normalized by the standard deviation of asperity surface height, *σ*_{s}. Term $\varphi (h\xafa)$ is the Gaussian (normal) probability distribution function^{4} for asperity summit height, $h\xafa$, which is normalized by *σ*_{s} to give domain $h\xafa\u2208$ [−3, 3] for ∼99.7% of asperities. Term $h\xafa\u22123+d\xaf$ is equal to the normalized contact interference, $\delta \xaf$, since $h\xafa\u22123$ represents the initial contact position of an asperity, assuming the Gaussian distribution. Terms *r* and *η* refer to asperity radius and density, respectively, where all asperities are assumed to have identical radius. Despite its simplicity, for contact problems emulating semi-infinite solids, the GW model exhibits strong agreement with experimental results for contact force and contact area, given prescribed displacements, prior to reaching half the yielding point [38]. Beyond half the yield point, however, agreement diverges quickly.

*p*

_{KE}, and asperity stiffness per unit area, $\kappa aKE$, are expressed as in Eqs. (15) and (16):

Note that to apply the KE method, the SM-FEM iteratively provides the unit contact force, *f*, and contact pressure, *p*. Due to the absence of an explicit analytical expression between the asperity stiffness and contact pressure, both the GW and KE solutions involve a simple Golden Section search [40] to determine $d\xaf$ and *κ*_{a}.

It is important to note that although both the GW and KE expressions assume asperity distributions to be Gaussian (normal), they can be adapted with any distribution to obtain the unit asperity foundation stiffness via Eqs. (8) and (9) as long as the probability distribution function for asperity height, $\varphi (h\xafa)$, is given.

### 2.2 Contact Interface Determination.

*y*

_{12}, between the axes of bodies 1 and 2 is equal to the sum of the relative displacement,

*v*

_{1}−

*v*

_{2}, and the corresponding relative initial coordinate position,

*y*

_{1c}−

*y*

_{2c}:

*f*(

*x*), the contact interferences of the work roll, sheet, and work roll asperities each have reciprocal ratios to their respective stiffnesses. Hence, the interferences,

*δ*

_{1}(

*x*),

*δ*

_{2}(

*x*), and,

*δ*

_{a}(

*x*), are readily obtained from the unit contact force,

*f*(

*x*), during iterative solution:

In contrast to explicitly including surface profile deviations in the interference calculation, as in Ref. [25], in the above formulation the roll surface profile is *geometrically uniform*, but includes micro-scale wavelength variations in the foundation stiffness along the roll axis direction to account for the roll roughness via Eq. (19).

### 2.3 Solution Method.

*m*identifies the load step and Δ

**is the load increment:**

*f*

*f*_{1}is assumed in the initial load step (i.e.,

*m*= 0) when computing estimated foundation contact stiffness terms in $K~F0$, which are added to the constant Timoshenko beam stiffness,

*K*_{T}, to obtain the initial estimated global stiffness matrix, $K~G0$. The next step’s estimated displacement increment, $\Delta u~$, is determined using the current stiffness and load increment:

**in Eq. (20) is the difference between the exact and estimated load, i.e.**

*e*## 3 Plane-Strain Contact Stiffness and Interference Results, and Comparison With Continuum Finite Element Analysis

As mentioned, prior to presenting 3D multi-scale contact mechanics results for a 4-high mill (in Sec. 4), the results detailing 2D plane-strain stiffness and interference behaviors with the multi-scale SM-FEM formulation are examined for a 4-high mill cross section to facilitate understanding of the relations between asperity deformations and bulk deformations. The contact stiffness and interference results are discussed in Secs. 3.1 and 3.2, respectively, considering two different surface roughness conditions on the upper work roll. Comparisons of the 2D model results to a solution obtained using a continuum finite element model with virtual material are provided in Sec. 3.3. Geometry for the 2D work roll and sheet contact problem is as seen earlier in Fig. 3 (Section *A-A* view). The 301 stainless steel sheet (body 1) has 0.233 mm thickness, and the work roll (body 2) diameter is 50.8 mm. Material for the work roll is carbon steel with elastic modulus 2.07 × 10^{11} Pa, Poisson ratio 0.29, and Brinell hardness 1.96 × 10^{9} Pa. For comparative purposes, Table 1 lists two different sets of surfaces topography (roughness) parameters that are applied to the work roll. Case 2 in Table 1 represents a rougher surface condition than Case 1. Thus, per unit contact area, Case 2 implies fewer asperities but increased asperity height. Under the same contact force, the asperities of Case 2 can be expected to undergo greater contact interference than those of Case 1. The sheet is assumed to be perfectly smooth.

Property | Case 1 | Case 2 |
---|---|---|

Surface height standard deviation, σ_{s} (μm) | 1.58 | 8.75 |

Asperity height standard deviation, σ_{a} (μm) | 1.15 | 8.61 |

Asperity radius, r (μm) | 6.34 | 3.27 |

Asperity density, η (μm^{−2}) | 7.32 | 2.51 |

Property | Case 1 | Case 2 |
---|---|---|

Surface height standard deviation, σ_{s} (μm) | 1.58 | 8.75 |

Asperity height standard deviation, σ_{a} (μm) | 1.15 | 8.61 |

Asperity radius, r (μm) | 6.34 | 3.27 |

Asperity density, η (μm^{−2}) | 7.32 | 2.51 |

### 3.1 Contact Stiffness Results for Section of 4-High Mill With Rough Upper Work Roll.

Figure 5 shows the resulting relationship between unit contact stiffness *k*_{f} (foundation stiffness) and unit contact force for the sheet (body 1), the work roll “bulk” as separate from the asperities (body 2), and the integrated work roll asperity stiffness along the roll-bite contact arc. Note that $kf1$ for the sheet is derived from either the elastic modulus (pre-yielding) or from the secant of the empirical unit rolling force versus plastic reduction relationship in Fig. 4 (post-yielding), depending on the unit force magnitude. As seen in Fig. 5(a), once the unit contact force reaches about 2.15 kN/mm, the sheet deformation transitions from elastic to plastic, and therefore the pre-yielding constant (elastic) $kf1$ is significantly decreased to a post-yielding secant value from the unit force versus reduction. This discontinuity in $kf1$ not only alters the overall equivalent foundation stiffness between the roll and sheet axes, it also shifts the greatest fractional contribution to the total contact interference displacement from the work roll and its asperities to the rolled sheet, as is shown later in Fig. 7. Referring again to Fig. 5, as anticipated, the fully elastic GW model predicts greater $kfa$ than the elastic-plastic KE model for both work roll roughness cases in Table 1. The enlarged view in Fig. 5(b) clearly indicates that, for both roughness cases, the asperity stiffness, $kfa$, is less than the bulk work roll stiffness, $kf2$, when the unit contact force is very small. Consequently, the total contact interference between the roll and sheet axes can be attributed more to the work roll’s asperities than to the body of the work roll in this small unit contact force region. This circumstance is also quantitatively shown in Fig. 7, which indicates the ratio of interference, *δ*, to total displacement *d* (between the roll and sheet axes) for the sheet, the work roll, and the asperities. Figure 5 also shows that, for both the GW and KE asperity theories, the Case 2 surface roughness leads to softer asperity stiffness, $kfa$, than for roughness Case 1. This is because Case 2 has a relatively rougher surface than Case 1.

Figure 6 shows a comparison of the equivalent unit foundation contact stiffness, *k*_{feq}, as a function of unit contact force with and without considering the integrated asperity stiffness. Hence, Fig. 6 compares the “original” perfectly smooth SM-FEM formulation result for *k*_{feq} with the new multi-scale, roughness-capable SM-FEM formulation for *k*_{feq}. With the new SM-FEM formulation, both GW and KE theory results are indicated in Figs. 6(a) and 6(b), which illustrate, respectively, roughness Cases 1 and 2. Since the surface roughness is included in the overall equivalent foundation stiffness using a series combination of the integrated parallel asperity stiffnesses between the work roll and the sheet (see again Fig. 3), the roughness effectively reduces the overall equivalent foundation in the contact interface. Note that $kfeqGW$ is greater than $kfeqKE$ due to the respective purely elastic and elastic-plastic material constitutive models. Therefore, $kfeqGW$ is closer to $kfeqoriginal$ and, furthermore, $kfeqGW$ approaches $kfeqoriginal$ faster with increasing contact force than does $kfeqKE$ because of the series stiffness combination governed by Eq. (4). In addition, with the smaller roughness of Case 1 shown in Fig. 6(a), both $kfeqGW$ and $kfeqKE$ approach $kfeqoriginal$ faster as the contact force increases than they do for the rougher Case 2 condition in Fig. 6(b). This is because, as seen in Fig. 5, for both GW and KE, the $kfa$ exceed the bulk roll stiffness, $kf2$, more quickly for Case 1 than for Case 2. Thus, the overall equivalent foundation stiffness, *k*_{feq}, is progressively influenced by the corresponding most rapidly increasing *k*_{f} contribution in Eq. (4). As long as the unit contact force exceeds that corresponding to the sheet’s elastic limit (yield point), the equivalent foundation stiffness rapidly converges, with or without the asperities. Again, the elastic to plastic discontinuity in the rolled sheet instantly reduces *k*_{feq} and thereafter changes the most significant relative contribution to the overall contact interference from the work roll to the sheet, as discussed next and seen in Fig. 7.

### 3.2 Contact Interference Results for Section of 4-High Mill With Rough Upper Work Roll.

Figure 7 shows both the total displacement between the sheet and work roll axes (right axes on each figure) as well as the ratio of individual contact interferences to the total displacement (left axes on each figure), both as a function of the unit contact force. Note that the total contact displacement, *d*, increases almost linearly before the sheet yielding point. At the yielding point, *d* has a dramatic and discontinuous increase, followed by further, quasi-linear, post-yielding increase at a rate larger than before yielding. The ratios of contact interferences to the total displacement, *δ*/*d*, decompose the significance of each component (asperities, sheet, work roll bulk) to the total contact interference. In each plot of Fig. 7, as the unit contact force increases from zero the asperity interference contribution, *δ*_{a}/*d*, drops rapidly because $kfa$ increases rapidly, especially for the Case 1 roughness with the GW model (Fig. 7(a)). In comparison to Fig. 7(d), which represents the Case 2 roughness with the KE model (and the lowest stiffness in Fig. 5), the decrease in the asperity interference contribution, *δ*_{a}/*d*, is less significant than for the other results. To understand the contribution of the relative bulk work roll interference, *δ*_{2}/*d*, note that it depends on the contact arc length given in Eq. (7); as the unit contact force increases, *δ*_{2}/*d* initially grows very rapidly, but then the rate of increase slows. This is because in the very low unit force region the contact force initially acts to deform the asperities, but soon afterward the force begins to deform the work roll body. Finally, as mentioned earlier, once the sheet yields, the largest relative contribution to the contact interference transitions to the sheet, and as a result, *δ*_{2}/*d* first abruptly drops and then continually decreases.

### 3.3 Comparison Between Multi-Scale SM-FEM Model and Continuum Finite Element Analysis With Virtual Material.

*µ*m, which corresponds to the mean asperity height based on a range of ±3 standard deviations with the Gaussian distribution. The virtual material elastic modulus used with the continuum FEA model is defined as

Although Eq. (23) is based on Xiao and Sun’s virtual material method [16], the total number of asperities is eliminated; Xiao and Sun obtained the elastic modulus by integrating a single asperity elastic modulus with the number of contacting asperities. As a result, for nearly smooth contact, Xiao and Sun's expression indicates *E*_{vm} approaches infinity since the number of contacting asperities also approaches infinity. However, it is noted that $h\xafa\sigma s$ and 3*σ*_{s} become zero, and $\delta =d\xaf\sigma s$. Thus, the integral term equates to *δ*^{−1/2} and hence a single asperity radius equates to the body radius so that the elastic modulus of the virtual material is the same as the body elastic modulus.

For the comparison to FEA, a unit force boundary condition of *f* = 2.29 kN/mm is applied as shown in the earlier Fig. 3 (recall from Fig. 5(a), the elastic-plastic transition corresponds to 2.15 kN/mm). In the FEA simulation, *E*_{vm} = 6.27 GPa from Eq. (23). Constitutive behavior of the sheet is based on the industry measurement of unit force versus plastic strain in Fig. 4. Elastic modulus of the work roll is *E*_{2} = 207 GPa. Figure 8(a) indicates the FEA mesh, which contains ∼1.45 million 3-node, linear, plane-strain triangular elements resulting from three re-meshing iterations to achieve convergence in the vertical displacement at the center of the roll/sheet contact interface. Using an i7 processor with 8 GB at 60% allocated memory, the FEA simulation requires approximately 38 min. In contrast, calculation time on the same processor for the improved SM-FEM is less than 4 s in matlab 2015.

Note that the primary advantage of the large-scale continuum FEA is its capability to attain detailed contact mechanics results. Figure 8(b) shows the vertical displacement contour. The enlarged view of the contact interface shows a discontinuity transition between classic Hertzian deformation and the virtual material deformation. Although *E*_{vm} is much smaller than *E*_{2}, the virtual material layer undergoes less displacement than the roll bulk because of its small thickness. Table 2 shows a comparison of the contact interferences between the FEA, the SM-FEM with GW, and the SM-FEM with KE. Because the GW statistical asperity model is fully elastic, the SM-FEM with GW predicts a much lower asperity interference, *δ*_{a}, than either SM-FEM with KE or the continuum FEA. The SM-FEM with KE predictions of *δ*_{1}, *δ*_{2}, and *δ*_{a} are all 15% larger than for FEA. This difference may stem from the assumption that contact width and roll-bite arc are equal; in reality the Hertz contact width is slightly smaller than the roll-bite contact arc. Equations (5)–(8) shows that an underpredicted contact width leads to underpredicted stiffness, perhaps suggesting why the SM-FEM with KE contact interference is lower than FEA in Table 3. Importantly, the comparisons of *δ*_{a} in Table 2 indicate that integrating the asperity stiffnesses along the contact arc to obtain $kfa$, then using it in the equivalent work roll/sheet foundation stiffness via series combination is appropriate with the SM-FEM when the elastic-plastic KE statistical asperity theory is incorporated, but not when the fully elastic GW asperity theory is used.

## 4 3D Multi-Scale Coupling Example in 4-High Mill With Rough Upper Work Roll

The improved, roughness-capable SM-FEM is now applied to illustrate 3D multi-scale effects in a 4-high mill stand. Shown in Fig. 9 are two “true” rough surface topographies to be applied to the upper work roll, which are identified by sampling probability distributions for the respective parameters of Cases 1 and 2 in Table 1. The asperity sampling applies a Gaussian probability distribution function, $\varphi (h\xafa)$, for the asperity summit height, $h\xafa$. Table 3 lists the post-sampling statistics characterizing these true surface roughness conditions, which are matched to those of Chang et al. [39] per their experimental measurements. Case 1 gives 1.237 *µ*m average roughness (Ra), while Case 2 results in 6.877 *µ*m Ra. In accordance with the Gaussian distribution, the asperity skewness is nearly zero and the kurtosis is close to 3. All sampled asperities are assumed to be independent and identically distributed over the roll-bite contact area with the upper work roll.

For the multi-scale modeling, only the more accurate elastic-plastic KE statistical asperity theory is included to account for asperity deformation. Note that in prior macro-scale sheet thickness profile comparisons between the original SM-FEM model and industry measurements for a CVC rolling mill [30], the anti-symmetric contact interface required a full-scale contact model (no symmetry). The study illustrated flaws in other approaches that exploited symmetry between the upper and lower mill sections. Likewise, in all of the finite element studies involving roll roughness described in the introduction, the symmetric models used presume *identical* roughness profiles on the top and bottom work rolls. For the 3D study here, the efficiency of the improved SM-FEM approach allows for a full model. Table 4 lists dimensional and mechanical properties of the 4-high mill and sheet, taken from a small production mill in industry that rolls mostly type 301 stainless steel. To more easily discern the inherent multi-scale effects, the 4-high mill model here does not include any applied roll bending, shifting, or other sheet flatness control mechanisms that can otherwise be readily applied. In addition, the rolls all have uniform nominal diameter, and the sheet has a rectangular entry thickness profile. Note that this absence of profile control measures will generate large variation in the unit contact force across the sheet width, but this is advantageous for the current study in that it helps reveal roughness effects that vary as function of the contact force. For the multi-scale analysis, 652 foundation-coupled, double-beam, foundation-coupled elements (see again Fig. 3), representing 10,432 deg of freedom comprising Eq. (2), are used with four and five Newton–Raphson iterations to convergence for Cases 1 and 2, respectively. Run times on an i7 processor with 8 GB RAM are <110 min (asperity integration dominates computation);

Work roll diameter (mm) | 50.8 |

Back-up roll diameter (mm) | 304.8 |

Roll length (mm) | 304.8 |

Top back-up roll screw down (mm) | 0.0325 |

Sheet entry thickness (mm) | 0.236 |

Sheet width (mm) | 203.2 |

Roll and sheet Poisson ratio | 0.28 |

Roll elastic modulus (GPa) | 207 |

Work roll diameter (mm) | 50.8 |

Back-up roll diameter (mm) | 304.8 |

Roll length (mm) | 304.8 |

Top back-up roll screw down (mm) | 0.0325 |

Sheet entry thickness (mm) | 0.236 |

Sheet width (mm) | 203.2 |

Roll and sheet Poisson ratio | 0.28 |

Roll elastic modulus (GPa) | 207 |

Figure 10 shows the resulting unit contact force distribution as a function of axial contact length for both WR roughness cases based on a 0.35 mm vertical displacement boundary condition (“screw down”) applied at the center of each neck of the top back-up roll (BUR). The smaller roughness, Ra = 1.237 *µ*m (Case 1), results in a total WR/sheet contact force of about 573 kN (obtained by integrating the WR/sheet unit force), while the larger roughness, Ra = 6.877 *µ*m (Case 2), produces about 548 kN total WR/sheet contact force. The smoother roll therefore experiences 4.56% greater total force. Maximum and minimum values of the unit force between the sheet and the work roll (labelled in Fig. 10) are also lower for the rougher work roll of Case 2. In addition, because of the static equilibrium requirement, the larger roughness also reduces the total contact force between work roll and the back-up roll. Thus, under the same vertical displacement boundary condition, the rougher surface can be said to “absorb” more of the contact load than the smoother surface, which corroborates the 2D rough surface contact behaviors and, more importantly, for the first time, quantitatively supports that which mill operators have anecdotally reported about roughness. Aside from the contact force value differences, the larger Case 2 roughness also exhibits greater amplitude in the “high frequency” force deviations than does the smaller Case 1 roughness, as seen in the enlarged views of the contact force near the sheet center. Note also that the contact force deviation is slightly more pronounced when the unit contact force is lower (i.e., near the sheet center), but there are some specific locations where particularly large asperities contradict this trend, as in Case 1 at −6 mm. Figure 11 illustrates the results for the rolled sheet thickness profile with work roll roughness Cases 1 and 2. Due to the different length scales of the asperities and the sheet thickness, it is not possible to directly visualize the effects of the roughness topographies on the thickness profile. Therefore, Fig. 11 also includes enlarged views for Cases 1 and 2 within the sheet width range −10 to 10 mm.

Because the rougher surface absorbs more contact energy, Case 2 with Ra = 6.877 *µ*m produces an average exit thickness of 0.2331 mm, while Case 1 with Ra = 1.237 *µ*m creates a reduced average thickness of 0.2329 mm. Thus, Case 1 undergoes 0.2 *µ*m or 6.90% more average reduction and achieves about 25 kN or 4.56% more total contact force than Case 2. Accordingly, the macro-scale sheet thickness profile of Case 1 results in greater relative sheet crown,^{5} quantitatively indicated as 3.261% for Case 1 versus 3.009% for Case 2 with reference to the sheet edges, and 0.408% for Case 1 versus 0.377% for Case 2 with reference to a standard 25 mm distance from the sheet edges. In this thin sheet, the crown difference implies greater tendency to induce wavy edge type flatness defects for Case 1 than for Case 2 [30]. At the micro-scale, the larger roughness of Case 2 (Ra = 6.877 *µ*m) results in greater transfer of roughness from the work roll to the sheet than does the smaller Case 1 roughness (Ra = 1.237 *µ*m). Note also that, from comparison between the Case 2 enlarged view for −10 to 10 mm with the Case 2 enlarged views of 10 to 30 mm, it can be seen that the high frequency thickness perturbation is not as pronounced at locations farther away from the sheet center. The above shows that the SM-FEM method reveals 3D coupling effects between the micro-scale work roll roughness effects and the macro-scale roll-stack deformation. Moreover, such quantitative connection between the roll roughness and the corresponding changes in the sheet flatness quality (which are strongly related to sheet crown ratio changes) can be used to gain new insights into multi-scale geometric quality trade-offs. Note that, as mentioned earlier, while the SM-FEM and KE asperity theories have already been separately validated elsewhere in the literature, carefully designed experimental cold rolling studies with suitable process control, roll roughness design, and roll/sheet surface metrology would provide for additional understanding to help bring about relevant quality improvements in thin metal sheet rolling.

## 5 Conclusion

A new multi-scale, simplified-mixed finite element method has been developed by integration of elastic-plastic statistical asperity distributions to represent an equivalent stiffness. When applied to thin sheet cold rolling, the method reveals new quantitative, coupled effects between rough work roll surfaces and the bulk roll-stack mechanics. Comparison with a previously validated, large-scale continuum finite element model using virtual roughness material indicates that the presented approach of integrating parallel asperity stiffnesses along the contact region, and then applying this equivalent asperity stiffness integral in series combination with the bulk elastic roll stiffness and the sheet elastic-plastic stiffness, is an appropriate model representation. Results for a full-scale, 4-high rolling mill with work roll surface roughness indicate that a relatively rougher roll “absorbs” more contact load than a smoother roll, which leads to decrease in both the contact force and the “crown” thickness profile of the rolled sheet for a prescribed back-up roll displacement. In addition, while the relatively rougher roll transfers more roughness to the sheet as expected, the coupled contact force distribution, which also depends on bulk bending, shearing, and Hertzian flattening, leads to roughness transfer that varies across the sheet width. The sheet thickness profile therefore ultimately realizes the coupled effects of micro-scale asperity deformations and macro-scale roll-stack deformations. In the future work, the roll-bite lubricant frictional contact can be considered explicitly (precluding the requirement for industry observed sheet constitutive behavior), since the SM-FEM can exploit parallel connectivity between asperity stiffness and lubrication film stiffness. Indeed, Xiao et al. [42] has verified parallel connection between asperity stiffness and lubrication film stiffness using the ultrasonically measured normal contact stiffness. The new SM-FEM method described can also be extended to provide the valuable supporting predictive capability to wear, friction, and fatigue (spalling) studies.

## Footnotes

Forward slip refers to relative sliding motion wherein the sheet speed is greater than the roll peripheral speed.

Both the GW and KE statistical asperity models applied Gaussian probability distributions, which often characterize roughness well, but other probability distributions can be used as determined, for example, by Fourier transformations with experimental data.

% Crown is sheet center thickness minus sheet edge reference thickness, normalized by the edge reference thickness (expressed as %) [30].

## Acknowledgment

The authors gratefully acknowledge financial support for this work from the National Science Foundation (CMMI-1555531). Any opinions, findings, or conclusions expressed in this paper are those of the authors and do not necessarily reflect the views of the National Science Foundation.

## Conflict of Interest

There are no conflicts of interest.

## Data Availability Statement

The datasets generated and supporting the findings of this article are obtainable from the corresponding author upon reasonable request. The authors attest that all data for this study are included in the paper.

### Appendix

*i*between bodies 1 and 2 respectively, the coupled global stiffness $KG1,2,i$ has a size 24 by 24 with all six translational and rotational degree of freedom per each of the four nodes (see two Timoshenko beam elements with intervening foundation in Fig. 3):

*l*

_{i}is the

*i*th element length. For bodies 1 and 2, and shape function matrix subset term $Npq$ (for $p$, $q$ ∈ [1,2]) in the foundation stiffness matrix $KF1,2,i$ is defined as: $pq$

**represents the full shape function representing both horizontal (**

*N**w*) and vertical (

*v*) displacements for conventional Timoshenko beams. Term

*θ*in Eq. (A2) is the angle of inclination between bodies 1 and 2, where, for example, $\theta =90deg$ for a vertically oriented 4-high mill as shown Fig. 1 [30].

In Eq. (A3), *X* represents axial translation, *S* is the twist stiffness coefficient of degree of freedom *θ*_{x}, and $\u03d1$ is a property of the shape and size of the cross section. For circular sections, $\u03d1$ is the same as *J*, the polar moment of inertia.