Finite Element Modelling of a Pair of Flexible Elements Contact Preloaded and Externally Loaded with an Arbitrary Force

In the paper, the numerical tests of the contact of a pair of flexible elements were presented. The non-linear phenomena occurring in the contact zone of two elements have been considered separately for the assembly and operational condition of the joint. The physical and mathematical models of the contact joint created using the finite element method were proposed. The contact joint is treated as a system composed of the subsystems: two of them are joined elements and the third one is the contact layer between the elements. The calculations results for the selected contact joint were given.


INTRODUCTION
In contact joints, the elements with contact surfaces that have previously undergone some mechanical treatment are most often joined [8]. Such surfaces do not adhere to each other over the entire nominal contact surface, but only on a small part of it, in separate randomly arranged micro-areas [37,47,48]. The actual contact surface depends on the shape and dimensions of micro-unevenness and the free space between them, the value of the clamping force, the mechanical properties of the contacting materials as well as the type of the contact joint and is a negligible fraction of the nominal contact surface [11,14,32,41,46,50]. As a result, significant deformations occur within the individual micro-contacts in the area of the nominal surface, which have a crucial impact on the distribution of interaction forces and mutual displacements of the contacting bodies.
The issues at the interface of two bodies are usually considered in two directions: normal and tangential [34,43]. For small values of pressure and deformations, for which the contact layer shows high plasticity, the physical relationships between normal displacements of points of this layer and normal pressure, as well as between tangential displacements and shear pressure are characterized by non-linear trajectories [19,26,29]. At the same time, both normal contact displacements and tangential displacements exhibit the elastic-plastic nature at the first loading, while at subsequent loading and unloading -elastic.
In order to describe the phenomena occurring in the contact zone of a pair of elements, an appropriate model of the rough surface of the contacting bodies and their contact should be adopted. There are two approaches to this issue found in the literature: analytical [6,22,23,35,36,39,40,49] and experimental [4,25,27,30]. Because the analytical research for physical relationships describing the contact displacements of interacting rough surfaces most often lead to very complex formulas, in this paper it was decided to use the relationships obtained in an experimental way to write up the mechanical characteristics of the contact. It was assumed that in the case of normal pressure and displacements, these characteristics can be presented with a good approximation by means of an exponential function, as in the papers [15,28]. However, in the case of shear pressure and tangential displacements, these characteristics can be described by classic non-linear curves, as in the paper [19].
Considering the issues of contact mechanics on a macroscopic scale, many contact details can be omitted. It is then purposeful to isolate the contact zone from two contacting bodies. The transition from the contact of these bodies to the physical model consists in cutting out of the surface layers with unevenness and inserting in their place a conventional third body in the form of a thin layer [45]. The conventional third body defined in this way has appropriate properties, including the features of the contact between real bodies, resulting directly from the empirical research. In addition, in the adopted contact model, there are, independently of each other, both the own stiffnesses of the connected bodies and the stiffness of the contact, determined by the stiffness assigned to the conventional layer separating these bodies. It also enables the use of numerical methods for modelling and calculating contact joints, among which the finite element method (FEM) is the most popular.
The analysis of the contact phenomena at the interface of a pair of joined elements using finite element systems has been carried out so far in many publications and there are many models of contact joints that have been described in them. These models can be divided into the following three groups: − models of the joint, in which the contact layer between the joined elements is not taken into consideration [1], − linear contact models with constant values of the contact stiffness coefficients over the entire nominal contact surface, commonly available in commercial finite element analysis systems [9,10,12,20,21,24,31,33,38,42], − non-linear contact models [5,7,13,51]. In the above-mentioned models, the possibility of taking into account the variable contact stiffness coefficients for each of the contact elements on the contact surface is not considered, which would be useful in the case of the contact analysis of the elements joined in a connection model with the experimentally defined stiffness characteristics of the contact layer between them [16]. In this paper, a method of modelling of the contact layer, by means of which it is possible to model the elastic properties of individual elements of the layer, is presented. The study was conducted on the example of a pair of flexible elements contact preloaded and externally loaded with an arbitrary force. The calculations were carried out in a Midas NFX 2017 R1 finite element system, but using some special procedures that supported this computer system.

