# Fractal Geometry-Based Hypergeometric Time Series Solution to the Hereditary Thermal Creep Model for the Contact of Rough Surfaces Using the Kelvin-Voigt Medium

Authors: Osama M. Abuzeid; Anas N. Al-Rabadi; Hashem S. Alkhaldi
Journal: Mathematical Problems in Engineering
 (2010)

---

## Abstract

This paper aims at constructing a continuous hereditary creep model for the thermoviscoelastic contact of a rough punch and a smooth surface of a rigid half-space. The used model considers the rough surface as a function of the applied load and temperatures. The material of the rough punch surface is assumed to behave as Kelvin-Voigt viscoelastic material. Such a model uses elastic springs and viscous dashpots in parallel. The fractal-based punch surface is modelled using a deterministic Cantor structure. An asymptotic power law, deduced using approximate iterative relations, is used to express the punch surface creep which is a time-dependent inelastic deformation. The suggested law utilized the hypergeometric time series to relate the variables of creep as a function of remote forces, body temperatures, and time. The model is valid when the approach of punch surface and half space is in the order of the size of the surface roughness. The closed-form results are obtained for selected values of the system parameters; the fractal surface roughness and various material properties. The obtained results show good agreement with published experimental results, and the methodology can be further extended to other structures such as the Kelvin-Voigt medium within electronic circuits and systems.


---

## Body

## 1. Introduction

