Modelling of elastoplastic, multi-scale and multi-contact problems: application to worm gears

. A nonconventional application of worm gears exploits the irreversibility of these power transmission devices in order to realize fast emergency braking. This application can be used to secure lifting devices. A limiting factor in the design of these instantaneous braking systems is the residual deformations of the worm/ wheel contacting teeth, due to the impact between them at each emergency stop. The prediction of these residual displacements requires solving of an elastic – plastic, multi-scale and multi-contact problem. Original numerical tools were developed in this study to solve the problem at global and local scales. The method has been validated by comparing the obtained results with 3D measurements on new and deformed worm/wheel pairs. In order to predict the issue of the worm gear after an impact, a criterion based on kinematic errors is proposed. Applying this criterion gives the maximal admissible torque for the braking system to be operational after the impact.


Introduction
Worm gears are one of the most important devices used for power transmission between spatial crossed axes. They are employed in applications requiring high reduction ratio, relatively low speed drive and compactness. They are indeed capable of transmitting a large reduction ratio with a single stage, consisting of a worm with helical thread and a toothed wheel, meaning small volume, low inertia and quiet operation compared to parallel axis gearing which would require multiple stages to achieve an equivalent ratio. Therefore, worm gears are widely employed in industrial applications such as automotive applications (mainly for actuators), steam turbines, lifts, conveyor belts and stringed musical instruments. However, they present some disadvantages compared to parallel axis gearing, such as delicate assembling, or high wear rate and lower efficiency due to the sliding between the contacting teeth. Another feature of worm gear is the self-locking effect that may occur when the wheel drives the system. In other terms, above some critical worm helix angle, the system becomes unable to reverse the direction of power transmission, because of the important friction induced between the worm/wheel pair.
A nonconventional application of worm gears exploits this latter feature in order to secure lifting devices during uncontrolled movements, such as over-speed, backwashing and gear failure. This feature is used in the Fastbrake© system developed by FOC Transmissions. An illustration of this system is shown in Figure 1. In the Fastbrake© system, the motion of the worm/wheel pair is synchronized in order to ensure a contactless rotation during normal operation of the lifting equipment. An electromechanical system is used to detect movement anomalies, and commands the stop of the worm servo-control if any such anomaly occurs. Only at this moment the wheel comes into contact with the worm, and because of the irreversibility of the worm gear system, the lifting device comes to stop. The response time of this braking system corresponds thus to the time needed for the wheel to take up the clearance and to come into contact with the worm. This working principle reduces considerably the response time of the braking system compared to standard disk brakes. A limitation of the system however is caused by the residual deformation of the worm and the wheel during impact following an emergency brake. Such residual deformation can actually reduce the clearance between wheel/worm teeth necessary for the contactless rotation of the worm gear during the normal operation of the lifting device. The question that arises is whether the braking system will be able to run after emergency brakes.
In order to predict the issue of the braking system, a problem with the following characteristics needs to be solved. First, it is a nonlinear elastic-plastic contact problem, therefore requiring an iterative approach adapted to contact situations [1]. Second, it is a multiscale problem since the residual deformation is a combination of global structural bending deformations and local contact deformations. Third, it is a multicontact problem since it requires the consideration of the contacts occurring simultaneously on several wheel/worm teeth [2]. In the literature, examples can be found on solving local plastic problems in gears. For instance, Huesmann et al. [3] employed an automatic incrementation technique to study the elastic-plastic contact between two teeth of a cylindrical gear. In this respect, they applied an algorithm of boundary element method, and used a fast iterative algorithm in order take into account the non-linear behavior of plastic deformation. He et al. [4] developed a model coupling an elastic-plastic finite element method together with continuous damage theory in order to investigate the contact fatigue crack initiation process of a heavy-duty wind turbine gear. They used this to model to determine the influence of plastic deformation on the fatigue life of this type of gears. However, to the best of the author's knowledge, no study has been conducted to examine the multi-scale plastic deformation of gears in general, and worm gears in particular.
In this paper, an efficient dialogue between complementary numerical tools is presented that allows overcoming the limitations of existing numerical tools solving elastic-plastic, multi-contact and multi-scale problems. Once a estimation of residual displacements is obtained, a criterion based on kinematic error is proposed in order to predict whether the braking system would be operational after an emergency stop or not.
2 Resolution of the elastic-plastic, multi-scale, multi-contact problem 2.1 General procedure Figure 2 is a diagram depicting the different steps involved in the dialogue between three numerical simulation software, namely: ROUVISLAM, ABAQUS and ISAAC, in order to predict the elastic-plastic, multi-scale and multi-contact displacement in worm gears.
Since the contact line location between the gears flanks is not known when an emergency stop occurs, the three simulation tools are used to perform the following computations at any kinematic position of the worm/ wheel pair.
-In step 1, the torque is applied on the worm/wheel pair in ROUVISLAM (developed by the LaMCoS Lab. and Mecalam Cie) in order to compute semi-analytically the quasi-static load sharing due to the applied torque. In other terms, it provides the nodal force field or equivalently the pressure field on the different worm and wheel teeth in contact, for a given kinematic position, assuming an elastic behavior. -In steps 2-3, the nodal force field is transferred to the standard ABAQUS software in order to compute by finite element method the residual teeth bending displacement field. -In steps 4-5, the pressure field is transferred to the ISAAC software (developed by the LaMCoS Lab.) in order to compute semi-analytically the residual contact displacement field. The choice of a semi-analytical solver for the resolution of the contact problem is justified by its faster convergence and its robustness compared to a finite element solver. In addition, the semi-analytical solver allows overcoming difficulties of finite element method in dealing with nodal forces. -In steps 6-8, the bending and the contact residual displacement fields resulting from the plastic deformation at the local (contact deformation) and global (tooth bending) scales are transferred to ROUVISLAM, in order to compute the kinematic errors, and to predict whether the worm gear will be able to continue to ensure the safety of the lifting system.
Several codes have been implemented in python and MATLAB in order to couple the different numerical simulation software. Namely, these codes allow for the transfer of the different fields (nodal force, pressure, bending residual displacement and contact residual displacement) between ROUVISLAM, ABAQUS and ISAAC software (steps 2, 4, 6 and 7).

