The main stage in the computation of structures made of elastomers by the displacement-based finite element method (FEM) is a derivation of the stiffness matrices. Their properties define the existence, stability and convergence of the FEM solutions, as well as the effectiveness of the method. At the same time, based on the FEM displacement-based methods for modelling elastomers often have a slow convergence, especially for massive bodies having complex curvilinear forms. The slow convergence is typical for the cases, where the approximation of shifts cannot be accurately modelled by considering displacements of the finite elements as a rigid whole. To solve this problem, the paper proposes variational relations for the tetrahedral finite element, which are developed on the base of the moment scheme of the FEM (MS-FEM). A numerical experiment shows that the results obtained by application of the MS-FEM outperform the solutions obtained by application of the conventional FEM scheme.
Rubber and rubber-like materials (elastomers) are widely used in various industries. Due to the widespread use, there is a need for the new and effective methods of design and computation of the elastomers structures [1].
Elastomer structures are used in the various fields of industry and normally have complex geometry, e.g. plates, discs, couplings, dampers, suspension brackets, bearings, rings, hinges, including the rubber seals for movable and fixed joints. Numerous examples of the use of thin-layer rubber-metal elements in technology are presented in the paper [2], as well as in [3; 4]. As a rule, modern design structures include various materials, where elastomers have undergone heavy loads. Therefore, the design should take into account the rigidity, strength, heat generation and cyclic deformation.
Rubber-like materials have a specific structure, based on multiple repeating of identical layers, where the length is tens of thousands times greater than the transverse dimensions. This causes the flexibility of the molecular chains, which leads to the appearance of highly elastic properties [5]. As a result, elastomers have the following features:
1. An ability under the influence of external (constant or changing in time) loading to experience significant deformations without destruction.
2. Weak compressibility of the elastomer, which causes specific methods of the computations in comparison with conventional materials.
3. In the case of deformation and highly elastic state, the equilibrium between force and displacement is established over a period of time.
Although physical studies of elastomers started over a hundred years ago, the computation methods of the stress and strain states are still under development. This is caused by the complexity of nonlinear differential equations that describe the solution [6].
The implementation of numerical methods needs improvement of the computation schemes and development of new effective algorithms. One of the problems arouses in the boundary values of the Poisson coefficient, which leads to the degeneration of the matrix of the system of equations. There are different directions to find a solution.
The first direction is characterized by the analysis of nonlinear problems (geometric nonlinearity, large deformations) based on the combination of the penalty methods and FEM [7-9]. The second direction is based on a reduced integration [10], where the displacement fields and values responsible for the poor compressibility of the elastomers are approximated by various functions. As a rule, the degree of the polynomial for the second function is less (on one unit) as for the first function. The third direction develops mixed variational principles, in which independent displacement and strain fields or displacement and stress fields were approximated [8]. The basis of this method is the principle of variation of the components of displacement and average stress. This variation principle is widely used in the design of the structures, made of elastomers.
A. Sakharov proposed so-called moment scheme of the finite element method [11]. The approximating function is decomposed into a Taylor series and the members that respond to the displacements and dummy shifts in deformation are subsequently rejected. This allows users to take into account the basic properties of rigid displacements for isoparametric and curvilinear finite elements of isotropic elastic bodies. However, the exact equations of deformation and displacement are replaced by approximate ones.
In the mechanics of elastomers, there is also a problem of choosing an optimum computational scheme, based on the specific methods of computational mathematics [17, 18]. To check the effectiveness of a particular computational scheme, the obtained intermediate and final results should be investigated for compliance with the mechanical sense of the problem [19, 20]. This is a necessary part of the method, as rounding errors and instability of computing algorithms can significantly alter the result.
Thus, analysis of the published works indicates that there are still open issues related to the problems of the mechanics of elastomers [12-15]. Existing approaches reduce the methods to a system of simplified hypotheses (considering a threedimensional problem as a two-dimensional; assuming linearity of deformation; incompressible or weakly compressible material etc.) and often present the method in a form not efficient for computations [16]. Significantly, these issues are related to the spatial representation of the FE. To solve this problem, the paper proposes variational relations for the tetrahedral finite element, which are developed on the base of the moment scheme of the FEM. The proposed method is based on our previous work [12]. 2
Derivation of variational relations for the tetrahedral FE
Let’s us derive the formulas of the stiffness matrix for the tetrahedral finite element
(Fig. 1).
(
N i ( x, y, z) = i 0 + i1 x + i 2 y + i 3 z
The formulas for the second and third directions (
4 v(x, y, z) = vi Ni (x, y, z) , Ni (x, y, z) = i 0 + i1x + i 2 y + i 3 z ; i=1 4 w(x, y, z) = wi Ni (x, y, z) , Ni (x, y, z) = i 0 + i1 x + i 2 y + i 3 z .
i=1
The shape functions for each face of tetrahedral FE (Fig. 1), given in the basic coordinate system, are determined by the formulas (5)
N1 ( x, y, z) = 1 − x − y − z , N 2 ( x, y, z) = x ,
N 3 ( x, y, z) = y , N 4 ( x, y, z) = z Formulas for the deformations of tetrahedral FE are defined as follows (6)
Functions that define the geometry for a similar FE in the basic coordinate system
have the form (8)
(
i i i i=1 i=1 i=1
The Jacobean transition from the basic to the local coordinate system is defined by the formula (9) where N1 x J = N1 y N1 z N2 x N2 y N2 z ... ... ...
Ni
x Ni z1i z 2i z 3i y Ni y z1i T = z11 , z12 , z13 , z14 T , z 2i T = z 21 , z 22 , z 23 , z 24 T , z 3i T = z 31 , z 32 , z 33 , z 34 T , pqr = x p yq z r p! q! r! (0 p + q + r 1). (8) (9) (10) (11) (12) of the ith node; z3i - applications of the ith node.
where i is the number of nodes; z1i - abscissae of the ith node; z 2i - ordinates Let’s define the coordinates of the nodes of the linear tetrahedral element: х1= 0, y1 = 0, z1 = 0; x2 = 1, y2 = 0, z2 = 0; x3 = 0, y3 = 1, z3 = 0; x4 = 0, y4 = 0, z4 = 1 The approximating function of displacement for a tetrahedral finite element can be represented as (11)
uk' = k' 000 + wk' 100 100 + wk' 010 010 + wk' 001 001 where wk' pqr - decomposition coefficients, pqr - a set of power coordinate functions, defined by the formula (12)
The components of the deformation tensor are decomposed into the Maclaurin series in the vicinity of the origin (13)
ij ij = eij (stg ) (stg )
stg (13)
In the decomposition of the deformation components, along with the coefficients of the deformations, there are also the coefficients of the rigid turns. This is a common reason for the slow convergence of the FE. In the proposed approach, we will dismiss these members of the series. After transformation of a given finite element, the strain tensors will have the following form (14) 11 = e10100 = k10'0b1k0'0 ; 22 = e20200 = k0'10b0k1'0 ; 33 = e30300 = k0'01b0k0'1 ; The transformation matrix A is the formula (17) (15) (16) Matrix Fijs has the following form (18) submatrices of the matrix Fijsэ are defined as (19) F1s1' = 0 0 0
The matrix of the power functions for the formula of displacements H ijkl lets present in the form (20) 2
0 0 H ijkl = J 0 0 0 0 2 0 0 0 0 2 0 0 0 0 0
The next step is to derive a formula to find specific energy of deformation for the volume change function. The matrix of connections of nodal displacements and the (21) (22) (23) (24) F1 F s' = F 2
F 3
, where 0 F 1 = 0 0 0
The matrix of power functions for the volume change formula H lets present in the form (23) K st = AT F s T H ijkl F t A+ AT F s T H F t A
ij kl
Finally, substituting the formula (24) by the formulas (17-23), we can obtain the values of the stiffness matrices for the tetrahedral finite element using the constants of ij , , and coordinates zk' .
Experimental Study: Compression of a Rubber Sprocket The following example compares two methods: the classical FEM and proposed MS-FEM, based on the moment scheme. For the modelling, FORTU-FEM system was applied [21-23].
Let’s consider the elastic cam coupling with a sprocket, intended for coaxial connection of shafts of mechanisms, e.g. a reducer and the electric motor (Fig. 2).
This coupling is made of two half couplings, on the inside is equipped with a hubcam, between which the rubber sprocket is enclosed. Sprocket teeth work on compression.
When the torque is transmitted in each direction, half of the teeth work. The efficiency of the rubber sprocket is determined by the magnitude of the buckling stress.
Sprockets for elastic cam couplings are designed to connect coaxial cylindrical shafts in torque transmission from 2.5 to 400 N/m and to reduce dynamic loads. The sprocket parameters for the computational experiment are presented in Fig. 3.
The object will be calculated with a torque direction counterclockwise.
The load force is applied to each tooth at points as far away from the center of the sprocket.
Characteristics of the material were taken from the standard ТМКЩ-С 7338-90.
The finite element discretization in this study was done by tetrahedral and parallelepiped elements (Fig. 4). Fig. 4. FE discretization of a rubber sprocket
a) tetrahedral FE; b) parallelepiped FE Fig. 5. shows a distribution of displacements in different directions.
To check the efficiency of the proposed computation scheme, it was compared with other methods in the FORTU-FEM system (Fig. 6).
A numerical experiment shows that the results obtained from the FEM using the Lagrange variational principle, with a Poisson coefficient varying from 0.470 to 0.4999, outperform the solutions obtained using the conventional FEM scheme. The application of MS-FEM for the computation of structures made of compressible materials gives the numerical results close to the analytical solutions.
Conclusion
The paper present variational relations for the tetrahedral FE based on the moment scheme of the finite element method. The results allow us to take into account the basic properties of rigid displacements for isoparametric and curvilinear finite elements, where the deformation components depend not only on derivatives of rigid rotations but also on translational and rotational displacements of the whole elements. In future work, the proposed method will be used for computation of elastomer structures in the stress and strain states. 5 4. H. Lee S, J. Shin K, S. Msolli et al. “Prediction of dynamic equivalent stiffness for rubber bushing using the finite element method and empirical modeling,” International Journal of Mechanics and Materials in Design (2019). DOI: 10.1007/s10999-017-9400-7. 5. F. Zhao, “Continuum constitutive modeling for isotropic hyperelastic materials,” Advances in Pure Mathematics (2016). DOI: 10.4236/apm.2016.69046 6. Zihan Zhao, Xihui Mu, Fengpo Du, Modeling and Verification of a New Hyperelastic
Model for Rubber-Like Materials (2019). DOI: 10.1155/2019/2832059. 7. Y. Chandra, S. Adhikari, E.I. Saavedra Flores, Ł. Figiel. Advances in finite element modelling of graphene and associated nanostructures. Materials Science and Engineering. (2020). DOI: 10.1016/j.mser.2020.100544. 8. Grebenyuk S.N. The shear modulus of a composite material with a transversely isotropic matrix and a fiber. Journal of Applied Mathematics and Mechanics. (2014).
DOI.org/10.1016/j.jappmathmech.2014.07.012. 9. FabioGalbusera, FrankNiemeyer, Mathematical and Finite Element Modeling, Biomechanics of Spine, Basic Concepts, Spinal Disorders and Treatments (2018), DOI: 10.1016/B978-0-12-812851-0.00014-8. 10. Erik R. Denlinger, Thermo-Mechanical Modeling of Large Electron Beam Builds, Thermo-Mechanical Modeling of Additive Manufacturing (2018). DOI: 10.1016/B978-0-12811820-7.00012-4. 11. Bazhenov, B.A., Sakharov, A.S., Solovei, N.A. et al. Moment scheme of the finite-element method in problems of the strength and stability of flexible shells subjected to the action of forces and thermal factors. Strength Mater (1999). https://doi.org/10.1007/BF02511170 12. V. Mezhuyev, V. Lavrik, Improved Finite Element Approach for Modeling ThreeDimensional Linear-Elastic Bodies, Indian Journal of Science and Technology (2015).
DOI:10.17485/ijst/2015/v8i30/57727/ 13. S.Kasas, T.Gmur, G.Dietler, Chapter Eleven - Finite-Element Analysis of Microbiological Structures, The World of Nano-Biomechanics (Second Edition). (2017).
DOI:10.1016/b978-0-444-63686-7.00011-0. 14. Naman Saklani, Zhenhua Wei, Alain Giorla, Benjamin Spencer, Subramaniam Rajan, Gaurav Sant, Narayanan Neithalath, Finite element simulation of restrained shrinkage cracking of cementitious materials: Considering moisture diffusion, ageing viscoelasticity, aleatory uncertainty, and the effects of soft/stiff inclusions, Finite Elements in Analysis and Design (2020), 103390. DOI :10.1016/j.finel.2020.103390. 15. G.D.Huynh, X.Zhuang, HGBui, G.Meschke, H.Nguyen-Xuan, Elasto-plastic large deformation analysis of multi-patch thin shells by isogeometric approach, Finite Elements in Analysis and Design (2020), DOI:103389. 10.1016/j.finel.2020.103389. 16. Subrato Sarkar, I.V.Singh, B.K.Mishra, A.S.Shedbale, LHPoh, A comparative study and ABAQUS implementation of conventional and localizing gradient enhanced damage models, Finite Elements in Analysis and Design (2019), DOI : 10.1016/j.finel.2019.04.001. 17. Marek Klimczak, Witold Cecot, Towards asphalt concrete modeling by the multiscale finite element method, Finite Elements in Analysis and Design (2020). 103367, DOI: 10.1016/j.finel.2019.103367. 18. Aurélien DOItrand, Eric Martin, Dominiqu, Leguillon. Numerical implementation of coupled criterion: Matched asymptotic and full finite element approach, Finite Elements in Analysis and Design (2020), 103344. DOI: 10.1016/j.finel.2019.103344 19. T.R.Walker, C.J.Bennett, T.L.Lee, A.T.Clare, A novel numerical method to predict the transient track geometry and thermomechanical effects through in situ modification of process parameters in Direct Energy Deposition, Finite Elements in Analysis and Design, (2020), 103347. DOI:10.1016/j.finel.2019.103347. 20. M.R.Javanmardi, Mahmoud R. Maheri Extended finite element method and anisotropic damage plasticity for crack propagation modeling in concrete, Finite Elements in Analysis and Design, (2019), DOI: 10.1016/j.istruc.2017.09.009. 21. Vitaliy Mezhuyev, Vladimir Lavrik, Ravi Samikannu. Development and Application of the Problem-Oriented Language FORTU for the Design of Non-Standard Mechanical Constructions. Journal of the Serbian Society for Computational Mechanics (2015) DOI:10.1007/978-3-319-07674-4_97. 22. Vitaliy Mezhuyev, Sergey Homenyuk, Vladimir Lavrik. Computation of elastomers properties using FORTU-FEM CAD system. ARPN Journal of Engineering and Applied Sciences (2015) 23. Vitaliy Mezhuyev, Vladimir Lavrik. Development and application of FORTU-FEM Computer-Aided Design System. 4th World Congress on Information and Communication Technologies (2014). DOI: 10.1109/WICT.2014.7077292