Study of flexible couplings non-linear dynamics using bond graphs

The aim of this work is to model the dynamics of flexible couplings. On the basis of a non-linear mathematical model solved by bond graph, the ranges of excitation frequency were determined, in which the movement of the couplings is chaotic. For three couplings, the 3D distributions of the largest Lyapunov exponent and correlation dimension diagram (CDD) were plotted. The proposed diagram (CDD) illustrates how the geometric structure of the attractor changes when the conditions of excitation change. The classic Poincare cross-section, completed by us with the density of points distribution, significantly enhances information about geometrical structures of strange attractors. It has been shown that in relation to large ranges of changes in the control parameter, the geometric structure of the strange attractor is stretched and curved. The areas with the highest densification of the Poincaré cross section are most often located in places where the chaotic attractor is curved.


Introduction
The subject of the model research included in this paper is the tire flexible coupling. The basic task for power transmission components, including mechanical couplings, is to connect two shafts to transfer the torque. From a theoretical point of view, the optimum mechanical properties of the flexible coupling to the dedicated drive system should not only be based on the transmitted torque (as in most cases) but also mainly, based on the viscoelastic properties of the flexible connector present in the flexible coupling. These properties are identified based on static or dynamic characteristics recorded under experimental conditions [1][2][3][4][5]. Limited access to this type of research results in the fact that the designers during the design of the couplings do not take these parameters into account when fitting the coupling to the drive system. The viscoelastic properties of the tire couplings are mainly determined by the characteristics of a rubber connector that joins the active and passive coupling member [6,7]. The basic task of the connector is to suppress the excited mechanical vibration and to minimize, to some extent, the dynamic surplus that occur in the transient states of the drive system. The beneficial effect of a flexible coupling to the operation of the drive system is significantly dependent on the correct choice of torsional stiffness [8,9]. Its incorrect selection may result in increased dynamic forces or in extreme situations lead to the appearance of resonance. Depending on the coupling connector material used, couplings with linear and non-linear mechanical characteristics are distinguished [8,10,11].
Knowledge of mechanical characteristics plays a significant role in modelling [12,13]. While the results of analyzes based on linear models are predictable and easy to interpret they are relevant for small amplitude oscillations. In contrary, nonlinear models for moderate excitations can hardly be said for predictability because of the irregular phenomena occurring in the system [14,15]. It is considered that chaotic phenomena are a typical feature of nonlinear systems, and systems of this type are highly sensitive to initial conditions.

Formulation of a mathematical model
In this paper, numerical and laboratory tests of flexible couplings were performed. The discussed system was modelled as a two-degree of freedom system in which rotating discs J1 and J2 were coupled via a non-linear elastic element cC and a linear viscous damper bC. Whereas the spring properties of the connector are described by the third-degree polynomial function, which maps the progressive mechanical characteristics of the coupling. The movement of the system is caused by an external harmonic torque MN, which was applied to the active member. On the other hand, the load torque on the passive member is constant MO (Fig. 1).
To derive differential equations of motion, a modified method of bond graphs was used. In the original form proposed by Paynter [16] and developed by Karnopp and Rosenberg [17], mathematical models derived from bond graphs are obtained in the form of state equations. The modification we apply allows direct derivation of differential equations of motion. This is achieved by introducing additional edges representing fictitious sources of flow type variables: Sf W P ' 1 ; Sf W P ' 2 : Their task is to assign to incidental 1-junction, flow type variables, representing the displacement speeds of the coupling discs. This solution allows the use of differential causality of kinetic energy storage elements: I W J 1 ; I W J 2 ; as a result of which there are no  conflicts of causality in the 1-jounction nodes. Topological structure of the so-called Lagrang's bond graph, describing the dynamics of the flexible coupling, is shown in the figure (Fig. 2). It maps cause-and-effect relations occurring between particular elements of the graph.
The relations occurring in the 1-junction and 0-junction of the bond graph constitute a formal basis for deriving the differential equations of motion. From the point of view of computer simulations implemented via bond graphs, an explicit representation of motion equations is not required. Nevertheless, bearing in mind quantitative and qualitative research, focusing in particular on the description of chaotic phenomena, knowledge of the mathematical model becomes indispensable.
In the analysed graph, there are three 1-junction nodes, with only two representing the generalized speeds of the mathematical model. These are incidental nodes with elements that store kinetic energy I. An incident 1-junction with an element that dissipates energy R and stores the potential energy C, has the task of summing disruptive and elastic force. The equations of motion are derived based on cause and effect relations occurring in the nodes of the 1-junction graph: With regard to the 0-junction, the cause-and-effect relations occurring in it are depicted in dependencies: Additionally, the relations occurring in the 1-junction with elements R and C should be taken into account: Equations from (1) to (3) form the formal basis for deriving the system of motion equations. Nevertheless, it is necessary to define the cause-effect relationships that characterize the individual edges of the graph: Including in (1) relations (2) to (4), finally a system of differential equations of motion is obtained: After entering a new variable u, defined as the difference of angular displacement of the clutch discs, the system of Eq. 5 reduces to one differential equation. In addition, we assume that the change of the driving torque, forcing the coupling movement, is described by the harmonic function MN = cos(ωW1t). However, the value of the torque MO loading the coupling remains constant over time. Taking into account the above assumptions and the nonlinear characteristics of the flexible element of the coupling leads finally to the differential equation: where:

Laboratory tests
Laboratory tests were carried out in relation to three couplings, differing in the size of the flexible connector. The laboratory stand for dynamic testing of flexible couplings consists of the cooperation of two M1 and M2 induction motors connected with each other by means of the tested flexible coupling, assuming that the M1 motor maintains a constant rotational speed (Fig. 3). Selected parameters of the tested couplings are shown in Table 1. Three sizes of a typical coupling with flexible connectors used in the drive systems of heavy working machines were selected.
The mechanical characteristics were identified at the laboratory stand and then approximated by the third-degree polynomial function. The obtained characteristics are shown in Fig. 4.
It was shown that the characteristics of the tested couplings are non-linear and depend on the size of the connector used. The largest values of force dependence on the torsion angle occur for the largest coupling.

Results of model simulations
Model tests were performed for a clutch with increased susceptibility of flexible connectors. Couplings of this type are designed to carry heavy loads with relatively small dimensions and compact design. The movement of the tested mechanical system is caused by a quadripolar induction motor with rated power P = 7.5 kW and speed ns = 1450 min -1 . During computer simulations, perfectly rigid shafts connecting the clutch with the electric motor and other components of the drive system were adopted. The adoption of such model assumptions was aimed at eliminating the influence of dynamic effects of the propulsion system, which are caused by the elastic deformations of the drive system components. As a result of such procedure, the influence of parameters characterizing the mechanical properties of the clutch on its dynamics can be assessed. In addition, during the numerical experiments, it was assumed that the value of the torque applied to the coupling is 80% of the nominal torque value. Substitute viscous damping coefficient δ, captures total energy losses in a flexible connector and resistance to motion in bearings.
Research on solutions of non-linear mathematical models is accomplished by using many numerical tools. Their usefulness is unprecedented because it allows to present and interpret the obtained results in an intuitive way. In the conducted research, only three of many were used. The chaotic motion zones were identified through the Lyapunov exponent. Poincaré maps [18] were used to characterize and visualize strange attractors, while the correlation dimension was used to characterize the observed geometric structures.
In the first stage of the model tests, areas, where the coupling movement is chaotic, were identified. In order to implement such a formulated research task, the maximum Lyapunov exponents were calculated, which were depicted in the form of three-dimensional graphs (Fig. 5). Its basic purpose is to distinguish between unpredictable chaotic behaviour and predictable periodic behaviour. The estimation of the character of the solution requires the calculation of at least the maximum Lyapunov exponent [19,20] or n of Lyapunov exponents. At the same time, the obtained multicoloured maps depend to a large extent on the assumed a priori distance ε, measured between two zero initial conditions of the trajectory. In the conducted numerical experiments, this distance was assumed ε = 10 -5 . To estimate the largest Lyapunov exponent, the numerical algorithm proposed by Wolf was used [21]. Chaotic movement takes place when the maximum Lyapunov exponent λ > 0, takes positive values. Bearing in mind the exact distinction of the chaotic motion areas, a zero plane was introduced, which in the graphs represents a transparent surface in a bright green colour. Areas (bands) where the Lyapunov exponent takes positive values are limited by curves, which can be approximated with a polynomial function of the third degree with high accuracy. At the same time, these curves are determined by the static characteristic reflecting the mechanical properties of the clutch. Comparison of three-dimensional maps of the largest Lyapunov exponent (Fig. 5) indicates that the increase in the static slope of the clutch reduces the area in which the movement is irregular. Note that the topology of 3D corrugations is similar for the considered coupling types with maximizing of the magnitudes for the D316 case. Subsequently, correlation di- Fig. 7 The distribution density of the Poincaré cross-section points obtained with zero initial conditions u .0/ = 0, P u .0/ = 0 for different coupling types (D375, D360, and D316) mension diagrams (CDD) were generated, based on which the attractor's complexity can be assessed (Fig. 6).
With respect to the widest bands presented in the CCD diagrams, the values of the frequency ω have been determined, for which the chaotic attractor assumes approximate to the maximum values of the correlation dimension. On its basis, at the given value p1 = 1, the width of the zones in which the couplings movement is irregular is determined. The smallest D316 coupling has three basic bands: first {0.77 < ω < 0.89}, second {0.99 < ω < 1.18} and third {1.65 < ω < 1.69}. For the second of tested couplings D360 there are two bands, the first one {0.84 < ω < 1} and second {1.39 < ω < 1.41}, which are located in a cross-section p1 = 1. For the third coupling D375, there are also two bands in which the movement occurs irregularly, first {0.77 < ω < 0.92} and sec-ond{1.29 < ω < 1.31}. However, in both the D316 and D375 couplings, there are narrow zones in the first bands ω 0.02, where the dynamics of the clutch is regular.
On the basis of plotted images of Poincaré cross-sections [22,23], which were additionally supplemented by analysing the density of point distribution on the phase plane, they were found to exhibit geometrical similarity. This is particularly noticeable with respect to the attractors of D360 and D375 couplings (Fig. 7).
The observed similarity should be related to the fact that the frequency values ω are located in the main variation bands of the largest Lyapunov exponent (Fig. 5). The reference of chaotic attractors to the characteristics (Fig. 4) indicates that its inclination affects proportionally to the displacement range u. In addition, the increase in the slope of the mechanical characteristics of the coupling causes the reduction of the attractor's correlation dimension. In individual bands, where the coupling work takes the form of a chaotic movement, various geometrical structures of chaotic attractors are observed.