ROUVISLAM: elastic, multi-contact problem
The resolution of a multi-contact problem is complex by nature, and it cannot be achieved analytically. Instead, the quasi-static computation of pressure distribution among the contacts occurring simultaneously between worm/ wheel pair is realized in the following steps: -In the 1st step, the precise geometry of the worm/wheel pair is generated through the simulation of the cutting processes [5]. -In the 2nd step, an unloaded kinematic simulation is performed in order to delimit the contact zone between the worm and the wheel [2]. A slipped mesh constituted of rectangles of constant size is then constructed following the unloaded contact line. (see Fig. 3a).
-In the 3rd step, the worm gear under load is simulated by solving simultaneously, over the slipped mesh covering the contact zone, a set of two equations, i.e. the equation of compatibility of displacements and the equation of load balance. The resolution is performed under quasistatic and linear elastic behavior assumptions.
The equation of compatibility of displacements [6] gives the angular gap after loading, y ang , at two zones (see Fig. 3b): Inside the contact zone: p i ≥ 0 and y ang Outside the contact zone: p i ¼ 0 and y ang i ≥ 0; ð2Þ with y ang where U ang 1 and U ang 2 are the angular displacements of the two surfaces in contact, ei ang is the angular gap between the 2 surfaces before loading, and a ang is the global body adjustment. These angular quantities (gap or displacements) are all given around the axis of the wheel. The angular displacements U ang i1 and U ang i2 are expressed as linear functions of the pressure: where U i1 and U i2 are the linear displacement of the two surfaces in contact, R i is the proportionality constant between the angular displacements, U ang i1 and U ang i2 , and the linear displacements, U i1 and U i2 . The coefficients C ij are the terms of the influence matrix expressing the displacements at any point i as a result of a pressure p j corresponding to a unit force f j applied at any point j. This matrix C ij is the sum of the different contributions, namely: the bending deformations of the worm and the wheel, C b;worm ij and C b;wheel ij , and the contact deformation of the surfaces in the contact zone, C c ij [7]. The bending influence matrix, C b;worm ij and C b;wheel ij , are calculated with an automatic finite element model in CATIA. In this model, the wheel's bore and the worm's shaft are clamped. Unit forces f j are then applied successively at each point j and the resulting displacements are then computed on the set of points i. In order to get only the contribution of the bending deformation and to avoid any error in proximity of the loaded nodes, the difference between computation results of the two distinct configurations is taken [8]. In the first configuration, the flanks opposite to those subjected to pressure are left free, while in the second case, these flanks are clamped.
The contact influence matrix C c ij is computed using Kalker's semi-analytical approach based on Boussinesq's equation [9,10]. This equation provides contact displacements U c i as a function of the pressure field p j : where E 1 and E 2 are the Young's modulus, and n 1 and n 2 are the Poisson's coefficients of the worm's and the wheel's materials.
On the other hand, the equation of load balance imposes that the pressure p must compensate for the driving torque C driving (defined in the frame of the worm): where s i is a rectangular surface around the point i, and n i is the normal direction to the surface s i , M i is the vector defining the position of the point i in the frame of the worm.
Combining the equation of compatibility of displacement and the equation of load balance gives rise to a recursive formula for the pressure p i at any point i: This formula is numerically solved using the fixed-point method, whose convergence is based on the global body adjustment a ang which must be constant in the contact zone. Load sharing can also be expressed in terms of nodal force field f i given by:

Force field transfer from ROUVISLAM to ABAQUS
The force field f i is calculated in the local slipped mesh constructed in ROUVISLAM, and the nodes of this mesh do not necessarily coincide with those of ABAQUS mesh. Therefore, an automatic procedure is developed to carry out the transfer of the nodal force field from ROUVISLAM slipped mesh to ABAQUS mesh (Fig. 4).
In this procedure, a nodal force f i , applied on the node i of the ROUVISLAM mesh, is assumed to have a corresponding set of nodal forces F iJ , applied on the nodes J of ABAQUS mesh. Furthermore, the nodal force f i is assumed to be a weighted average of the set of forces F iJ : where w iJ is the weight associated to a force F iJ . It can be easily shown that this assumption preserves the sum of nodal forces after field transfer: X The greater the weight w iJ is, the greater will be the contribution of the force applied on the node J in ABAQUS mesh to the force applied on the node i in ROUVISLAM mesh. More precisely, w iJ is chosen such that it is maximized when the distance r iJ between points i and J tends to zero when r iJ =∞: With this choice, we will get: with a being fixed to 0.001. The distance r 0 is manually adjusted in order to have a nodal force field F iJ on the ABAQUS mesh visually representative of the nodal force field f i on ROUVISLAM mesh.

ABAQUS: Elastoplastic global bending problem
An automatic procedure creates the ABAQUS model, including the importation of the CATIA model (generated via ROUVISLAM software), the application of the nodal forces as well as the application of boundary conditions. ABAQUS software is customized in order to carry out the calculation of the residual bending displacements. In this model, the wheel's bore and the worm's shaft are clamped. In addition, quadratic prism elements are also used to mesh the geometry. In order to account only for the residual displacements, the computations are carried in 2 steps. In the 1st step (loading step), the geometry is loaded so the elastoplastic displacements can be determined. In the 2nd step (unloading step), the geometry is unloaded so to restore the elastic displacements and get pure residual displacements.
On the other hand, in order to account only for bending displacements, each of the worm and the wheel geometries is partitioned into two parts. The 1st part comprises the contact zones and it is modeled using an elastic constitutive law, and the 2nd part comprises the rest of the geometry and it is modeled with an elasticplastic constitutive law (see Fig. 5). An automatic procedure is implemented to calculate the partitioning parameters that minimize the size of the elastic part. With such geometry partitioning, while the contact displacements will be present in the loading step, they will be completely eliminated at the end of the unloading step, so only displacements due to the bending of the foot of the teeth will be retained. Figure 6 shows an example of nodal force field and the resulting residual bending displacement of a worm.