Surface topography plays a significant role in tribology, that is, in problems of friction, wear, lubrication, and contact [1]. Therefore, the problem of analysis of rough surfaces attracts the attention of engineers and applied mathematicians. Historically, the following engineering parameters, statistical in nature, were used for the characterization of surface roughness: (1) the root mean square of the heights, σ, (2) the root mean square of slopes, σm2, (3) and the root mean square of curvatures, σk2. However, it was realized that the topography of engineered surfaces is too complex to be described completely by a few statistical parameters. Thus, it was found that roughness has a multiscale nature and requires sophisticated mathematical techniques for its description.First attempts to model the distribution of heights of surface asperities utilize the classical random field theory which assumed that the functions of surface model are differentiable. In particular, this implies that limiting values forσm2 and σk2 should exist as the sample interval tends to “0” [1]. However, it turned out that such limiting behavior is in contradiction with the results of advanced investigations of surfaces. For example, the exponential behavior of the autocorrelation function implies that the engineering parameters should tend to infinity rather than to constant values when the sampling interval is infinitely reduced [2]. Furthermore, it was shown that the profiles of a large number of both natural and artificial surfaces have the following form of the spectral density function G(ω)~1/ωυ where υ≈2, and ω is the spatial frequency (see (1.13)). It follows from this that all wavelengths are equally represented in the profile and that there exists no characteristic scale; in other words, after arbitrary magnification roughness looks like before. Moreover, it was found that the values of engineering parameters depend on the measurement scale, that is, these parameters are scale dependent [2, 3].The fractal approach was introduced as an attempt to give a scale-invariant characterization of surface topography. The idea of fractality of roughness was experimentally verified on real surfaces as well as when applied to mathematically simulated profiles [4]. Figure 1 shows a picture of popular fractals, that is, the middle-third Cantor set, the von Koch curve, graph of the Weierstrass-Mandelbrot function C in the range 0≤x≤3 (p=1.5 and γ=0.5 where p and γ are two numerical parameters (see (1.17)) where the trend of the function is ~x2, and trajectories of a fractional Brownian process for different Hurst index H and fractal dimension D [5–7]. The index H is the key parameter of the fractal surface which describes the smoothness of the surface.Common fractals: (a) the middle-third Cantor set, (b) the von Koch curve, (c) the Weierstrass-Mandelbrot functionC in the range 0≤x≤3 (p=1.5 and γ=0.5), and (d) trajectories of a fractional Brownian process for different H and D.
(a)(b)(c)(d)Evidently, roughness of the surface of a body has a great influence on stress fields that arise when two deformable bodies are pressed together. Analysis of the effect of roughness on the contact interaction of solids has attracted wide attention [8].One of the most popular models for studying contact of rough bodies is the Greenwood and Williamson (GW) model based on the use of the Hertz theory [9], where it is important to mention that (GW) model is a nonscale-invariant [10]. Currently, the development of models of contact between nominally flat fractal rough surfaces presented for the Cantor profile is an active area of research [11]. Various contact problems utilizing Cantor profile were considered [12–17]. All these models consider the one-level Cantor profile. It is argued that such profile is simple for analytical analysis. However, it has a minor drawback: all asperities of the profile have one-level character, while all real roughness has a hierarchical structure [17].It is accepted that fractal dimension is not a compressive geometric parameter that could characterize alone the behavior of contacting rough bodies [1]. Moreover, the employment of the fractal approach in the study of surfaces has several drawbacks. The proposed model can be both fractal and nonfractal depending on values of the structural parameters. Regardless of this, the model profile remains rough and possesses certain self-affine properties. The iterative regular construction of the profile allows us to analyze its prestructures, that is, prefractals, of arbitrary generation.In this introduction, important and relevant definitions and methods that are attributed to fractal geometry with the application to the modeling of rough surfaces will be fully presented. Furthermore, the important differences between mathematical and utilized physical fractals will be explicitly highlighted.


### 1.1. Mathematical Definition of Fractals

Mandelbrot stated that a set in a metric space is called a fractal set if the Hausdorff-Besicovitch dimension of the set is greater than its topological dimension [18]. Let X be a compact metric space and O be the totality of open balls in X. The Hausdorff s-measure of a subset S⊂X which is defined for s≥0 as the following limit:(1.1)mH(S,s)=limσ→0infG∈O{∑V∈G(dimV)s:S⊆⋃V∈GV,diamV≤δ}.
Here G is finite or denumerable subset of O. It was proven that there exists a value s0 such that
(1.2)mH(S,s)={∞,fors<s0,0,fors>s0.The Hausdorff dimension of the setS, denoted by dimHS, is the number s0 such that (1.2) holds. Unfortunately, the calculation of the Hausdorff dimension of mathematical objects often demands a lot of effort. Even to find some estimations of the dimension, it is necessary to overcome a number of rather complex mathematical difficulties [19]. This issue called for the use of other definitions of dimension which are useful in applied mathematics for the characterization of fractal objects. One such alternative is the box dimension [1]. The analytical calculation of the box dimension is usually easier since the corresponding definition of this dimension involves coverings by spheres of equal radii.LetE be the Euclidean dimension of the space in which a set S is embedded. For δ>0, let N(δ) be the smallest number of E-dimensional balls or cubes of diameter d needed to cover the set S. The box counting dimension or box dimension, denoted by dimBS, can be defined if the following limit exists:
(1.3)dimBS=limδ→0+logN(δ)-logδ.It can be proven thatdimBS does not change if one takes N(δ) to be either the smallest number of δ-cubes that cover S: the number of δ-mesh cubes that intersect S; or the smallest number of sets of diameter at most d that cover S; or the largest number of disjoint δ-balls with centers in S. Unfortunately, the box dimension is not always equal to the Hausdorff dimension. For example, the set S={0,1,1/2,…,1/n…} has unequal values for the Hausdorff and box dimensions for dimHS≠dimBS=1/2. However, it can be proven that dimHS≤dimBS.As a simple alternative to the Hausdorff measure, we can introduce thes-measure ms of a set as the following limit:
(1.4)ms(S)=limδ→0+N(δ)δs
and define the box dimension as the value s=D such that ms(S) has a jump from 0 to ∞ similar to the behavior of mH(S,s) in (1.2).On the other hand, the difficulties involved with calculating the Hausdorff dimension are the reason for the opinion that the Hausdorff dimension is not generally used in applications in the study of fractal and non fractal curves that are originated in other sciences such as in biology, engineering, physics, quantum physics and computing [20–27].


### 1.2. Physical Concept of Fractals

Evidently, it is impossible to carry out the scaling procedure for any real physical object down to infinitely small scales. Hence, the mathematical concept of the Hausdorff measure is applicable only to mathematical models of objects rather than to the objects themselves and, of course, the Hausdorff dimension cannot be obtained by experimental procedures. In this sense there are no actual fractal objects in nature. For physical objects, the box dimension cannot be calculated analytically but it is estimated by experimental or numerical calculations. However, various errors can arise during such numerical calculations. There is no canonical definition of physical fractals and there are numerous methods for the practical estimation of the fractal dimension of an object. The cluster fractal dimension is taken as the first example of a physical fractal dimension definition.Let a whole cluster be imagined as consisting of elementary parts of the sizeδ* [1]. An object can be modeled as a fractal cluster with dimension D when the model considers scales R such that δ*<R<Δ*, where δ* and Δ*are the upper and lower cutoffs for the fractal representation. To get the value D of the dimension, the considered region is discretized into cubes with side length δ*. Then the smallest number of E-dimensional cubes needed to cover the cluster (N(δ*)) is counted. One says that the cluster is fractal if the numbers N(δ*) satisfy the so-called number-radius relation for different sizes of the considered region of the cluster R as follows:
(1.5)N(δ*)≈(Rδ*)D,δ*<R<Δ*.The value ofD is estimated as the slope of linear growth of ln(N(δ*)) plotted against ln(R). The power D is usually called the cluster dimension or mass dimension.In literature, various methods were utilized to estimate the fractal dimension of a physical object. However, the notion of fractal dimension is not well-defined in that the relative value does depend on the approach used. Indeed, only for the mathematical box dimension of a fractal setS it is proven that dimBS is the same when using various specific schemes of covering [28], while for physical fractals the estimations of the fractal dimension inevitably involve various techniques, distinct scale ranges, and various computation rules. Therefore, the obtained values can differ strongly and it is unlikely that they could be fruitfully compared for distinct objects. Thus, even in the case of physical objects of a similar nature, it would be wrong to consider the fractal dimension of these objects as their specific property without referring to the estimation technique involved.


### 1.3. Self-Similarity and Self-Affinity of Surfaces

Let us recall that a one-to-one mappingM of a plane π onto a plane π′ is called a similarity mapping with coefficient λ>0, or simply a similarity, when the following property holds: if {A,B} are any two points of π, and {A′,B′} are their images under M, then |A′B′|=λ|AB| [29]. It is known that any similarity transformation of a plane is a homogeneous (isotropic) dilation of coordinates {x′=λx,z′=λz} up to a rotation and translation. A set S is called statistically self-similar if under homogeneous scaling with the coefficient λ, where 1>λ>0, it is identical from the statistical point of view to the set S′=λS.In practice, it is impossible to verify that all statistical moments of the two distributions are identical. Frequently, a setS is said to be self-similar if only a few moments do not change under scaling [30]. A one-to-one mapping M of a plane π onto a plane π′ is called an affine mapping, if the images of any three collinear points are collinear in turn [29]. In general, an affine transformation of a plane may be given in any coordinate system as a nondegenerative linear transformation. In practical studies of rough surfaces, one often considers a particular affine mapping, with anisotropic scaling, that is given coordinate wise by x′=λx and z′=λHz. Here z is a graph of a surface profile and H is some scaling exponent.One says that a fractal is self-affine if it is invariant from the statistical point of view under quasihomogeneous (anisotropic) scaling. It is possible to show that usually a quasihomogeneous transformation is a particular case of Lipschitz homeomorphism [1, 17]. The Hausdorff dimension of a set S does not change under the action of the Lipschitz homeomorphism L as follows:
(1.6)dimHS=dimHL(S).The ideas of self-similarity and self-affinity are very popular in studying surface roughness because experimental investigations show that usually profiles of vertical sections of real surfaces are statistically similar to themselves under repeatedly magnifications; however, the profiles should be scaled differently in the direction of nominal surface plane and in the vertical direction. The self-affine fractals were used in a number of papers as a tool for description of rough surfaces [3, 31, 32]. Two standard examples of self-affine fractals are the trace of the fractional Brownian motion and the Weierstrass function. The former is a statistical fractal while the latter is a deterministic fractal.


### 1.4. Brownian Surfaces and Random Fractals

Fractional Brownian processes are widely used in creating computer-generated surfaces, in particular landscapes. For example, a profile can be constructed as a graph of 1D(fBm)VH(x) of index H, where x is taken as the time and z is the random variable of the single valued function VH(x) with the following property:
(1.7)〈[VH(x+δ)-VH(x)]2〉~δ2H,0<H<1,
where 〈〉 denotes averaging over the ensemble, and H is the Hurst index. The scaling behavior of the different traces, VH(x), is characterized by a particular H which relate the typical change in Δz(x), where z(x)=VH(x), is the trace of thefBm, and the change in the spatial coordinate Δx by the simple scaling law [30, 33, 34]:
(1.8)Δz(x)~ΔxH.It is known that, with probability equal to “1”, the following holds [28]:
(1.9)dimHVH(x)=dimBVH(x)=2-H.The autocorrelation function is one of the main tools for studying statistical models of rough surfaces. The autocorrelation functionR(δ) of the profile is(1.10)R(δ)=limT→∞12T∫-TT[z(x+δ)-z̅][z(x)-z̅]dx=〈[z(x+δ)-z̅][z(x)-z̅]〉
or
(1.11)R(δ)=limT→∞12T∫-TTz(x+δ)z(x)dx-(z̅)2,
where z̅ is the average value of the profile function z(x).Another tool for the characterization of surfaces is the spectral density functionG(ω) which is the Fourier transform of R(δ):
(1.12)G(ω)=2π∫0∞R(δ)cosωδdδ,R(δ)=2π∫0∞G(ω)cosωδdω.In general, it is accepted in fractional Brownian motion that [14]:(i)
if the autocorrelation functionR(δ) of the profile z(x) satisfies R(0)-R(δ)~δ2(2-s), then it is reasonable to expect that the box dimension of the graph z(x) is equal to s, note that one can find R(0)-R(δ)~δ2H for the fBm defined by (1.7).(ii)
if the profilez(x) has spectral density:
(1.13)G(ω)~1ωυ,
then it is reasonable to expect that the box dimension of the graph z(x) is equal to (5-υ)/2 [1]. The above conclusions are valid for mathematical models of the profile, for which the relation 2(5-s)=υ-1orυ=5-2s holds. The exponent υ varies typically between 0 and 2. Usually, it is assumed that the same conclusions concerning the box dimension are valid for physical fractals as well. It is shown that real surfaces approximately satisfy the property in (1.13) in wide range of scales [35]. The moments mn of the spectral density G(ω) provide a useful description of the surface roughness:
(1.14)mn=∫ω0∞ωnG(ω)dω,
where ω0=2π/λ0 is the wave number corresponding to the profile length λ0. It is possible to show that m0 is the variance of heights (rms height) of the surface, m1 is the variance of slopes (rms slope) and m3 is the variance of curvatures (rms curvature) [36].


### 1.5. Weierstrass-Type Functions and Modeling of Rough Surfaces

A number of researchers have used the Weierstrass-type functions for fractal modeling of surface roughness [3, 31, 32] and fractal modeling applications such as in quantum computing [20, 21]. The real Weierstrass-type function can be defined as:(1.15)W(x;p)=∑n=0∞p-γnh(pnx),p>1,0<γ<1,
where h is a bounded Höِlder function of order greater than β. The following complex generalization of the W(x;p)was considered:(1.16)W̃(x;p)=∑n=-∞∞p-γn[(1-eipnt)eiΦn],p>1,0<γ<1,
where Φn are arbitrary phases [29].The Weierstrass-type functions are continuous everywhere and differentiable nowhere. In addition, their graphs are curves whose fractal dimension exceeds one. Fractal properties of these functions including the Weierstrass-Mandelbrot (WM) functionC and the Takagi-Hopson function T:(1.17)C(x;p)=∑n=-∞∞p-γn(1-cospnx),p>1,0<γ<1,(1.18)T(x;p)=∑n=-∞∞p-γn|pnx-[pnx+12]|,p>1,0<γ<1,
have been studied in numerous papers [2, 14, 29, 31]. By direct calculations, one may obtain:
(1.19)|W̃(x+δ;p)-W̃(x;p)|~δγ,
which is similar to the behavior of (1.7) of fractional Brownian motion. The box dimension of the Weierstrass function graphs is D=2-γ and it is believed that their Hausdorff dimension is the same [28, 37]. Currently, the only known bounds for the Hausdorff dimensions are D-(c/logp)≤dimHgraphC≤D, provided that p is large and constant c is large enough [19]. It is possible to calculate the spectral density of the WM function W̃(x;p) as follows:(1.20)G(ω)=∑n=-∞∞δ(ω-pn)p2(2-D)n,
where δ is the Dirac delta function. Some arguments for approximating this discrete spectral density by a continuous spectral density G̃(ω)~1/ω5-2D, whose exponent (5-2D) is in agreement with (1.13) with respect to the box dimension were suggested [2]. The following truncated WM function
(1.21)W̃1(x;p)=A(D-1)∑n=n1∞p(D-2)ncos2πpnx
is often used for fractal characterization of the surface topography [3, 31, 32]. Here n1 is an integer number, which corresponds to the low cutoff frequency of the profile, and A is the so-called characteristic length scale of the profile. The number n1 depends on the length L of the sample and is given by pn1=1/L and the parameter A determines the position of the spectral density along the logG axis. It was stated that both parameters A and D of the function W1 are scale-invariant characteristics of the roughness. However, the extensive experimental studies of this fractal characterization model showed that the values of parameters A and D are not unique and depend on instruments or resolution of a given instrument [1].Evidently, the functionC(x;p) is not homogeneous. Nevertheless, it exhibits the property C(pkx;p)=pkγC(x;p), with k∈Z where Z is the set of all integers which looks similar to the definition of a homogeneous function hd of degree d, that is, hd(λx)=λdhd(x) for λ>0.Thus, the graph of the functionC(x;p) near any point x0 is repeated in scaling form near all points pkx0,k∈Z. This scaling (self-affine) property was often attributed to fractal features of the graph. However, this discrete scaling property is the main property of the so-called parametric-homogeneous (PH) functions introduced [1, 17] which strictly satisfy the equation bd(pkx;p)=pkdbd(x;p), with k∈Z where d is degree of homogeneity. As examples of 1-dimensional fractal PH-curves we can consider the graphs of functions b1 and b2 with degrees d=1 and d=2, respectively:
(1.22)b0(x;p)=x-γC(x;p),b1(x;p)=xb0(x;p),b2(x;p)=x2b0(x;p).Because of (1.6), these functions have the same Hausdorff dimension as the WM function C(x;p) whose box-dimension is D.Another consequence is that the WM functionC(x;p), with C(x;p)~x2-D can be used only as an example of fractal profile and it cannot be considered as the general fractal functional model for simulations of the rough surface profiles. The assumption that the WM function represents the general fractal properties of rough profiles can lead to wrong conclusions concerning surface roughness parameters and their distributions.The solution to the problem of mechanical contact between elastically deforming solids was obtained by Hertz [8]. Subsequently, several approaches were used to analyze the contact interaction between the soft layer and the indenting object surface [38–42]. These methods are based upon Radok's technique of replacing the elastic constants in the elastic solution by the corresponding integral or differential operators, which appear in the stress-strain relations for linear viscoelastic materials. Furthermore, these studies assumed that the surfaces of contacting solids are smooth, excluding from consideration all real solids, which have a certain degree of roughness and waviness regardless of how fine their finish is [43].Various models for the approach of the fractal punches were considered [11–16]. In the previously cited works, different constitutive relations were considered: (1) linear elastic material [11], (2) rigid-perfectly plastic material [13], (3) elastic-perfectly plastic material [12], (4) linear viscoelastic creep model via Maxwell medium [14], (5) linear viscoelastic creep model via standard linear solid (SLS) material [15], and (6) linear thermoviscoelastic relaxation model via Maxwell medium [16].The objective of this work is to introduce an alternative approach, using fractal geometry, to study the deformation of a viscoelastic surface as a function of the force applied and the bulk temperature. In this model friction force effect is assumed to be negligible. The development of the fractal model of the rough surface is carried out using fractional Brownian motion in conjunction with Cantor set. The Radok's technique [44] is then used to derive the thermoviscoelastic model from the corresponding elastic model. The main contribution of this work is a mathematical model for the time-dependent-creep of a rough surface (cf. (6.7)). This model relates the creep to time, temperature, external applied load, fractal dimension of the rough surface, and various material properties.Section2 presents the fractal model, where the Cantor structure is built and its fractal dimension is presented. Section 3 presents the discrete and continuous elastic model. In Sections 4 and 5, the effect of temperature on the viscoelastic behavior is presented and the Arrhenius's relation is introduced. In Section 6, the elastic viscoelastic correspondence is presented which consists of replacing the elastic constant in the elastic solution by the corresponding integral or differential operators from the viscoelastic stress-strain relations. Also, in Section 6 a new continuous model for the creep contact of a thermovisco-elastic punch is presented. In Section 7, the results obtained from the new model is presentd and compared with an experimental results obtained from literature. In Section 8, conclusions and future work are presented.


## 1.1. Mathematical Definition of Fractals

Mandelbrot stated that a set in a metric space is called a fractal set if the Hausdorff-Besicovitch dimension of the set is greater than its topological dimension [18]. Let X be a compact metric space and O be the totality of open balls in X. The Hausdorff s-measure of a subset S⊂X which is defined for s≥0 as the following limit:(1.1)mH(S,s)=limσ→0infG∈O{∑V∈G(dimV)s:S⊆⋃V∈GV,diamV≤δ}.
Here G is finite or denumerable subset of O. It was proven that there exists a value s0 such that
(1.2)mH(S,s)={∞,fors<s0,0,fors>s0.The Hausdorff dimension of the setS, denoted by dimHS, is the number s0 such that (1.2) holds. Unfortunately, the calculation of the Hausdorff dimension of mathematical objects often demands a lot of effort. Even to find some estimations of the dimension, it is necessary to overcome a number of rather complex mathematical difficulties [19]. This issue called for the use of other definitions of dimension which are useful in applied mathematics for the characterization of fractal objects. One such alternative is the box dimension [1]. The analytical calculation of the box dimension is usually easier since the corresponding definition of this dimension involves coverings by spheres of equal radii.LetE be the Euclidean dimension of the space in which a set S is embedded. For δ>0, let N(δ) be the smallest number of E-dimensional balls or cubes of diameter d needed to cover the set S. The box counting dimension or box dimension, denoted by dimBS, can be defined if the following limit exists:
(1.3)dimBS=limδ→0+logN(δ)-logδ.It can be proven thatdimBS does not change if one takes N(δ) to be either the smallest number of δ-cubes that cover S: the number of δ-mesh cubes that intersect S; or the smallest number of sets of diameter at most d that cover S; or the largest number of disjoint δ-balls with centers in S. Unfortunately, the box dimension is not always equal to the Hausdorff dimension. For example, the set S={0,1,1/2,…,1/n…} has unequal values for the Hausdorff and box dimensions for dimHS≠dimBS=1/2. However, it can be proven that dimHS≤dimBS.As a simple alternative to the Hausdorff measure, we can introduce thes-measure ms of a set as the following limit:
(1.4)ms(S)=limδ→0+N(δ)δs
and define the box dimension as the value s=D such that ms(S) has a jump from 0 to ∞ similar to the behavior of mH(S,s) in (1.2).On the other hand, the difficulties involved with calculating the Hausdorff dimension are the reason for the opinion that the Hausdorff dimension is not generally used in applications in the study of fractal and non fractal curves that are originated in other sciences such as in biology, engineering, physics, quantum physics and computing [20–27].


## 1.2. Physical Concept of Fractals

Evidently, it is impossible to carry out the scaling procedure for any real physical object down to infinitely small scales. Hence, the mathematical concept of the Hausdorff measure is applicable only to mathematical models of objects rather than to the objects themselves and, of course, the Hausdorff dimension cannot be obtained by experimental procedures. In this sense there are no actual fractal objects in nature. For physical objects, the box dimension cannot be calculated analytically but it is estimated by experimental or numerical calculations. However, various errors can arise during such numerical calculations. There is no canonical definition of physical fractals and there are numerous methods for the practical estimation of the fractal dimension of an object. The cluster fractal dimension is taken as the first example of a physical fractal dimension definition.Let a whole cluster be imagined as consisting of elementary parts of the sizeδ* [1]. An object can be modeled as a fractal cluster with dimension D when the model considers scales R such that δ*<R<Δ*, where δ* and Δ*are the upper and lower cutoffs for the fractal representation. To get the value D of the dimension, the considered region is discretized into cubes with side length δ*. Then the smallest number of E-dimensional cubes needed to cover the cluster (N(δ*)) is counted. One says that the cluster is fractal if the numbers N(δ*) satisfy the so-called number-radius relation for different sizes of the considered region of the cluster R as follows:
(1.5)N(δ*)≈(Rδ*)D,δ*<R<Δ*.The value ofD is estimated as the slope of linear growth of ln(N(δ*)) plotted against ln(R). The power D is usually called the cluster dimension or mass dimension.In literature, various methods were utilized to estimate the fractal dimension of a physical object. However, the notion of fractal dimension is not well-defined in that the relative value does depend on the approach used. Indeed, only for the mathematical box dimension of a fractal setS it is proven that dimBS is the same when using various specific schemes of covering [28], while for physical fractals the estimations of the fractal dimension inevitably involve various techniques, distinct scale ranges, and various computation rules. Therefore, the obtained values can differ strongly and it is unlikely that they could be fruitfully compared for distinct objects. Thus, even in the case of physical objects of a similar nature, it would be wrong to consider the fractal dimension of these objects as their specific property without referring to the estimation technique involved.


## 1.3. Self-Similarity and Self-Affinity of Surfaces

Let us recall that a one-to-one mappingM of a plane π onto a plane π′ is called a similarity mapping with coefficient λ>0, or simply a similarity, when the following property holds: if {A,B} are any two points of π, and {A′,B′} are their images under M, then |A′B′|=λ|AB| [29]. It is known that any similarity transformation of a plane is a homogeneous (isotropic) dilation of coordinates {x′=λx,z′=λz} up to a rotation and translation. A set S is called statistically self-similar if under homogeneous scaling with the coefficient λ, where 1>λ>0, it is identical from the statistical point of view to the set S′=λS.In practice, it is impossible to verify that all statistical moments of the two distributions are identical. Frequently, a setS is said to be self-similar if only a few moments do not change under scaling [30]. A one-to-one mapping M of a plane π onto a plane π′ is called an affine mapping, if the images of any three collinear points are collinear in turn [29]. In general, an affine transformation of a plane may be given in any coordinate system as a nondegenerative linear transformation. In practical studies of rough surfaces, one often considers a particular affine mapping, with anisotropic scaling, that is given coordinate wise by x′=λx and z′=λHz. Here z is a graph of a surface profile and H is some scaling exponent.One says that a fractal is self-affine if it is invariant from the statistical point of view under quasihomogeneous (anisotropic) scaling. It is possible to show that usually a quasihomogeneous transformation is a particular case of Lipschitz homeomorphism [1, 17]. The Hausdorff dimension of a set S does not change under the action of the Lipschitz homeomorphism L as follows:
(1.6)dimHS=dimHL(S).The ideas of self-similarity and self-affinity are very popular in studying surface roughness because experimental investigations show that usually profiles of vertical sections of real surfaces are statistically similar to themselves under repeatedly magnifications; however, the profiles should be scaled differently in the direction of nominal surface plane and in the vertical direction. The self-affine fractals were used in a number of papers as a tool for description of rough surfaces [3, 31, 32]. Two standard examples of self-affine fractals are the trace of the fractional Brownian motion and the Weierstrass function. The former is a statistical fractal while the latter is a deterministic fractal.


## 1.4. Brownian Surfaces and Random Fractals

Fractional Brownian processes are widely used in creating computer-generated surfaces, in particular landscapes. For example, a profile can be constructed as a graph of 1D(fBm)VH(x) of index H, where x is taken as the time and z is the random variable of the single valued function VH(x) with the following property:
(1.7)〈[VH(x+δ)-VH(x)]2〉~δ2H,0<H<1,
where 〈〉 denotes averaging over the ensemble, and H is the Hurst index. The scaling behavior of the different traces, VH(x), is characterized by a particular H which relate the typical change in Δz(x), where z(x)=VH(x), is the trace of thefBm, and the change in the spatial coordinate Δx by the simple scaling law [30, 33, 34]:
(1.8)Δz(x)~ΔxH.It is known that, with probability equal to “1”, the following holds [28]:
(1.9)dimHVH(x)=dimBVH(x)=2-H.The autocorrelation function is one of the main tools for studying statistical models of rough surfaces. The autocorrelation functionR(δ) of the profile is(1.10)R(δ)=limT→∞12T∫-TT[z(x+δ)-z̅][z(x)-z̅]dx=〈[z(x+δ)-z̅][z(x)-z̅]〉
or
(1.11)R(δ)=limT→∞12T∫-TTz(x+δ)z(x)dx-(z̅)2,
where z̅ is the average value of the profile function z(x).Another tool for the characterization of surfaces is the spectral density functionG(ω) which is the Fourier transform of R(δ):
(1.12)G(ω)=2π∫0∞R(δ)cosωδdδ,R(δ)=2π∫0∞G(ω)cosωδdω.In general, it is accepted in fractional Brownian motion that [14]:(i)
if the autocorrelation functionR(δ) of the profile z(x) satisfies R(0)-R(δ)~δ2(2-s), then it is reasonable to expect that the box dimension of the graph z(x) is equal to s, note that one can find R(0)-R(δ)~δ2H for the fBm defined by (1.7).(ii)
if the profilez(x) has spectral density:
(1.13)G(ω)~1ωυ,
then it is reasonable to expect that the box dimension of the graph z(x) is equal to (5-υ)/2 [1]. The above conclusions are valid for mathematical models of the profile, for which the relation 2(5-s)=υ-1orυ=5-2s holds. The exponent υ varies typically between 0 and 2. Usually, it is assumed that the same conclusions concerning the box dimension are valid for physical fractals as well. It is shown that real surfaces approximately satisfy the property in (1.13) in wide range of scales [35]. The moments mn of the spectral density G(ω) provide a useful description of the surface roughness:
(1.14)mn=∫ω0∞ωnG(ω)dω,
where ω0=2π/λ0 is the wave number corresponding to the profile length λ0. It is possible to show that m0 is the variance of heights (rms height) of the surface, m1 is the variance of slopes (rms slope) and m3 is the variance of curvatures (rms curvature) [36].


## 1.5. Weierstrass-Type Functions and Modeling of Rough Surfaces

A number of researchers have used the Weierstrass-type functions for fractal modeling of surface roughness [3, 31, 32] and fractal modeling applications such as in quantum computing [20, 21]. The real Weierstrass-type function can be defined as:(1.15)W(x;p)=∑n=0∞p-γnh(pnx),p>1,0<γ<1,
where h is a bounded Höِlder function of order greater than β. The following complex generalization of the W(x;p)was considered:(1.16)W̃(x;p)=∑n=-∞∞p-γn[(1-eipnt)eiΦn],p>1,0<γ<1,
where Φn are arbitrary phases [29].The Weierstrass-type functions are continuous everywhere and differentiable nowhere. In addition, their graphs are curves whose fractal dimension exceeds one. Fractal properties of these functions including the Weierstrass-Mandelbrot (WM) functionC and the Takagi-Hopson function T:(1.17)C(x;p)=∑n=-∞∞p-γn(1-cospnx),p>1,0<γ<1,(1.18)T(x;p)=∑n=-∞∞p-γn|pnx-[pnx+12]|,p>1,0<γ<1,
have been studied in numerous papers [2, 14, 29, 31]. By direct calculations, one may obtain:
(1.19)|W̃(x+δ;p)-W̃(x;p)|~δγ,
which is similar to the behavior of (1.7) of fractional Brownian motion. The box dimension of the Weierstrass function graphs is D=2-γ and it is believed that their Hausdorff dimension is the same [28, 37]. Currently, the only known bounds for the Hausdorff dimensions are D-(c/logp)≤dimHgraphC≤D, provided that p is large and constant c is large enough [19]. It is possible to calculate the spectral density of the WM function W̃(x;p) as follows:(1.20)G(ω)=∑n=-∞∞δ(ω-pn)p2(2-D)n,
where δ is the Dirac delta function. Some arguments for approximating this discrete spectral density by a continuous spectral density G̃(ω)~1/ω5-2D, whose exponent (5-2D) is in agreement with (1.13) with respect to the box dimension were suggested [2]. The following truncated WM function
(1.21)W̃1(x;p)=A(D-1)∑n=n1∞p(D-2)ncos2πpnx
is often used for fractal characterization of the surface topography [3, 31, 32]. Here n1 is an integer number, which corresponds to the low cutoff frequency of the profile, and A is the so-called characteristic length scale of the profile. The number n1 depends on the length L of the sample and is given by pn1=1/L and the parameter A determines the position of the spectral density along the logG axis. It was stated that both parameters A and D of the function W1 are scale-invariant characteristics of the roughness. However, the extensive experimental studies of this fractal characterization model showed that the values of parameters A and D are not unique and depend on instruments or resolution of a given instrument [1].Evidently, the functionC(x;p) is not homogeneous. Nevertheless, it exhibits the property C(pkx;p)=pkγC(x;p), with k∈Z where Z is the set of all integers which looks similar to the definition of a homogeneous function hd of degree d, that is, hd(λx)=λdhd(x) for λ>0.Thus, the graph of the functionC(x;p) near any point x0 is repeated in scaling form near all points pkx0,k∈Z. This scaling (self-affine) property was often attributed to fractal features of the graph. However, this discrete scaling property is the main property of the so-called parametric-homogeneous (PH) functions introduced [1, 17] which strictly satisfy the equation bd(pkx;p)=pkdbd(x;p), with k∈Z where d is degree of homogeneity. As examples of 1-dimensional fractal PH-curves we can consider the graphs of functions b1 and b2 with degrees d=1 and d=2, respectively:
(1.22)b0(x;p)=x-γC(x;p),b1(x;p)=xb0(x;p),b2(x;p)=x2b0(x;p).Because of (1.6), these functions have the same Hausdorff dimension as the WM function C(x;p) whose box-dimension is D.Another consequence is that the WM functionC(x;p), with C(x;p)~x2-D can be used only as an example of fractal profile and it cannot be considered as the general fractal functional model for simulations of the rough surface profiles. The assumption that the WM function represents the general fractal properties of rough profiles can lead to wrong conclusions concerning surface roughness parameters and their distributions.The solution to the problem of mechanical contact between elastically deforming solids was obtained by Hertz [8]. Subsequently, several approaches were used to analyze the contact interaction between the soft layer and the indenting object surface [38–42]. These methods are based upon Radok's technique of replacing the elastic constants in the elastic solution by the corresponding integral or differential operators, which appear in the stress-strain relations for linear viscoelastic materials. Furthermore, these studies assumed that the surfaces of contacting solids are smooth, excluding from consideration all real solids, which have a certain degree of roughness and waviness regardless of how fine their finish is [43].Various models for the approach of the fractal punches were considered [11–16]. In the previously cited works, different constitutive relations were considered: (1) linear elastic material [11], (2) rigid-perfectly plastic material [13], (3) elastic-perfectly plastic material [12], (4) linear viscoelastic creep model via Maxwell medium [14], (5) linear viscoelastic creep model via standard linear solid (SLS) material [15], and (6) linear thermoviscoelastic relaxation model via Maxwell medium [16].The objective of this work is to introduce an alternative approach, using fractal geometry, to study the deformation of a viscoelastic surface as a function of the force applied and the bulk temperature. In this model friction force effect is assumed to be negligible. The development of the fractal model of the rough surface is carried out using fractional Brownian motion in conjunction with Cantor set. The Radok's technique [44] is then used to derive the thermoviscoelastic model from the corresponding elastic model. The main contribution of this work is a mathematical model for the time-dependent-creep of a rough surface (cf. (6.7)). This model relates the creep to time, temperature, external applied load, fractal dimension of the rough surface, and various material properties.Section2 presents the fractal model, where the Cantor structure is built and its fractal dimension is presented. Section 3 presents the discrete and continuous elastic model. In Sections 4 and 5, the effect of temperature on the viscoelastic behavior is presented and the Arrhenius's relation is introduced. In Section 6, the elastic viscoelastic correspondence is presented which consists of replacing the elastic constant in the elastic solution by the corresponding integral or differential operators from the viscoelastic stress-strain relations. Also, in Section 6 a new continuous model for the creep contact of a thermovisco-elastic punch is presented. In Section 7, the results obtained from the new model is presentd and compared with an experimental results obtained from literature. In Section 8, conclusions and future work are presented.


## 2. Fractal Model

The surface profile of the punch, in contact with a rigid half-space, will be constructed on the basis of Cantor set [11]. The contacting surface is constructed by joining the segments obtained at successive stages of the construction of a Cantor set to one another, Figure 2, where L0 correspond to the profile nominal length, and h0 is equal to the twice rms height of the roughness.Figure 2
The fractal middle-third Cantor structure, whereE0 is the initiator step, E1, E2, and E3 are the other generated step of cantor structure, L's are the lengths of the E's steps, h's are the heights of E's steps, and F is the applied load.At each stage of profile construction, the middle section of each initial segment is discarded so that the total length of the remaining segments is1/a times the length of the initial segment, where a>1. The depth of the recesses (measured from the last step) at the (i+1)th construction step of the fractal surface is 1/b times less than the depth of the ith step, where b>1. From this it can easily be shown that the horizontal length and recess depth of the (i+1)th step are, respectively(2.1)Li+1=a-1Li=a-(i+1)L0,(2.2)hi+1=b-1hi=b-(i+1)h0,
where it is assumed that the surface is smooth in a direction perpendicular to the plane of the page. This restriction is not expected to have a significant effect since it is possible to construct a fractal Cantor surface perpendicular to the plane of the page [11].At theith generation, the Cantor structure contains N=2i segments, each of length δi=(2a)-iL0 [11]. The profile of the surface in Figure 2 can be considered as a certain graph of a step function.It can be seen that, during an iterative step in constructing the surface, scaling in the horizontal direction isΔχi+1=(2a)-1Δχi, while in the vertical direction, the corresponding fluctuations Δzi at the ith generation can be defined by considering the probability of obtaining the value, zi=b-ih0.The fluctuationΔzi at the ith generation can be obtained by assuming the Δzi scales as the expected value ziP(zi) in which Δzi∝ziP(zi) [8], where P(zi) is the probability of obtaining the value zi, that is, P(zi)=(Li-Li+1)/L0, and it is found that P(zi)=a-i(1-1/a).Thus, the expected value of the fluctuation at the(i+1)th generation is related to the expected value of the fluctuation at the ith generation through zi+1P(zi+1)=(ab)-iziP(zi). Hence Δzi+1=(ab)-iΔzi, and thus (Δzi+1/Δzi)=(Δχi+1/Δχi)2-D, from which the self-affine fractal dimension for the contour of the Cantor structure is derived as:
(2.3)D=1+ln2ln2a-lnbln2a=1+Dc-lnbln2afor(1<D<2),
where Dc is the fractal dimension of the Cantor set (0<Dc<1). Equation (2.3) will be used in the next section in the development of the approach-force model.


## 3. The Continuous Elastic Model

Qualitatively, two size scales are manifested in the contact problem [8]:(1)
the bulk scale, for which the elastic compression would be calculated by the Hertz theory and its limitations,(2)
the roughness scale, where the asperities act like a compliant layer on the surface, and so all the deformations are limited in a surface layer which represents all the asperities;bh0 in Figure 2, and their deformation is assumed to be linear elastic [45].In this paper, the approach of the punch of Cantor structure surface and lengthL0 will be considered. It is to be noted that the obtained relation may be applied for all problems with surfaces having the same fractal dimension. The contact between two rough surfaces can be modeled as the contact of an effective surface with a rigid flat surface [10]. Hence, a solution for the deformation of an equivalent surface generated using the Cantor structure can be modified for the problem at hand. The bodies treated in this work will be assumed to be isotropic and homogeneous, and obey linear force-displacement laws. The yield strength σy, the modulus of elasticity E, and coefficient of thermal expansionα, are all assumed to be independent of temperature.Furthermore, it is assumed, with reference to Figure2, that there exists a series of one-dimensional elastic bars, distributed in a way such that the distance from the initiator step E0 to the generated step E3 is indicated by h0, from E1 to E3 is indicated by h1, from E2 to E3 is indicated by h2, and so forth, [14–16]. By letting F3 be the force required to compress E3 until E2, F2 be the force require to compress E3 and E2 until E1, and F1 be the force required to compress E3, E2, and E1 until E0, and assuming unit depth, one obtains the following discrete force-displacement relations:
(3.1)F3=h2k3,F2=(h1-h2)k2+h1k3,F1=(h0-h1)k1+(h0-h2)k2+h0k3,
where ki=ELi/bh0 is the stiffness of the ith step, E is the modulus of elasticity of the material used, bh0 could be understood from Figure 2, and hi and Li can be calculated using (2.1) and (2.2), respectively.It is to be noted that thermal forces and deflections may arise in heated body either because of a nonuniform temperature distribution, or external constraints, or as a combination of these causes. The problem is assumed to be a steady state one with no internal heat source.Next, by lettingΔFi+1=Fi-Fi+1, then from (3.1) one can conclude the general equation for any number of steps as follows:
(3.2)ΔFi+1=EL0b(b-1)b-ia-i.In order to find a recursive relation as in (3.2) for the approach u,one lets Fi+1 be the limit force for protrusion of the (i+1)th generation. It is assumed that when the limit load is reached, the punch approaches a distance Δui+1, equals to the difference between the heights protrusion of ith and (i+1)th generations. Consequently the second generation E2 deflects a distance u2=h1, E1 deflects a distance u1=h0, and E0 deflects a distance u0=bh0, so Δu1=u0-u1=h0(b-1), and Δu2=h0b-1(b-1). Accordingly:(3.3)Δui+1=h0(b-1)b-i.The above-mentioned assumptions are sufficient to determine the dependence of the limit loadF on the approach u. The effects of the remote load and the bulk temperature will be first studied separately and then superimposed. Using the fact that, when the limit load increases from Fi+1 to  Fi, the punch is approached by an amount Δui+1, and by utilizing (3.2) and (3.3), the remote load effect is given by the discrete force displacement relation:
(3.4)ΔFi+1Δui+1=EL0bh0a-i.Asi→∞, (3.4) yields the following asymptotic behavior for the strain ε:
(3.5)ε=(bχEFL0)1/(χ+1),
where the strain ε could be defined as u/h0 and χ=lna/lnb [14–16]. It is to be noted that, the effect of the applied external load should not exceed a limiting yield load Fy, where Fy=σyL0.Since the interest of this work is to consider the viscoelastic behavior, the principle of correspondence [44] will be used in the next section to obtain a viscoelastic model corresponding to the elastic model presented in (3.5).


## 4. Effect of Temperature

Temperature has a dramatic influence on rates of viscoelastic response, and in practical work it is often necessary to adjust a viscoelastic analysis for varying temperature. This strong dependence of temperature can also be useful in experimental characterization, for example, if a viscoelastic transition occurs too quickly at room temperature, for easy measurement, the experimenter can lower the temperature to slow things down and vice versa.In some viscoelastic materials, the relation between time and temperature can be described by correspondingly simple models. Such materials are termed “thermorheologically simple” [46]. For such simple materials, the effect of lowering the temperature is simply to shift the viscoelastic response (plotted against log-time) to the right without change in shape. This is equivalent to increasing the relaxation time τ, without changing the relaxation modulus.A time-temperature shift factor(aT) can be defined as the horizontal shift that must be applied to a response curve, measured at an arbitrary temperature T in order to move it to the curve measured at some reference temperature Tref.If the creep time obeys an Arrhenius relation, the shift factor can be shown to be [47]:(4.1)logaT=Q2.303R(1T-1T0),
where Q is the activation energy (J/mol), R is the gas constant (J/mol·K), and T is the temperature (K).


## 5. Arrhenius Relation

The creep properties of materials are usually described by reference to the dependence of the creep (ε) on the applied stress, time and temperature, which may be written as:
(5.1)ε=f(σ,t,T).One way to simplify this function is to make it to be separable into three functions of stress, time and temperature as follows [38]:
(5.2)ε=f1(σ)f2(t)f3(T).Temperature has a significant effect on the creep of materials. In some steels, it is found that the temperature has a pronounced effect than the strain rate [48].Thermal forces and deflections may arise in a heated body either because of a nonuniform temperature distribution, or external constraints. The problem is assumed to be a steady state one with no internal heat. Arrhenius relation is a simple, but remarkably accurate, formula for the temperature dependence, where according to Arrhenius law, the temperature dependence is given as [49]:
(5.3)f3(T)=Bexp(QRT),
where B is a constant, Q is the activation energy and R is the ideal gas constant. The functions f1(σ) and f2(t) will be included in (6.4) in Section 6.It is clear that at the reference temperature,T0, the temperature function f3(T) is equal to unity, and the creep will be a function of the stress and time, as it was shown in (5.1). From (5.3) at the reference temperature T0, it can be easily shown that the value of the constant B will be:
(5.4)B=exp(-QRT0).


## 6. Elastic Viscoelastic Correspondence

The simplest approach to this problem consists of replacing the elastic constant in the elastic solution by the corresponding integral or differential operators from the viscoelastic stress-strain relations [44]. This approach can be applied to the contact problem provided that the loading program is such that the contact area is increasing throughout [38].A Kelvin-Voigt linear model is employed to describe the viscoelastic behavior of the compliant layer. Such a model is an arrangement of spring and dashpot in parallel, as shown in Figure3, in which ηand Ee are the Newtonian viscosity and the elastic modulus, respectively. The time-dependent force-displacement relation could be written in the operator form using the linear differential time operator ∂t≡∂/∂t as shown in (6.1) [50]:
(6.1)ε=σEe+η∂t.Figure 3
The linear Kelvin-Voigt model, whereηis the Newtonian viscosity, Ee is the elastic modulus, F is the applied load and u indicates the points to be displaced.It is clear that simple the constant of proportionality between stress and strain does no longer exist. The viscoelastic operator corresponding to the modulusE in (3.5) could be written as [50]:
(6.2)1E→1Ee11+τ∂t,
where (τ≡η/Ee) is a characteristics parameter with units of time called the retardation time.Creep test is a widely used standard test, wherein a forceP0 is suddenly applied at time t=0 on the viscoelastic model and then maintained constant thereafter, while measuring the approach as a function of time. The applied force can be expressed as a function of time with the aid of the unit step function[U(t)]. Thus
(6.3)F=P0[U(t)].By substituting (6.1)–(6.3) in (3.5) one obtains(6.4)f1(σ)f2(t)=(u(t)h0)1=(bχτEeP0L0)1/(χ+1)t1/(1+χ)Γ((2+χ)/(1+χ))1F1[11+χ;2+χ1+χ;-tτ]Tables of Laplace transforms [51] were utilized to obtain (6.4), where  Γ  represent the gamma function, and 1F1[c;d;x] is the Kummer's confluent hypergeometric function which could be expressed as [52]
(6.5)1F1[c;d;x]=1+cdx+c(c+1)x2d(d+1)2!+c(c+1)(c+2)x3d(d+1)(d+2)3!+⋯
or
(6.6)1F1[c;d;x]=∑n=0∞(c)n(d)nxnn!.By substituting the value of the constantB from (5.4) into (5.3), and substituting (3.5) and (5.3) into (5.2), the creep stain as a function of stress, time, and temperature, is obtained as follows:
(6.7)ε=f1(σ)f2(t)f3(T)=u(t)h0=(bχτEeP0L0)1/(χ+1)t1/(1+χ)Γ((2+χ)/(1+χ))1F1[11+χ;2+χ1+χ;-tτ]×exp(QR(1T-1T0)),
where u(t)/h0 is the approach, P0/L0 is the applied stress per unit depth, T is the bulk temperature, T0 is the reference temperature, τ≡η/Ee is the retardation time, ηis the Newtonian viscosity, Ee is the elastic modulus, Q is the activation energy, and R is the ideal gas constant. The effect of the fractal dimension D appears through the constant χ which combines the two scaling parameters {a,b}, that is, χ=lna/lnb (cf. (2.3)). The results obtained using this analytical model are presented and discussed in the next section.


## 7. Results and Discussion

A new continuous model for the creep contact of a thermovisco-elastic punch has been presented in (6.7). The model presents an approximate closed form solution for the approach u(t)/h0 of the fractal surface as a function of the applied load P0/L0. In order to use (6.7), values for the system parameters Ee,  η,  a and b are needed. The constant parameters a and b, in this equation, characterize the Cantor structure of the rough surface and are related through the fractal dimension D of the rough surface.For the Kelvin-Voigt model to be capable of describing the experimental results of various viscoelastic materials, the viscous coefficients, η(GPa·sec), and the elastic modulus Ee (GPa), should be selected properly. In this study the modulus Ee is selected to be 130 GPa and the retardation time, τ≡η/Ee, is selected to be of the order “1” [53]. The value of n is taken large enough for the result to converge to an accepted accuracy; it is assumed that n=550.The Cantor structure, shown in Figure2, is built from the middle-third Cantor set where the parameter a=1.5 is held fixed, giving a Cantor set dimension Dc=0.63093 [18]. Two different values of the parameter b (b=1.155 and b=1.29) yield two different dimensions of the Cantor structure; D=1.5 and D=1.4, respectively, which are used to verify the proposed analytical model. It was pointed out [11] that only for b≤2 the profile of the surface of the contacting body is fractal.In order to examine the validity of the presented model, results obtained using this model for selected values of the system parameters were compared with those obtained experimentally by the results in [12, 54]. These results displayed the approach-force relation between a flat rough surface and an ideally smooth and rigid counter surface. The specimens used in these experiments were made of carbon steel (0.45 percent carbon). Their surface roughness resulted from different finishing processes; face turning, grinding, and bead-blasting. The experiments were conducted using MATLAB for a wide range of the nominal load, up to 600 MPa. The error in the experimental measurements was determined to be approximately ∓0.5μm for the approach, and ∓5 MPa for the load [54]. Furthermore, the fractal dimension, D, of a ground stainless steel surfaces is D=1.5 [12]. The value of h0, which corresponds to twice rms height for the ground surface, is taken as 6.6μm [54], and the value D=1.5 [37] is used to calculate the value of the parameter b,using (2.3), for a fixed value of the parameter a.Figure4 shows a set of numerical creep data obtained by applying the model which was presented in (6.7). In Figure 4, the strain u/h0 is plotted versus the nondimensional time record t/τ for ΔT=0, with two different fractal dimensions D=1.5 and D=1.4, and two applied stresses σ1=300 MPa and σ2=500 MPa. As might be expected, higher strain rates occur for higher temperatures for a constant stress.Figure 4
Non-dimensional time-strain curves for Kelvin-Voigt model whenΔT=0, with D=1.5 and D=1.4, and two applied Stresses σ1=300 MPa and σ2=500 MPa.Figure5 shows the strain u/h0 versus the nondimensional time record t/τ. It shows, also, a set of numerical creep data obtained applying (6.7) for three temperature differences ΔT=0, ΔT=100 and ΔT=200 with two different fractal dimensions D=1.5 and D=1.4, and constant applied stress σ1=300 MPa. The behavior is also similar to the known typical creep curves; it indicates that higher strain rates occur for higher temperatures.Figure 5
Non-dimensional time-strain curves for Kelvin-Voigt model when for three temperature differencesΔT=0, 100 and 200, with D=1.5 and D=1.4, and constant applied stress σ=300 MPa.Figure6 is similar to Figure 5 in all of its aspects except for the applied stress, where in this case σ=500 MPa.Figure 6
Non-dimensional time-strain curves for Kelvin-Voigt model when for three temperature differencesΔT=0, 100 and 200, with dimensions D=1.5 and D=1.4, and constant applied stress σ=500 MPa.Figures4–6 show that higher strain rates result for higher fractal dimensions, which could be explained from the definition of the fractal dimension itself. Lower fractal dimensions means less roughness and consequently the bulk material dominates, while higher fractal dimensions means excessive roughness and consequently the asperities deformation dominate.Figure7 presents the isochronous stress-strain creep curves accompanied by experimental results available in the literature [54] for nondimensional time durations t/τ≈44 (which is held constant), fractal dimension D=1.4, and for three temperature differences ΔT=0, ΔT=100 and ΔT=200.Figure 7
Experimental results and analytical isochronous Stress-Strain curves for constantt/τ≈44 for three temperature differences ΔT=0, ΔT=100 and ΔT=200, with D=1.4.Figure8 also presents the isochronous stress-strain creep curves accompanied by experimental results available in the literature [54] which are conducted at the room temperature. The analytical results are shown for the nondimensional time durations t/τ≈42 (which is held constant), fractal dimension D=1.5, and for temperature differences ΔT=0, ΔT=100 and ΔT=200.Figure 8
Experimental results and analytical isochronous Stress-Strain curves for constantt/τ≈42 for three temperature differences ΔT=0, ΔT=100 and ΔT=200, with D=1.5.Figures7 and 8 show good agreement between the presented proposed model and the experimental data. For D=1.4 the nondimensional time duration, t/τ, required to get an agreement between the experimental results and the isochronous curves is about 44.2, while it is about 44 for the fractal dimension D=1.5.It is clear that the relatively longer duration of agreement occurs for lower fractal dimensions which could be attributed to the same reason mentioned above, that is, lower fractal dimensions means less roughness and consequently the bulk material dominates, while higher fractal dimensions means excessive roughness and consequently the asperities deformation dominate. The mathematical model also shows instability whent/τ exceeds 45.


## 8. Conclusions and Future Work

As well known, creep analysis is a nonlinear time-dependent phenomenon. The model which is modified in this work presents a solution to the thermal creep-contact of rough surfaces as a hypergeometric time series. Fractal geometry, via Cantor set, is utilized to model roughness of the creeping contact surfaces. The results obtained by this model turn out not to be too far from reality, since tests, at room temperatures, on the actual contact area of ground metal surfaces show that they contain sets of parallel ragged-edged scratches of different depths.Since the construction of the Cantor structure is periodic in its nature, it undergoes the same construction procedure at each hierarchical level producing contact areas that are all of the same size. Therefore, the presented analytical model provides an approximate, not exact, simulation of the approach of the viscoelastic rough surfaces, where the presented model shows a fairly good agreement with the available experimental results. As linearity is an inherent assumption, it is not expected from this model to be able to describe exactly the real material behavior; roughness and deformation. A nonlinear thermoviscoelastic stress-strain relation is required for the reproduction of real material behavior. It is also clear that the specific character of the fractal model has little effect on the asymptotic behavior of the process, and the fractal dimensionD which provides a measure of the rate at which a surface is changing is of most importance. The solution obtained in this work provides further insight into the effect that surface structure has on the deformation process, and it also provides indications of the effect that different surface forming processes may have on the subsequent surface deformation. Furthermore, in the averaged sense, the Cantor structure model appears to provide fairly reasonable results.For future work, it is intended to extend the methodology which is used in this paper for the application to the contact-surfaces within electronic and electrical devices and circuits such as resistors, capacitors and inductors that arise in electronic manufacturing systems. Furthermore, it is intended to further investigate the applications of fractals in the emerging quantum computing domain.



---

*Source: 652306-2010-07-11.xml*