Analysis of Continuous and Discontinuous Cases of a Contact Problem Using Analytical Method and FEM

In this paper, continuous and discontinuous cases of a contact problem for two elastic layers supported by a Winkler foundation are analyzed using both analytical method and finite element method. In the analyses, it is assumed that all surfaces are frictionless, and only compressive normal tractions can be transmitted through the contact areas. Moreover, body forces are taken into consideration only for layers. Firstly, the problem is solved analytically using theory of elasticity and integral transform techniques. Then, the finite element analysis of the problem is carried out using ANSYS software program. Initial separation distances between layers for continuous contact case and the size of the separation areas for discontinuous contact case are obtained for various dimensionless quantities using both solutions. In addition, the normalized contact pressure distributions are calculated for both cases. The analytical results are verified by comparison with finite element results. Finally, conclusions are presented.


INTRODUCTION
Since contact problems have different application areas in structural mechanics such as railways, foundation grillages, rolling mills, pavements of highway and airfield etc. (see, for example, Garrido and Lorenzana, 1998;Comez et al., 2004;Kahya et al. 2007), there has been an increasing attention on the contact problems. There is a large body of literature associated with contact problems both analytically (Keer et al., 1972;Weitsman, 1972;Ratwani and Erdogan, 1973;Gladwell, 1976 to Layer 2. In the analyses, it is assumed that all surfaces are frictionless, and only normal tractions can be transmitted across the contact surfaces. Since every material has a weight in nature and including body forces approaches the solution to reality, body forces of the elastic layers are taken into account. Including body forces also converts the type of the problem from receding contact to continuous or discontinuous contact problem. Due to symmetry about y-axis, it is sufficient to consider only one-half of the problem geometry. Where applicable, the germane quantities are reckoned per unit length in the z-direction. Analytical solution of the problem is obtained using theory of elasticity and integral transform technique. In the absence of body forces, the stress and the displacement components for Layer 1 and Layer 2 may be obtained as (Birinci and Erdol, 2003)       y y ih i i i i 0 2 u x, y A B y e C D y e sin x d where subscript h indicates the case without body forces. u=u(x,y) and v=v(x,y) are the x-and ycomponents of the displacement vector, respectively. x (x, y)  , y (x, y)  and xy (x, y)  are the stress components of the layers. i  is an elastic constant and equals . The subscripts 1 and 2 refer to the Layer 1 and Layer 2, respectively. A i , B i , C i and D i (i=1, 2) are unknown coefficients which will be determined from boundary conditions of the problem. For the case in which body forces exist, the particular parts of y  stress components for Layer 1 and Layer 2 are obtained as (Birinci and Erdol, 2003)   p(x) gh    is sufficiently small, then the contact between elastic layers, y=h 1 , 0<x<∞ , will be continuous, and A i , B i , C i and D i (i=1,2) must be determined from the following boundary conditions: where 0 k is the stiffness of the Winkler foundation. By making use of boundary conditions (9-16), unknown A i , B i , C i and D i (i=1,2) constants can be determined. By substituting these constants into Eq. (8), after some simple manipulations, one may obtain the contact pressure along the interface 1 y h  as (Birinci and Erdol, 2003) In this study, three different loading cases are considered. For each loading case, j P (j=1,2,3) can be determined as in Eq. (17) to zero after replacing , the expression in which cr x (initial separation distance) will be obtained can be written as where j f equals to corresponding loading j P (j=1,2,3) given in Eqs. (22-24), cr  is critical load factor and it can be written as After substituting j f ( ,a / h)  (j=1,2,3) values into Eq. (25) and solving the integral numerically, x cr values, for which initial separation takes place may be obtained for various dimensionless quantities.

Discontinuous Contact Case
where b and c are unknowns, and functions of  . Note that the unknown function (x)  in Eq.
(35) may be replaced by Utilizing the boundary conditions (27-33) and (35), the new constants A i , B i , C i and D i (i=1,2) which appear in Eqs.
(1-5) may be obtained in terms of the unknown function (x)  . Then, Eq. (34) gives the following singular integral equation which (x)  is the unknown function: One may notice that because of the smooth contact at the end points b and c, the function (x)  equals zero at the ends and the index of the integral Eq. (38) is equal to -1 (Erdogan and Gupta, 1972). The consistency condition of the integral Eq. (38) defined as Defining the following dimensionless quantities, and substituting Eqs. (43-45) into the integral Eq. (38), the single-value condition (37), and the consistency condition (42), the following equations are obtained: To insure smooth contact at the end points of the separation area, let Using the appropriate Gauss-Chebyshev integration formula (Erdogan and Gupta), Eqs. (46) and (47) become Eqs. (53) and (54) give (n+2)

FINITE ELEMENT ANALYSIS OF THE PROBLEM
In this section, the problem has been studied by the finite element method (FEM) using a commercial package program ANSYS. ANSYS Multiphysics is a powerful interactive environment for modeling and solving all kinds of scientific and engineering problems based on partial differential equations (PDEs). To solve the PDEs, ANSYS Multiphysics uses the proven finite element method (FEM). The software runs the finite element analysis together with adaptive meshing and error control using a variety of numerical solvers (Biswas and Banerjee, 2013). The problem is considered as a two-dimensional contact problem and the material of the layers are assumed elastic and isotropic. The physical system under consideration exhibits symmetry in geometry, material properties and loading. Due to the symmetry of the problem, only one half of the geometry of the problem is to be modeled. The geometry and the applied load are shown schematically in Fig. 2 and Finite Element model of the problem before analysis is shown in Fig. 3. In the study, two dimensional solid elements (PLANE 183) are used to model the layers. The PLANE 183 element is defined by six nodes having two degrees of freedom at each node: translations in the nodal x and y directions. The element may be used as a plane element (plane stress, plane strain and generalized plane strain). Winkler foundation is modelled by linear spring element (COMBIN 14). The COMBIN14 element or the longitudinal element spring-damper option is a uniaxial tension-compression element with up to two degrees of freedom at each node: translations in the nodal x and y directions. No bending or torsion is considered. If damping is neglected, the spring element will simply represent the linear Winkler model of one parameter which is the simple model to represent soil (Al-Azzawi et al., 2010).
The contact region is meshed by Surface-to-Surface CONTA172 and TARGE169 contact elements. CONTA172 is used to represent that of the mechanical contact analysis. The target surface, defined by TARGE169, was therefore used to represent 2-D ''target'' surfaces for the associated contact elements CONTA172. Plane strain finite elements are used for the meshing of the entire geometry. Frictionless surface-to-surface contact elements are used to model the interaction between the contact surfaces and Augmented Lagrangian method is used as the contact algorithm.

COMPARISON AND VERIFICATION OF ANALYTICAL AND FINITE ELEMENT APPROACHES
In this section, by using the described methods in the previous sections, initial separation distances between layers for continuous contact case, the size of the separation areas for discontinuous contact case, and the normalized contact pressure distributions for both contact cases are given for various dimensionless quantities such as h / h and  . Additionally, the results obtained using FEM have been compared and verified with analytical results.

Latin American Journal of Solids and Structures 12 (2015) 1771-1789
The effects of 2 1 /   and 0 1 k k /   on the initial separation distance for continuous contact case are presented in Table 1. It demonstrates that the initial separation distance increases with increasing of 2 1 /   . This means that as Layer 2 becomes more rigid compared to Layer 1, the initial separation occurs at a more distant point from the axis of symmetry. Moreover, increasing stiffness of the Winkler foundation ( 0 1 k k /   ) reduces the value of initial separation distance cr x and the initial separation occurs at a point closer to the axis of symmetry. As the change of initial separation distance cr x dwindles away and cr x approaches to a constant value.  Figure 4 illustrates that the effect of stiffness of the Winkler foundation on the normalized contact pressure between layers for continuous contact case. The results of this figure reveal that maximum normalized contact pressure occurs at the axis of symmetry and maximum normalized contact pressure and initial separation distance decrease with increasing of 0 1 k k /   . Table 2 depicts that the effect of different loading cases on the initial separation distance cr x . Examination of Table 2 indicates that maximum initial separation distance is obtained from loading case 1 [ The normalized contact pressure distributions for continuous contact case ( cr    ) and discontinuous ( cr    ) contact case are depicted in Fig. 6. As seen in Fig. 6, there are three regions in the case of discontinuous contact. These are continuous contact region, separation zone and also continuous contact region where the effect of the external load decreases and disappears infinitely. Deformation shapes after finite element analysis related to Fig. 6 for continuous and discontinuous cases are shown in Fig. 7. Finally, the examination of all tables and figures shows that finite element solution indicates a good agreement with analytical solution.

CONCLUSIONS
The main purpose of this paper is to compare and verify results obtained using analytical method and FEM for a contact problem in the continuous and discontinuous contact cases. So, the problem is solved analytically using linear elasticity theory and integral transform technique. Then, initial finite element model of the problem is created by ANSYS software and finite element analysis is performed. According to this study, the followings can be deducted:  Initial separation distance increases as (  Initial separation distance increases with increasing of 2 1 /   . This means that as Layer 2 becomes more rigid as against Layer 1, the initial separation occurs at more distant point from the axis of symmetry.  Maximum normalized contact pressure occurs at the axis of symmetry and maximum normalized contact pressure and initial separation distance decrease with increasing of 0 1 k k /   .  Maximum initial separation distance is obtained from loading case 1 [ 0 p(x) p  ].
 Maximum normalized contact pressure between layers is obtained from loading case 3  Increasing the value of 1 h / h results in an decrease of separation distances (b and c) between the layers for both loading cases  There are three regions in the case of discontinuous contact. These are continuous contact region, separation zone and also continuous contact region where the effect of the external load decreases and disappears infinitely  It is observed that the results obtained using FEM is in a good agreement with the analytical results.