Pressure field transfer from ROUVISLAM to ISAAC
An automatic procedure is developed to carry out the transfer of the nodal force field f i from ROUVISLAM slipped mesh to ABAQUS mesh. The field transfer is carried out by a simple interpolation: the pressure on an element k of ISAAC is obtained by interpolation of the values of the pressures on the neighboring elements i (the first 6 neighbors) of ROUVISLAM mesh.

ISAAC: Elastoplastic local contact problem
The computation of the contact residual displacements is performed with ISAAC which is far more accurate in solving contact problems than the standard ABAQUS software with equivalent computation time. ISAAC is a 3D elastic-plastic contact code based on a semi-analytical solver, initially developed by Jacq et al. [1]. It integrates two principal modules: an elastic contact solver, and a plasticity loop. The elastic contact solver module takes the external load W as input, and gives the contact surface G C , the pressure distribution p and the elastic surface displacement u e as output. This module solves simultaneously a set of three equations.
The 1st equation is the load balance between the total external load W and the pressure p applied at any point x = (x 1 , x 2 ) on the contact surface G C : The 2nd equation expresses the surface separation h as the sum of the initial surface separation h i , the rigid body displacement d, and the elastic normal displacements u e 3 of the two bodies in contact: The elastic normal displacement u e 3 can be expressed using the Boussinesq relation: In this expression, U (x, j) is the influence matrix (also called Green's function) which contains the displacement at any point x as a result of a unit load at a point j.
The 3rd equation imposes the contact condition, i.e. the non-interpenetration between the two bodies. The surface separation h must be larger than or equal to zero, or more specifically: In the contact solver, the contact surface G C is initially unknown. Therefore, the resolution algorithm is iterative so it allows the determination of both contact surface G C and residual contact displacement field. To reduce the computing time significantly, a conjugate gradient method (CGM) developed by Polonsky and Keeris is used [11]. In order to reduce further the computation time, the costly multiplication operations in each CGM iteration are performed using discrete convolution fast Fourier transform (DC-FFT) as presented by Liu et al. [12].
Then comes the plasticity loop module that takes the pressure distribution p (as determined by the contact solver) as input, and gives the residual normal displacement u r and the sub-surface residual stress s r as output. The plasticity loop requires the knowledge of the constitutive law of the plastic material, and it is performed in the framework of Betti's reciprocal theorem [13]. This latter theorem allows the calculation of total displacements u and stress s as the summation of the contributions of the elastic and the residual states: where Ais any point of the surface G of the contacting bodies, and B is any point within the volume V (or more specifically the finite plastic strain nuclei V p ) of the bodies in contact.
In order to fully solve the elastoplastic contact problem, it is necessary to account for the geometry modification as a result of the residual surface displacements u r (determined in the plasticity loop) in the elastic contact solver. Actually, an iteration between the two modules is necessary to get a converging solution of the elastoplastic problem.
ISAAC software has been used together with its various modules to solve not only plastic contact problems [1,[14][15][16][17], but also to account for the effect of the viscoelasticity [18,19], the anisotropy [20], the presence of sub-surface heterogeneities or surface coating [21][22][23] as well as the wear of the surfaces in contact [17,24]. The current work considers only the plasticity loop, whereas the contact surface G C , the pressure distribution p and the elastic surface displacement u e are computed in the contact solver of ROUVISLAM. For displacement calculation in ISAAC, the teeth surfaces were meshed using constant-size rectangular elements (element size: 2 mm Â 0.1 mm). The different teeth are then simulated one by one. Figure 7 shows a pressure field and the resulting residual contact displacement on one the worm's teeth.

Displacement field transfer to ROUVISLAM
Once the bending residual displacements (ABAQUS) and the contact residual displacements (ISAAC) are determined, ROUVISLAM can estimate the influence of these displacements on the operation of the system. This is done by integrating the residual displacements in the calculation of the kinematic error used as a criterion for the operation of the Fastbrake© system (see below).

Operation criteria based on kinematic error
In this study, we consider only the plastic deformation to determine whether the worm gear will continue to run after a first emergency stop. Another issue is related to the fracture limit of the materials used in the fabrication of the worm and the wheel. The pressures values obtained in this study are far below these limits.
Assume an assembled worm/wheel pair, with the teeth of the wheel positioned midway between the teeth of the worm at a given kinematic position. If the worm is fixed, then the clockwise and anticlockwise kinematic errors are the maximum (algebraic) displacements that the wheel can make while rotating in, respectively, the clockwise or anticlockwise direction before coming into contact with the worm. In order to ensure the operation of a meshing worm gear, the difference between the clockwise and the anticlockwise kinematic errors, also called the backlash, must be always positive. Otherwise, the teeth of the worm and the wheel would interpenetrate, meaning that the system would stop running. Things become more complicated in the Fastbrake© system after an emergency stop. First because in the Fastbrake© system, the worm and the wheel are not meshing; they actually run independently in order to ensure their contactless rotation of the worm gear. Second because after an emergency stop, some of the worm's and wheel's teeth will be deformed. So if the worm gear is to run again, then different configurations are encountered, depending on whether the worm and the wheel's teeth facing each other are deformed or not, namely: (i) non-deformed worm's zone against non-deformed wheel's zone, (ii) non-deformed worm's zone against deformed wheel's zone, (iii) deformed worm's zone against non-deformed wheel's zone, and (iv) deformed worm's zone against deformed wheel's zone. Now if the worm's teeth is initially positioned in the middle of the available space between the wheel's teeth in one of the configurations (for instance non-deformed worm's zone in front of non-deformed wheel's zone), then the positive backlash in this initial configuration avoids the penetration of the worm's and the wheel's teeth. However after a while, a new configuration (for instance nondeformed worm's zone in front of deformed wheel's zone) might be encountered, and the worm/wheel positioning made at the initial configuration will not necessarily avoid the contact between worm and wheel teeth in the new configuration. This is true even when the backlash in this new configuration is positive since the initial positioning of the worm's and wheel's teeth were made at the initial configuration, not at the new one. Therefore, in the case of the emergency brake system, another criterion based on a proper alignment of the clockwise and anti-clockwise curves is proposed, as shown in Figure 8, which guarantees the contactless operation of the worm gear. Indeed, if the clockwise kinematic error of the new configuration is larger than the anti-clockwise kinematic error of the initial configuration, then contactless rotation is a guarantee of the non-meshing worm gear at the new configuration. The initial and the new configurations are however arbitrary. Therefore, the contactless operation of the non-meshing worm gear can be guaranteed in all configurations, only with a positive "gap" between the smallest clockwise kinematic error and the largest anti-clockwise kinematic error, computed for the different configurations. Otherwise, if instead of the positive gap, the clockwise kinematic error at some configuration becomes smaller than the kinematic error in another configuration, then the worm/wheel pair may interpenetrate.
It is possible to know in advance the configurations corresponding to the minimum clockwise kinematic error and the maximum anti-clockwise kinematic error. Indeed, as shown in Figure 9, the deformation of the wheel involves a shift in the kinematic error, whereas the deformation of the  worm involves only a peak. In fact, when scanning the kinematic positions, the deformed zone always remains in contact for the wheel, while it only makes a brief passage for the worm. With the two deformations involved, we observe a shifted peak in the clockwise kinematic error meaning that it is the less favorable configuration as one could expect. On the other hand, the undeformed worm gear configuration involves no change in the anti-clockwise kinematic error change, meaning that it is the most favorable configuration, as one could expect as well. Therefore, in order to verify the contactless operation after an emergency brake, the less favorable configuration, i.e. deformed worm-wheel pair, must be compared with the most favorable configuration, i.e. undeformed worm-wheel pair.

Validation
Numerical simulations are compared to measurements in order to validate the method presented above to solve the elastic-plastic, multi-scale and multi-contact problem encountered in the worm gear used in the Fastbrake© system. For the numerical simulations, the geometrical parameters were precisely extracted from 3D point cloud as obtained with a Coordinate Measuring Machine (CMM) (Zeiss MC 550, Germany). A portion of the wheel containing 1/5 of the teeth was used for the different simulations in order to reduce computation time. In addition, the constitutive elastic-plastic laws of the material of the wheel and the one of the worm were identified via a tensile test machine (ZwickRoell, France). The measured data were used directly in ABAQUS model or fitted with a swift law for ISAAC simulations.
The contact tests data were provided by FOC Transmissions. These tests were conducted with the help of a pressing device (see Fig. 10a). In this device, the worm gears were subjected to impacts of various forces, by pressing the worm against the wheel, with the latter being clamped inside a housing. Figure 10b shows an example of the impacted wheel. The geometrical parameters were measured again for the impacted worm gears using the CMM and compared to those obtained previously on the non-impacted worm gears, in order to determine the residual displacements.
The kinematic position of the worm/wheel pair in which the impact occurred was identified by comparing the position and the shape of the impact-induced marks with those of the resulting pressure zones as obtained with ROUVISLAM software (see Fig. 10b). After the comparison, the kinematic position of impact of the examined worm/wheel pair turned out to be equal to 70 mrad. Figure 11 shows the comparison between numerical and experimental values of the tooth space of the worm and the wheel subjected to a load of 392400N (corresponding to a torque on the worm of 750 Nm) at the kinematic position = 70 mrad. The close proximity between simulation and measurements at such low scale deformations demonstrates the robustness of the developed numerical model.

Parametric study: maximum displacement and limit of operation
Using the different tools developed in the study, it is possible to study the influence of the different geometrical and the operational parameters on the value of the residual displacements of the worm/wheel pair after an emergency stop. This study can be particularly used to predict the maximum admissible load (or torque) to which the worm gear can be subjected so the Fastbrake© system can run again after the emergency brake.
In what follows, the influence of the load and the kinematic position is considered. Figure 12 shows the maximal total displacement obtained at different loads and kinematic positions. Several observations can be made. First, larger displacements are observed in the case of the worm. This may be explained by the fact that the elliptical form of the worm's thread rigidifies this latter against bending. Second, a net increase of displacements is observed with increasing load as one may expect, and the increase depends on the kinematic position. Figure 13 shows the maximal clockwise kinematic error obtained at different loads as a function of the kinematic position ranging from À70 to +70 mrad. The clockwise kinematic errors are compared in the same figure with minimal anti-clockwise error (obtained with non-deformed worm and wheel's teeth). This curve suggests a maximal admissible torque value somewhere between 900 and 1050 Nm. Above this maximal admissible value, the system is not guaranteed to continue to run.    In this work, an efficient dialogue procedure is developed in order to solve a complex contact problem, i.e. elasticplastic, multi-scale and multi-contact problem. Finite element method is used to compute the influence coefficients of the worm and the wheel (ROUVISLAM), as well as the bending of these with the applied nodal force field (ABAQUS). Semi-analytical method is used to compute the load sharing due to the applied torque (ROUVISLAM), and to find out the contact residual displacements due to the applied pressure field (ISAAC). The method gives good agreement with measurements. A criteria based on the kinematic errors is proposed to predict whether the system will continue to run after a first emergency stop or not. A parametric study was carried out in order to find the maximum load to which the system can be subjected so it can continue to run after the emergency stop. The procedure can be further enhanced for large displacements, by including an iteration process between load sharing calculation (ROUVISLAM software) and residual displacements calculation (ABAQUS and ISAAC  softwares). Indeed, the computed load sharing is made under elastic linear assumption, which can be considered valid for small displacements. By including the effect of residual displacements, in particular the flattening of the surfaces under the contact, on the computation of the load sharing, the elastic-plastic effect can be modeled more precisely. The developed procedure can be also used as a tool for the design of worm gears. For this design study, the impact on the operation criterion of each of the worm/wheel geometrical and assembly parameters must be determined. Once done, the developed numerical procedure can be used to optimize these parameters in order to guarantee longer-life to worm gear systems for the fast braking application. Finally, the integration of contact aspects in the developed numerical tools (namely with ISAAC) opens the possibility to study the effect of surface treatment on the resistance of the worm gear system to fracture. Poisson's coefficient of the worm n 2 Poisson's coefficient of the wheel G Body surface G C Contact surface s e Elastic stress s r Residual stress V Body volume V p Plastic strain nuclei