**Mathematical Modelling of Drive System with an Elastic Coupling Based on Formal Analogy between the Transmission Shaft and the Electric Transmission Line**

#### **Andrzej Popenda 1,\*, Marek Lis 1, Marcjan Nowak <sup>1</sup> and Krzysztof Blecharz <sup>2</sup>**


Received: 14 January 2020; Accepted: 28 February 2020; Published: 4 March 2020

**Abstract:** In the paper, the kinematic structure of the transmission shaft between the driving motor and the working mechanism is studied. The analysis is based on electrical and mechanical similarities. The equivalent circuits, typical for electrical systems, are defined for the transmission shaft concerned. Modelling of the transmission shaft based on a formal analogy between the transmission shaft and the electric transmission line is also proposed. The results of a computer simulation and experimental test are presented. The results confirm the high conformity of the proposed mathematical model with the physical object.

**Keywords:** electric drive; transmission shaft; electric transmission line; electrical and mechanical similarities; kinematic structure; equivalent circuit; mathematical modelling

#### **1. Introduction**

Rapid development and diversity of electric drive systems require extensive theoretical and practical knowledge as well as the use of a wide range of theoretical and constructional solutions from the contemporary designers and users. Electric drive systems are used in all branches of industry; therefore, the trouble-free operation of these systems is of crucial importance. An analysis of the states of electric drive operation is often related to ensuring its safety. A key issue is to detect the mechanical resonance phenomenon or the phenomena close to resonance as the most unsafe ones for the system. Ignoring the analysis and the ongoing diagnostics of vibrations in electromechanical systems results in many failures. The example, including fatalities, is inter alia the failure of turbo generator set in the power plant Kostromska (Russia). In the commission view, the cause of the failure was ignoring torsional vibrations, which led to the resonance in the turbine and transmission shaft included in the turbo-generator set.

A very interesting case of the electric submersible pumps system with a multi-component system powered by a soft start system is shown in [1]. In this real case scenario, the incorrect control case led to resonance phenomena and a shaft breakdown. Therefore, comprehensive research is needed to analyse the operating states of electrical and electromechanical systems. One method to investigate the dynamic behaviour of electrical and electromechanical systems is an experimental test directly on the system. In many cases, the experimental tests on the system are not possible, due to technical reasons, costs of experiments or their duration. The alternative approach to the investigation into the system behaviour is the formulation of the model and carrying out the investigations into the model. The model of a system is the representation using mathematical relationships, a physical model or an analogue

model. The physical model of the real-world system includes smaller-scale components having the same physical nature as a real-world system. In the physical model, there are phenomena described by the same dependencies, differing only in the order of magnitude. The analogue model includes components having a different physical nature than the real-world system, but these components are easier to implement. The analogue models are based on physical analogy. An electrical circuit may be defined as the analogue model of the mechanical system. The mathematical model is the set of mathematical relationships, on the basis of which the behaviour of a system can be predicted. The mathematical model is most often formulated as differential equations describing the operation of a system. Each physical model has a corresponding mathematical model.

The mathematical model of a mechanical object is usually a set of partial differential equations. They are hard to solve both analytically and numerically. Discrete models of systems include ordinary differential equations; therefore, they are used in practice most often. Real-world mechanical systems are usually nonlinear, where the nonlinearity is determined by material properties, clearances, nonlinear nature of dissipation forces and characteristics of elastic elements. A limited possibility of analysis of nonlinear differential equations leads towards the use of linear models or the linearization of those. There are a lot of physical systems, which may be represented by linear models with accuracy acceptable for practice. In the process of creating the models of physical systems, the equations of motion are formulated.