Conclusions
Model tests were carried out for various values of the mathematical model parameters in order to demonstrate the complexity of the non-linear dynamics of the coupling. On the basis of the conducted research, the following conclusions can be specified: 1. Plotting the CDD diagram allows setting the value of the control parameter for which the chaotic attractor takes the most complex geometrical structure. 2. From the point of view of stable and regular coupling operation, during its design, it is necessary to aim for the flexible element of the coupling to be characterized by a high stiffness, which is confirmed by the three-dimensional maps of the largest Lyapunov exponent (Fig. 4). 3. In order to ensure predictable and regular work of flexible coupling, it is recommended to avoid areas in which the maximum Lyapunov exponent takes positive value. This can be achieved by modifying the mechanical properties of flexible couplings connectors. 4. It was shown that the characteristics of the tested couplings are non-linear and depend on the size of the flexible connector used.
The numerical results of the nonlinear model of the flexible tire coupling indicate the need for further model studies to assess the damping effect of the coupling connector on its dynamics. In addition, qualitative analyzes of the individual parameters influence the structure of the 3D Lyapunov exponent and the Poincare cross section should be performed. In this paper, we used the simplest choice of the nodal initial conditions. Therefore we have not investigated multiple solutions which usually appear in such nonlinear systems [24].
Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http:// creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.