PHYSICAL MODEL OF THE CONTACT JOINT
A general structure of the contact joint model for the j-th elementary contact area is illustrated in Fig. 1. The model is based on a fragment of a flexible flange element that is fastened to a flexible support, for example, as in the case of multi-bolted connections [18]. In the normal direction to the surfaces of these elements, a non-linear spring is introduced, similar to the Winkler contact model [2]. In addition, a pair of fictional springs are appended to the model, placed perpendicular to the spring in the normal direction.
The experimentally determined relationships between normal displacements of points of the contact layer u n and normal pressure p n , as well as between tangential displacements u t and shear pressure p t , can be presented by means of waveforms, as shown in Fig. 2. A description of these characteristics is provided in [3]. On its basis, the following formulas can be derived for the springs stiffness in the j-th elementary area of the contact model ( Fig. 1): (1) (2) where: c zj is the contact stiffness coefficient in the normal direction of the j-th element of the contact model, c xj , c yj are contact stiffness coefficients in the tangential direction of the j-th element of the contact model, p n is normal contact pressure, A j is the j-th elementary contact area, S j is the cross-sectional area of the tangential spring of the j-th element of the contact model, and c, m, R, S are the experimentally determined constants [3]. Structuring of the contact layer model is run in the following four steps: − division of the contact surface between the flange and the support into elementary contact areas, − assignation mesh nodes in the centres of gravity of elementary contact areas, − insertion of normal and tangential springs at the nodes identified in the previous step, − creation of the 2D finite element mesh on the contact surface. On the basis of the 2D finite element mesh on the contact surface between the flange and the support, a uniform 3D finite element mesh for all the volume of the flange and the support is generated.
The equilibrium equation for the contact joint model can be written in the form [18]: where: K is the stiffness matrix, q is the displacements vector, and p is the loads vector. In the assembly state, the contact joint is loaded only by normal forces. Then, as a result of solving Eq. (3), normal displacements of nonlinear springs u nj are calculated. By multiplying them by the stiffness of these springs determined with formula (1), one obtains the reactions in the springs R j . The linearization of the characteristics of non-linear springs is performed using the secant method (Fig. 3). The calculations are carried out in an iterative process in which the following condition is checked [17]: (4) where: R j is the reaction in the j-th non-linear spring, R' j is the reaction in the j-th nonlinear spring after the linearization, and e 1 is the relative error. In the operational state, the preloaded contact joint is loaded with an external force. The linearization of the characteristics of non-linear springs is carried out in a similar way as in the previous stage of calculations, starting from their operating points determined at that stage and verifying the condition analogous to (4). In addition, a second iterative process is introduced here, in which the following condition is checked [44]: (5) where: F j is the resultant tangential force for the j-th elementary contact area, T j is the friction force on the j-th elementary contact area, and e 2 is the relative error.

EXAMPLE CALCULATIONS
According to the presented method, the computations of a contact joint shown in Fig. 4 were realized. The dimensions of each of the joined elements were: 15 mm x 30 mm x 9 mm. The preload of every non-linear spring F nj was equivalent of 1 kN. All the springs were preloaded simultaneously. After the assembly process, the joint was loaded by a force F e with coordinates (200 N, 300 N, 250 N) as shown in Fig. 4.
The values of the constants in Eqs. (1) and (2) were taken as for the surfaces machined by peripheral grinding and collected in Table 1.
The obtained calculation results included, among other things, normal contact pressure distributions for the two tested conditions. They are shown in Fig. 5. On the basis of the graphs, it can be seen that after the external load was applied to the system, the values of normal pressure between the joined elements were reduced, but no looseness occurred at their contact.
As a result of calculations, the shear pressure values for individual elementary areas of the contact layer were also obtained. They are given in comparison with limit shear pressure in Table 2. On their basis it can be concluded that there was no slippage between the joined elements.

CONCLUSIONS
In the paper, a method for conducting the computations of a contact joint between a pair of  flexible elements with use of its experimental normal and shear characteristics was presented. The method is related to the assembly stage and the operational stage of the joint loaded with an arbitrary force, and the calculations were prepared taking into account the variability of stiffness coefficients for each of the elements of the contact layer. This was possible by creating special calculation procedures and using them in conjunction with the calculations carried out in the finite element system. On the basis of the calculations results, it can be stated that after the external load is applied to the system, the distribution of normal pressure between the joined elements changes, but at the appropriate preload of the system there is no looseness at their contact. Similarly, due to the appropriate system preload, no slippage occurs between the joined elements.
With the developed method, it is possible to model the multi-bolted connections treated as a system under the conditions of any external loads. The phenomena occurring in the contact joints between the elements joined in a multi-bolted connection are similar to the ones occurring in other structural connections, so the method can also be used in the analysis of these connections.