An electric motor, being part of an electric drive, is coupled with a working mechanism via a driving shaft that is an element of mechanical power transmission. Mechanical power transmissions can be a single-path or a multi-path and can also include gear trains and clutches [2]. The long driving shafts, defined as transmission shafts, are used first of all in the drive systems for the steel industry, mainly in the drive systems for rolling mills—the transmission shafts are over 10 m long and their diameters are of 0.5 to 0.8 m. Transmission shafts are also used in drive systems for polymerization reactors. The length of these shafts is from 4 to 7 m. Moreover, transmission shafts are used in hydro generator sets, ship drive systems, submarine drive systems, etc. Depending on the length and cross-section, transmission shafts can demonstrate different susceptibilities to the impact of the moment of torsion, as measured by a value of the angle of twist [2]. In the case of short mechanical coupling, the value of the angle of twist is insignificant and may be omitted by the assumption of rigid mechanical coupling, whereas in the case of longer mechanical coupling the value of the angle of twist cannot be ignored and such coupling should be considered as the elastic one [1–10].

In order to model the transmission shaft, classical methods based on the d'Alembert principle or Newton's second law as well as the variational method are used. The use of the variational method results in the distributed parameter model consisting of partial differential equations [3]. An example of a partial differential equation used to describe a long elastic element is the wave equation [11]. Partial differential equations may be solved analytically, but this is both tedious and time-consuming. Thus, these equations are usually transformed into ordinary differential equations by using finite differences [2]. The model of the long elastic element, based on the ordinary differential equations, corresponds with a multi-mass lumped parameter model, e.g., in [2,12,13], especially the two-mass model, e.g., in [6,8,14], which does not guarantee the accurate results of numerical analysis. The equations included in the multi-mass model are usually solved using numerical integration.

Research on stresses and oscillations of long drive shafts using various methods of analysis is undertaken in the literature [15–18]. In [18], the author presents a mechanical analysis of the propeller shaft on a warship in the event of a variable moment of shaft load. Modelling of the object is selected on a classic approach using a system of nonlinear equations. A similar approach, but using a different model and a different scale, with a three-mass system, is shown in [16]. The author shows comparative results between the analytical calculations of the two-mass system and the numerical calculations for the three-mass system. An extended approach to analysing the issue of three-mass and two-mass system including changes in the moment of inertia is shown in [18]. In all of the above described publications, the authors used a classical approach in different variations.

In this paper, the kinematic structure of the transmission shaft between the driving motor and the working mechanism was studied. The telegrapher's equations and their d'Alembert solution are proposed for the mathematical description of the transmission shaft. The proposed modelling is an alternative to the multi-mass model based on ordinary differential equations as well as a distributed parameter model based on partial differential equations. The advantage of the proposed model is its simplicity, due to the fact that the model is based on discrete equations, which do not require numerical integration in contrast to the models based on differential equations. The presented modelling is based on electrical and mechanical similarities. Identifying these similarities is very helpful for electricians in finding a relevant interpretation of mechanical systems, which is particularly important in the case of professionals dealing with electromechanical energy converters or drive systems [2]. The aim of the presented study was to achieve very high conformity of the mathematical model with the physical object, whilst keeping a short time for numerical calculations.

#### **2. Methods of Modelling of Transmission Shaft on the basis of Electrical and Mechanical Similarities**

The mass of the real transmission shaft (or any element) is distributed continuously. Representing such a transmission shaft by a model based on lumped parameters causes discrepancies in the results of the analysis in relation to accurate models. These discrepancies decrease with the number of points of concentration in the model. Dividing a transmission shaft with distributed mass into several shorter elements, described by lumped parameters such as mass, elasticity and damping of the *i*-th element, it is possible to obtain the results of computer simulation essentially not different from the results obtained if the shaft is divided into infinite number of elements, that corresponds with the wave model. The process of the abovementioned division is referred to as discretization of kinematic structure and corresponds with the multi-mass lumped parameter model of the transmission shaft.

In order to present the synthesis of the proposed models of the transmission shaft, the following assumptions have been adopted:


On the basis of the above-given assumptions, the kinematic structure of transmission shaft (Figure 1) transferred into *m* discrete elements was formulated [2], where *J*1, ... , *Jm*, *Cs*,12, ... , *Cs*,*m*−1,*m*, *Sc*,12, ... , *Sc*,*m*−1,*m*, *D*12, ... , *Dm*−1,*<sup>m</sup>* are the moments of inertia, coefficients of torsional stiffness, coefficients of torsional susceptibility (*Sc* = 1/*Cs*) and coefficients of viscous friction inside the respective elements of divided transmission shaft; *D*1, *Dm* are the coefficients of friction defined for bearings.

$$\underbrace{J\_{\begin{subarray}{c}\mathcal{I}\in\mathcal{I}\\O\_{1}\end{subarray}}}\_{D\_{1}}\underbrace{\begin{subarray}{c}J\_{1}\\\hline\\D\_{1}\end{subarray}}\_{D\_{12}}\underbrace{C\_{s,12}}\_{\begin{subarray}{c}\begin{subarray}{c}\mathcal{I}\_{s,12}\\D\_{12}\end{subarray}}\end{B\_{s,23}}\_{D\_{23}}\underbrace{\dots\dots\dots\dots\dots\right}\_{=\prod}\underbrace{\begin{subarray}{c}C\_{s,m-1,m}\\D\_{m-1,m}\end{subarray}}\underbrace{\begin{subarray}{c}J\_{m}\\\hline\\D\_{m}\end{subarray}}\_{D\_{m}}\underbrace{\begin{subarray}{c}\mathcal{I}\_{m}\\\hline\\D\_{m}\end{subarray}}\_{D\_{m}}$$

**Figure 1.** Kinematic structure of transmission shaft transferred into *m* discrete elements.

In accordance with the electromechanical energy conversion theory, the generalised coordinates are forces or torques in mechanical systems and voltages in electrical circuits, whereas the generalised velocities are the velocities or angular velocities in mechanical systems and currents in electrical systems. As a consequence, the parameters of transmission shaft are equivalent to those of the transmission line, i.e., the moment of inertia is equivalent to inductance, the torsional susceptibility is equivalent to capacitance, the coefficient of friction, if any, on the shaft surface is equivalent to the resistance and the coefficient of the viscous friction inside the shaft is equivalent to the inverse of the conductance. Adopting the abovementioned analogies and considering the specificity of transmission shaft, the equivalent circuits for the transmission shaft (Figure 2) has been created.

**Figure 2.** Equivalent circuit corresponding with the kinematic structure depicted in Figure 1.

On the basis of the above given equivalent circuit the *m* equations of moments and torques as well as *m*−1 equations of angular velocities, corresponding with the Kirchhoff's circuit laws typically applied in order to analyse electrical circuits, may be written for the mechanical coupling concerned [2]:

$$M\_1 - D\_1 \omega\_1 - \frac{d}{dt}(I\_1 \omega\_1) - D\_{12}(\omega\_1 - \omega\_2) - M\_{\varepsilon, 12} = 0,\tag{1}$$

$$-M\_{\varepsilon,k-1,k} + D\_{k-1,k}(\omega\_{k-1} - \omega\_k) - D\_{k,k+1}(\omega\_k - \omega\_{k+1}) - \frac{d}{dt}(l\_k \omega\_k) - M\_{\varepsilon,k,k+1} = 0 \\ \text{for} \quad k = 2, \dots, m-1,\tag{2}$$

$$D\_{\mathbf{c},m-1,m} + D\_{m-1,m}(\omega\_{m-1} - \omega\_m) - \frac{d}{dt}(f\_m \omega\_m) - D\_m \omega\_m - M\_m = 0,\tag{3}$$

$$
\omega\_k - \omega\_{k+1} = \frac{d}{dt} \Big( S\_{c,k,k+1} M\_{c,k,k+1} \Big) \qquad \text{for} \qquad k = 1, \ldots, m-1. \tag{4}
$$

On the basis of the equivalent circuit (Figure 2) it can be concluded that the transmission shaft may be analysed in a similar way as the electric transmission line [2]. In order to formulate the distributed parameter model of transmission shaft, the following equivalent circuit of a shaft section, represented by the two-port system, can be used (Figure 3). The adopted length of the shaft section is Δ*x*. It is assumed that the section Δ*x* is short enough to use, in respect of this section, the dependencies adequate for the mathematical description of the lumped parameters model.

**Figure 3.** Equivalent circuit of a transmission shaft section.

*Energies* **2020**, *13*, 1181

The equivalent circuit of transmission shaft section (Figure 3) is similar to that of the transmission line section [19]. There is a formal analogy between the lossless transmission shaft and the lossless transmission line.

Balancing moments of torsion around the closed-loop and angular velocities at the point leads to the following dependencies [2]:

$$M(\mathbf{x}, t) = f' \Delta \mathbf{x} \frac{\partial \omega(\mathbf{x}, t)}{\partial t} + M(\mathbf{x} + \Delta \mathbf{x}, t), \tag{5}$$

$$
\omega(\mathbf{x},t) = S\_c' \Delta \mathbf{x} \frac{\partial \mathcal{M}\_c(\mathbf{x} + \Delta \mathbf{x}, t)}{\partial t} + \omega(\mathbf{x} + \Delta \mathbf{x}, t), \tag{6}
$$

$$M(\mathbf{x} + \Delta \mathbf{x}, t) = D' S\_c' \frac{\partial M\_c(\mathbf{x} + \Delta \mathbf{x}, t)}{\partial t} + M\_c(\mathbf{x} + \Delta \mathbf{x}, t), \tag{7}$$

where *J* , *D* , *S <sup>c</sup>* are specific (per unit of length) parameters of the shaft, i.e., specific moment of inertia, specific coefficient of friction and specific coefficient of torsional susceptibility, respectively, *S <sup>c</sup>* = ρ/*GJ* , ρ is the mass density in kg/m3 and *G* is the shear module in GPa. Transforming Equations (5), (6) and (7) as well as assuming Δ*x* → 0 the following equations can be obtained:

$$-\frac{\partial M(\mathbf{x},t)}{\partial \mathbf{x}} = f' \frac{\partial \omega(\mathbf{x},t)}{\partial t},\tag{8}$$

$$-\frac{\partial \boldsymbol{\omega}(\mathbf{x},t)}{\partial \mathbf{x}} = S\_c' \frac{\partial \mathcal{M}\_c(\mathbf{x},t)}{\partial t},\tag{9}$$

$$M\_{\mathbf{c}}(\mathbf{x},t) = M(\mathbf{x},t) + D' \frac{\partial \omega(\mathbf{x},t)}{\partial \mathbf{x}}.\tag{10}$$

The further equations result from Equation (8) as well as Equations (9), (10) and (11):

$$-\frac{1}{l'}\frac{\partial^2 M(\mathbf{x},t)}{\partial \mathbf{x}^2} = \frac{\partial^2 \omega(\mathbf{x},t)}{\partial \mathbf{x} \partial t},\tag{11}$$

$$-\frac{\partial \boldsymbol{w}(\mathbf{x},t)}{\partial \mathbf{x}} = S\_c' \frac{\partial \mathcal{M}(\mathbf{x},t)}{\partial t} - \frac{S\_c' D'}{f'} \frac{\partial^2 \mathcal{M}(\mathbf{x},t)}{\partial \mathbf{x}^2}.\tag{12}$$

A lossless transmission shaft may also be considered, if *D* = 0:

$$-\frac{\partial M(\mathbf{x},t)}{\partial \mathbf{x}} = f' \frac{\partial \omega(\mathbf{x},t)}{\partial t}, \quad -\frac{\partial \omega(\mathbf{x},t)}{\partial \mathbf{x}} = S'\_{\varepsilon} \frac{\partial M(\mathbf{x},t)}{\partial t}.\tag{13}$$

If the shaft parameters do not change along the axis, Equation (13) may be decoupled in respect of moment of torsion and angular velocity:

$$-\frac{\partial^2 M}{\partial \mathbf{x}^2} = f' S\_c' \frac{\partial^2 M}{\partial t^2}, \quad -\frac{\partial^2 \omega}{\partial \mathbf{x}^2} = f' S\_c' \frac{\partial^2 \omega}{\partial t^2},\tag{14}$$

$$-\frac{\partial^2 M}{\partial \mathbf{x}^2} = \frac{1}{v^2} \frac{\partial^2 M}{\partial t^2}, \quad -\frac{\partial^2 \boldsymbol{w}}{\partial \mathbf{x}^2} = \frac{1}{v^2} \frac{\partial^2 \boldsymbol{w}}{\partial t^2}, \boldsymbol{v} = \frac{1}{\sqrt{\|'\boldsymbol{S}'\_c\|}} = \sqrt{\frac{\overline{\mathbf{G}}}{\rho}},\tag{15}$$

where *v* is phase velocity in m/s. The general solution of the equations of type (15) was given by d'Alembert [19]:

$$M(\mathbf{x}, t) = \frac{1}{2} \left[ M\_1 \left( t + \frac{\mathbf{x}}{v} \right) + M\_1 \left( t - \frac{\mathbf{x}}{v} \right) \right] + \frac{v}{2} \int\_{t - \frac{\mathbf{x}}{v}}^{t + \frac{\mathbf{x}}{v}} \frac{\partial M(0, \tau)}{\partial \mathbf{x}} d\tau,\tag{16}$$

*Energies* **2020**, *13*, 1181

with the boundary conditions:

$$M\_1(t) = M(\mathbf{x}, t) \Big|\_{\mathbf{x} = 0^\prime} \quad \frac{\partial M(0, t)}{\partial \mathbf{x}} = \left. \frac{\partial M(\mathbf{x}, t)}{\partial \mathbf{x}} \right|\_{\mathbf{x} = 0}. \tag{17}$$

On the basis of the first Equation (13) and the Equation (16) the following Equations (18), (19) and (20) may be written:

$$\frac{\partial M(0,t)}{\partial \mathbf{x}} = -f' \frac{\partial \omega(0,t)}{\partial \mathbf{t}} = -f' \frac{\partial \omega\_1(t)}{\partial \mathbf{t}},\tag{18}$$

$$\frac{\upsilon}{2} \int\_{t-\frac{\pi}{\upsilon}}^{t+\frac{\pi}{\upsilon}} \frac{\partial M(0,\tau)}{\partial \mathbf{x}} d\tau = -\frac{1}{2} \sqrt{\frac{l'}{S'\_c}} \int\_{t-\frac{\pi}{\upsilon}}^{t+\frac{\pi}{\upsilon}} d\omega\_1(\tau) = -\frac{\upsilon\_\upsilon}{2} \Big[\omega\_1(t+\frac{\pi}{\upsilon}) - \omega\_1(t-\frac{\pi}{\upsilon})\Big],\tag{19}$$

$$M(\mathbf{x}, t) = \frac{1}{2} \left[ M\_1 \left( t + \frac{\mathbf{x}}{\upsilon} \right) + M\_1 \left( t - \frac{\mathbf{x}}{\upsilon} \right) \right] - \frac{z\_{\upsilon}}{2} \left[ \omega\_1 \left( t + \frac{\mathbf{x}}{\upsilon} \right) - \omega\_1 \left( t - \frac{\mathbf{x}}{\upsilon} \right) \right] \tag{20}$$

where *zv* is wave impedance, *zv* = - *J*/*S <sup>c</sup>* = *vJ* . The analogous formula for calculation of angular velocity may also be written as follows:

$$
\omega(\mathbf{x},t) = \frac{1}{2} \Big[ \omega\_1 \Big(t + \frac{\mathbf{x}}{\upsilon} \Big) + \omega\_1 \Big(t - \frac{\mathbf{x}}{\upsilon} \Big) \Big] - \frac{1}{2z\_\upsilon} \Big[ M\_1 \Big(t + \frac{\mathbf{x}}{\upsilon} \Big) - M\_1 \Big(t - \frac{\mathbf{x}}{\upsilon} \Big) \Big].\tag{21}
$$

On the basis of the Equations (20) and (21) the moments of torsion at the shaft beginning (input) and the shaft end (output) may be expressed as follows:

$$M\_1(t) = z\_\upsilon \omega\_1(t) - z\_\upsilon \omega\_m \left(t - \frac{l}{\upsilon}\right) + M\_m \left(t - \frac{l}{\upsilon}\right) \tag{22}$$

$$M\_{\rm mf}(t) = -z\_{\rm \overline{\nu}} \omega\_{\rm mf}(t) + z\_{\overline{\nu}} \omega\_1 \left(t - \frac{l}{\upsilon}\right) + M\_1 \left(t - \frac{l}{\upsilon}\right) \tag{23}$$

where *l* is the shaft length. The abovementioned moments of torsion may also be expressed in the discrete form:

$$M\_1(j) = z\_{\overline{v}}\omega\_1(j) - z\_{\overline{v}}\omega\_{\mathfrak{m}}(j-n) + M\_{\mathfrak{m}}(j-n),\tag{24}$$

$$M\_{\rm mf}(j) = -z\_{\rm \upsilon}\omega\_{\rm mf}(j) + z\_{\rm \upsilon}\omega\_1(j-n) + M\_1(j-n),\tag{25}$$

with the initial conditions *M*1(0), *Mm*(0), ω1(0) and ω*m*(0), where: *n* is the number of numerical calculation points, *n* = *l*/*vh*, *h* is the width of the step size between the points, *j* = 0,1, ... . The number *n* is also a discrete-time of moving the mechanical wave from one end of the shaft to the other. The continuous-time of moving the wave along the shaft axis, corresponding with number *n*, is given in the Equations (22) and (23) as *l*/*v*. The optimum number of numerical calculation points is the number, the increase of which does not cause significant differences in simulation results, whereas the decrease of this number causes significant differences in simulation results. The terms from the previous steps ω(*j–n*) and *M*(*j*–*n*) in the Equations (24) and (25) may be determined using shift registers. For *j* ≤ *n*:

$$M\_1(j) = z\_v \omega\_1(j),\ M\_m(j) = -z\_v \omega\_m(j). \tag{26}$$

The angular velocities at the beginning (input) and the end (output) of the transmission shaft may be calculated on the basis of the equations of motion for the rotors of both an electric motor and working machine:

$$J\_c \frac{d\omega\_1(t)}{dt} = M\_\varepsilon(t) - D\_1 \omega\_1(t) - M\_1(t),\tag{27}$$

$$J\_L \frac{d\omega\_m(t)}{dt} = M\_m(t) - D\_m \omega\_m(t) - M\_L(t),\tag{28}$$

where *Me*(*t*) and *ML*(*t*) are motor output torque and a load torque of working machine, *Je* and *JL* are moments of rotor inertia for motor and working machine, respectively, *D*<sup>1</sup> and *Dm* are coefficients of friction defined for bearings.

The model with two lumped damping parameters *D*1*<sup>m</sup>* and *Dm*<sup>1</sup> may be adopted in order to represent the viscous friction inside the shaft. Then, the Equations (27) and (28) take the following form:

$$J\_c \frac{d\omega\_1(t)}{dt} = M\_c(t) - D\_1\omega\_1(t) - D\_{1m}(\omega\_1(t) - \omega\_m(t)) - M\_1(t),\tag{29}$$

$$J\_L \frac{d\omega\_m(t)}{dt} = M\_m(t) - D\_m \omega\_m(t) - D\_{m1}(\omega\_m(t) - \omega\_1(t)) - M\_L(t). \tag{30}$$

Such an extension of the model allows for maintaining the advantageous features resulting from the simplicity of numerical calculations.

#### **3. Simulations and Experimental Results**

#### *Research Implementation Methods*

To verify the proposed solution of the long driveshaft modelling algorithm by the use of an electric transmission line, simulation and experimental research were performed. In the first step, the algorithm computer program in C++ programming language was written by the research team. Open and free software was used to write the simulation program. The development of a simulation program independently enables the optimization of the algorithm code that is to be implemented in the future in an inverter with a digital signal processor.

The simulation contains the brushless dc motor (BLDC) model with a full control system, long driveshaft model according to Equations (24) and (25) and a direct current generator model. In the simulation, the parameters of 4 kW BLDC motor and steel transmission shaft of length 0.66 m and diameter 0.008 m, according to the real laboratory setup elements, were taken. Other parameters of the simulation model and reference units are shown in Table 1.


**Table 1.** Model parameters.

In order to determine the optimum number of numerical calculation points, the numerical simulations of the transmission shaft behaviour were carried out, in which various numbers of numerical calculation points *n* were taken into account, respectively, 10, 15, 20 and 25 points. Figure 4 shows angular velocities for beginning (input) and end (output) of the transmission shaft for different calculations variants. The end of the transmission shaft was loaded by a rotating mass *JL* <sup>=</sup> 0.11 kg·m2 and a mechanical torque of rated value.

**Figure 4.** Angular velocity at shaft input (**a**) and output (**b**) vs. time.

Figure 5; Figure 6 show respective the differences in angular velocity and angle of twist for the transmission shaft for different calculations variants.

**Figure 5.** Difference in angular velocity vs. time.

**Figure 6.** Angle of twist vs. time.

If a lower number of numerical calculation points is chosen (10 and 15) then the effect similar to the shortening of the shaft occurs (slightly less amplitude of oscillation and slightly less natural frequency as well as the decrease of steady angle of twist after 8 s from starting the system at the time *t* = 1 s).

For 20 and 25 points, the differences in amplitude of oscillation, natural frequency and steady angle of twist are not significant. It can be read from the graphs. For example, the difference in amplitude of oscillation and difference in steady angle of twist for 10 and 25 points, respectively, is in a range of 1 degree, which is around 5% in respect to the steady angle of twist, whereas for 20 and 25 points, respectively, this difference is a fraction of a percentage. The differences can disappear for an infinite number of points; therefore, a balance between the accuracy of numerical calculations and the speed of numerical calculations should be found.

In the next step, the simulation results were verified by the experiment. Experimental investigations were made using a test stand including a BLDC motor coupled with a DC generator by a steel long driveshaft. The parameters of the shaft (length, diameter) in the experimental setup was designed in accordance with the power of experimental setup components. This setup allows for connecting the shafts of various diameters, which allow for a good representation of the parameters of real-world systems. The laboratory setup is shown in Figure 7. In addition, in the figure, a simplified diagram of an experimental stand layout corresponding to the actual configuration was presented. The test stand is fully configurable; therefore, it is possible to mount additional steel discs at the shaft on the DC generator side that allow to change the moment of inertia of the load. Two rotary incremental encoders with a resolution of 3600 pulse/rev installed on the transmission shaft ends, are used to measure angular displacement and angular velocities. Hall effects transducers were used to measure motor phase currents. The measuring signals from the transducers are sent to a laboratory computer equipped with two multifunctional I/O devices PCIe-6351 from National Instruments.

Laboratory Virtual Instrument Engineering Workbench (LabVIEW) was used for data collection, whereas a computer program to visualise the measurement datasets was created by the authors. The test stand is fully automated and allows testing of the driveshaft with different geometries and made of different materials.

All presented time waveforms correspond to starting the motor and the operation under rated load after 5.5 s. The characteristic time waveforms for angular velocities at shaft input and output, the difference in angular velocity, angle of twist and phase current of BLDC motor, are shown in Figure 8, Figure 9, Figure 10, Figure 11, Figure 12.

**Figure 7.** Stand for testing a drive system with elastic coupling.

**Figure 8.** Angular velocity at shaft input vs. time.

**Figure 9.** Angular velocity at shaft output vs. time.

**Figure 10.** Difference in angular velocity vs. time.

**Figure 12.** Phase current of BLDC motor under load vs. time.

High similarity of both the experimental results confirms good accuracy of the proposed model.

#### **4. Discussion**

The results of the computer simulation and experimental tests presented in the paper deal with the modelling of the transmission shaft being a mechanical component of an electric drive. The high similarity of both experimental and simulation results confirms the high conformity of the proposed mathematical model, based on electrical and mechanical similarities with the physical object. This conformity is crucial from a practical point of view for design activities.

The telegrapher's equations and their d'Alembert solution are proposed in the article for the mathematical description of the transmission shaft. If the boundary conditions for the two points at the ends of the shaft are adopted, the solution is determined only for these points, without the need to track the processes within the shaft. If the shaft is part of a complex system, including an electric motor, power converter, mechanical coupling and working machine, then only the variables at the shaft ends can be taken into account, without the need to estimate the phenomena for any other value of shaft axis coordinate x. This assumption leads to the very simple equations of the discrete model of the lossless transmission shaft. This approach is known as Bergeron diagram method. The solution of the abovementioned equations does not require numerical integration, which can cause a problem with the stability of the integration method and it is often carried out in more than one stage. As a consequence, a shorter time for carrying out the numerical calculations is required if the proposed model is used than in the case of the models based on differential equations. Thus, high accuracy of the proposed model goes hand in hand with its simplicity and, as a consequence, with a short time of numerical calculations in contrast to the modelling based on numerical integration of differential equations or analytical solution of partial differential equations. For example, in the case of the multi-mass model of the transmission shaft divided into *m* discrete elements, 2*m* − 1 ordinary differential equations have to be integrated numerically *k* times for each step, where *k* is the order of numerical integration method. The use of the proposed model, based on the d'Alembert solution for telegrapher's equations, requires for each step the single numerical calculation of 4*n* − 2 discrete equations of type X [j + 1] = X [j] and 2 discrete equations given in the article as (24) and (25), where *n* depends on the shaft length, phase velocity and width of the step size in accordance with dependency given in the article. In order to ensure high accuracy of numerical integration of ordinary differential equations, which is comparable with the accuracy of numerical calculation of the discrete equations included in the proposed model, the width of the step size for the numerical integration should be tens of times shorter than for the numerical calculation of the abovementioned discrete equations. This means that using the proposed model leads to less numerical calculations in the same unit time than for the models based on differential equations.

#### **5. Conclusions**

Long elastic elements, especially transmission shafts, may be described mathematically by multi-mass lumped parameter models based on ordinary differential equations or distributed parameter models based on partial differential equations, which may be solved analytically. However, the analytical solution is both tedious and time-consuming. Thus, partial differential equations are usually transformed into ordinary differential equations by using finite differences. In both the abovementioned cases, the ordinary differential equations are usually integrated numerically.

In the paper, the kinematic structure of the transmission shaft between the driving motor and the working mechanism is studied. The analysis is based on electrical and mechanical similarities. The equations of moments and angular velocities, analogous to the Kirchhoff's law-based equations typically applied in the mathematical analysis of branched electrical circuits, are defined. Modelling of transmission shaft based on a formal analogy between the transmission shaft and the electric transmission line is also proposed. The high conformity of the proposed mathematical model with the physical object has been confirmed by the comparison of the results of computer simulation and experimental tests.

The advantage of the proposed model is its simplicity, due to the fact that the model bases on discrete equations, which do not require numerical integration in contrast to the models based on differential equations. The presented approach to modelling of the entire electromechanical system can be easily developed towards the system containing more than one elastic coupling, electric motor or working machine.

The proposed mathematical model may be used to analyse the impact of control and load on an electromechanical system as well as the failure states associated with or close to the resonance phenomenon. Simulation results may be used to design electric drives as well as to optimise electromechanical systems during operation.

The proposed transmission shaft model based on discrete equations is an open-loop model, but nothing precludes its use in a closed-loop control system in an advanced electric drive. This will be the subject of further research work of the authors team.

**Author Contributions:** Conceptualization, A.P. and K.B.; Methodology, A.P.; Software, A.P. and K.B.; Validation, A.P., M.L., M.N. and K.B.; Formal analysis, A.P. and K.B.; Investigation, A.P., M.L., M.N. and K.B.; Resources, A.P., M.L., M.N. and K.B.; Data curation, A.P., M.L., M.N. and K.B.; Writing—Original draft preparation, A.P. and K.B.; Writing—Review and editing, A.P. and K.B.; Visualization, A.P., M.L., M.N. and K.B.; Supervision, A.P. and K.B.; Project administration, A.P. and K.B.; Funding acquisition, M.L. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **Nomenclature of All Variables**


#### **References**


© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Energies* Editorial Office E-mail: energies@mdpi.com www.mdpi.com/journal/energies

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18