An E ﬃ cient Analytical Approach to Investigate the Dynamics of a Misaligned Multirotor System

: This paper deals with the analytical investigation of a multirotor system connected by a ﬂexible coupling subjected to dynamic angular misalignment. A novel analytical approach is employed for this purpose, namely the optimal auxiliary functions method, which proved to be successfully used to obtain explicit analytical solutions to a system of strongly nonlinear di ﬀ erential equations with variable coe ﬃ cients, useful in dynamic analysis of the considered multirotor system. The proposed procedure proved to be very e ﬃ cient in practice for solving complicated nonlinear vibration problems.


Introduction
Misalignment is one of the most commonly observed faults in rotating mechanical systems [1][2][3], being often a major cause of machinery vibration and representing a great concern to designers and maintenance engineers, since perfect alignment cannot be achieved in practice.
Even if a perfect alignment is initially realized, in many cases this condition cannot not be maintained for long time during the working regime due to many disturbing factors, such as base foundation movements, improper machine assembly, and thermal distortion and lubrication effects [4].
Beside many other possible defects in rotating machines (such as rotor unbalance, rotor bends, cracks, rubs), misalignments generate important dynamical problems which should be avoided, especially when higher speeds and higher loads are envisaged. Obviously, the importance of these aspects is increased for rotating machines having a high capital cost. In principle, one could find two types of coupling misalignments, which are parallel misalignment and angular misalignments, but often a combination of these types is common in the practice of rotating machinery.
The problem of misalignment in rotating machines is both analytically and experimentally investigated. Analytical predictions are often useful in order to have a deeper insightinto the essence of the phenomenon, especially in the stage of design, while experimental investigations are performed in practice on existent mechanical systems in order to diagnose possible defects.
Both analytical and experimental approaches are useful in order to provide better functioning conditions from dynamical point of view. Basically, the analytical and numerical approaches are the most inexpensive ways to analyze and predict misalignment phenomena. However, the experimental approach remains the way to validate these predictions. A combination of these approaches could be sometimes useful.
The need for a complete understanding of the phenomena related to misalignment has generated a large body of literature in the last few decades. Al-Hussain [5] examined two rotors which are found in connection with a flexible mechanical coupling having an angular misalignment. In another work, Al-Hussain and Redmond [6] examined, from a theoretical and numerical point of view, the lateral and torsional responses to pure parallel misalignment of two rotating shafts revealing that the natural frequencies are excited due to parallel misalignment. Sekhar and Prabhu [7] proposed a higher order FEM model for a rotor-coupling-bearing system which incorporates coupling misalignment reaction forces and moments, which were identified in both parallel and angular misalignments. By using the proposed model, the vibration responses can be predicted.
Beside analytical and numerical approaches, abundant experimental investigations of misalignments using various techniques were reported. Pérez et al. [8] characterized the parallel misalignment in rotating machines presenting the discrete time interval measurement system (DTIMS) as an effective experimental approach to identify and to measure the angular vibration produced in case of a parallel misalignment between coupled shafts of a rotating machine.
Sinha et al. [9] illustrated a reliable method to determine the misalignment and the rotor unbalance from a single machine run-down. This approach is proved by means of experimental data corresponding to a machine with two bearings and a flexible coupling to the motor.
Experimental investigations on vibration response of misaligned rotors were conducted by Patel and Darpe [10], whoefficiently used full spectra and orbit plots to emphasize misalignment fault, proposing new misalignment diagnostics recommendations. In the existing literature, there are few publications concerned with analytical investigation of misalignment, mainly because of the complexity of the governing equations, which are difficult to be solved analytically without adopting harsh simplifying hypotheses. An efficient approach proved to be the harmonic balance method which is employed in [11] for analyzing a misaligned rotor system.
There are available in the literature many analytical methods intended to provide approximate analytical solutions to nonlinear problems-such as the multiple time scales method [12], the differential transform method (DTM) [13], the Adomiandecomposition method (ADM) [14], the variationaliteration method (VIM) [15], the homotopyanalysis method (HAM) [16], the homotopyperturbation method (HPM) [17], the optimal homotopyasymptotic method (OHAM) [18], and the optimal parametric iteration method (OPIM) [19]-some of them being applicable also for strongly nonlinear problems as is the case of misalignment problems under study in this work. Such analytical approximate approaches are widely used since exact solutions for nonlinear problems are rare [20].
Many scientists concerned with complex nonlinear problems use numerical or semi-analytical approaches to the study of advanced nonlinear problems, since these are easier to implement, but an analytical approach, even not so easy to implement, would offer a deeper insight into the physical nature of the problem and allows further developments, such as the study of stability. Analytical developments represent a constant concern of scientists and the method presented in this paper is a new analytical approximation method which provides explicit and accurate analytical solutions for complex nonlinear problems.
The present work is devoted toward an application of a new approach, the optimal auxiliary functions method (OAFM) to investigate dynamic phenomena generated by angular misalignment of two rigid rotors connected by a flexible coupling, the bearings being isotropic in dynamic performances [21]. This system, which is numerically investigated in [21], generates a strongly nonlinear problem, which is very hard to solve through classical analytical techniques (without considering small parameters).
To the knowledge of the authors, no study has presented explicit analytical solutions to such a complex problem without using simplifying hypotheses. In this work, we will develop a dynamic analysis using no simplifying hypothesis and avoiding the assumption of small parameters. Highly accurate explicit analytical solutions will be obtained for a very complicated nonlinear dynamic system, which will prove the efficiency of the proposed approach.

Governing Equations
In Figure 1 is shown the model of the multirotor system with flexible couplings which transmits torque between the two rotors. In Figure 1 is shown the model of the multirotor system with flexible couplings which transmits torque between the two rotors.
We consider that the two rotors are rigid and the bearings are isotropic, m1 and m2 are the lumped masses, k1 and k2 are bearings' stiffness, and ktis the angular stiffness of the flexible couple, the coordinates of the mass centers of discs 1 and 2 are O1(x1,y1,z1) and O2(x2,y2,z2), referring to the coordinates system Oxyz which is set up in the static equilibrium position (Figure 1). If x and y represent the displacement of the first mass at the geometric center, θ is the angle between rotors, Ω is the speed of rotation, γ is the initial phase angle, and e1 and e2 are the eccentricities of the unbalanced masses, then we have The terms related to e2 could be neglected since the offset h of disc 2 is much larger than e2. In these conditions, using the Lagrange's equation after computing the kinetic and potential energy, and considering that m1=m2=m, k1=k2=k, e1=e, then the governing equations can be written as [21] 0 t  cos  ml  t  cos  sin  kl  kx  2   t  cos  sin  ml  t  sin  cos  ml  2  t  cos  sin  ml  t  cos  cos  ml  x  m  2 0  t  sin  ml  t  sin  sin  kl  ky  2   t  sin  sin  ml  t  cos  cos  ml  2  t  sin  sin  ml  t  sin  cos  ml  y  m   We consider that the two rotors are rigid and the bearings are isotropic, m1 and m2 are the lumped masses, k1 and k2 are bearings' stiffness, and ktis the angular stiffness of the flexible couple, the coordinates of the mass centers of discs 1 and 2 are O1(x1,y1,z1) and O2(x2,y2,z2), referring to the coordinates system Oxyz which is set up in the static equilibrium position ( Figure 1). If x and y represent the displacement of the first mass at the geometric center, θ is the angle between rotors, Ω is the speed of rotation, γ is the initial phase angle, and e1 and e2 are the eccentricities of the unbalanced masses, then we have x 1 = x + e 1 cos Ωt y 1 = y + e 1 sin Ωt (1) and x 2 = x + l sin θ cos Ωt + e 2 cos θ cos(Ωt + γ) y 2 = y + l sin θ sin Ωt + e 2 cos θ sin(Ωt + γ) The terms related to e2 could be neglected since the offset h of disc 2 is much larger than e2. In these conditions, using the Lagrange's equation after computing the kinetic and potential energy, and considering that m1 = m2 = m, k1 = k2 = k, e1 = e, then the governing equations can be written as [21] 2m ..
In this work, we present a new procedure, namely the optimal auxiliary functions method to study the nonlinear vibrations of the multirotor system. The accuracy of the results obtained using through the proposed procedure is illustrated by numerical simulations developed to validate the obtained analytical results. This approach is independent of any small or large parameters.

Basics of the Optimal Auxiliary Functions Method
In order to show the basics of OAFM [22], we consider a nonlinear equation in the general form L[u k (τ)] + N[u k (τ)] = 0 , k = 1, 2, . . . , n (13) in which L represents the linear differential operator and N is the nonlinear differential operator, τ represents the independent variable and u k (τ) is an unknown function (k = 1, 2, ..., n). The initial/boundary conditions are It is well-known that an exact solution for strongly nonlinear differential equations of type (13) and (14) is very hard to be identified. To obtain an approximate solution u k (τ) we consider that this solution can be expressed in the form where u k0 (τ) and u k1 (τ) which represent the initial approximation and the first approximation will be determined as described below.
Substituting (15) into (13) one obtains where Ci, i = 1, 2, ..., s are unknown parameters at this moment. The initial approximation can be obtained from the linear equations and the first approximation from the Equations (16) and (17) The nonlinear term in Equation (18) is expanded in the form however, to avoid the difficulties appearing in solving the nonlinear Equation (19) and implicitly in finding the approximate solution u k (τ, C i ), instead of the last term arising in Equation (18), we propose another expression, such that Equation (18) can be rewritten in the form where A1 and A2 are the arbitrary auxiliary functions which depend on u k0 (τ) and on the convergence-control parameters Cj and Cl, j = 1, 2, ..., p, l = p + 1, p + 2, ... s, and F(N[u k0 (τ)]) are functions depending only on some expressions which appear in the component of N[u k0 (τ)]. It should be emphasized that the linear operator L, the auxiliary functions A1 and A2 (namely optimal auxiliary functions) and F(N[u k0 (τ)]) are not unique, but A1 and A2 are of the same form like u k0 (τ). More explicitly, for example, if u k0 (τ) is a polynomial function, then A 1 u k0 (τ), C j and A 2 u k0 (τ), C q are sums of polynomial functions. In case of u k0 (τ) is an exponential function then A1 and A2 are sums of exponential functions. In the case when u k0 (τ) is a trigonometric function then the auxiliary functions are sums of trigonometric functions, and so on. If N[u k0 (τ)] = 0, one obtains u k0 (τ) as an exact solution of Equation (1), but this rarely happens for nonlinear equations.
The values of the convergence-control parameters Cj and Cq may be optimally determined by means of different methods such as the least square method, collocation method, Galerkin method, Ritz method, but the preferred one should be minimizing the square residual error. If then using the operator in which the approximate solution u k (τ, C i ) is given by Equation (15), the unknown parameters C1, C2,...,Cs can be obtained from the system of algebraic equations In this way, after finding the optimal values of the initially unknown parameters Ci, the approximate solution is well-determined. Hence, our procedure contains the auxiliary functions A1 and A2 which ensure a rigorous way to adjust and control the convergence of the approximate solutions u k (τ, C i ).
Here, one should underline the importance of a properly choosing functions A1 and A2 which appear in the construction of the solution u k1 (τ, C i ).

Numerical Results
In order to show the validity of the OAFM, Equations (10)- (12) have been analyzed for the following cases:

Case 1
For the first case we consider A = 0.06; α = 0.07; β = 0.08; ω0 = 1.2; ωt = 0.02; L = 20; E = 0.05. In this case, using the procedure of minimizing the residual error [18], one obtains the optimal values of the convergence-control parameters and approximate frequency as It is important to mention that the value of the frequency obtained through numerical integration is ω NUM =0.662478, and should be emphasized that the OAFM lead to high precision results. Figure 2 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach. Figure 3 shows the residual (21) corresponding to the first approximate solution θ. It is important to mention that the value of the frequency obtained through numerical integration is ωNUM=0.662478, and should be emphasized that the OAFM lead to high precision results. Figure 2 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach. Figure 3 shows the residual (21) corresponding to the first approximate solution θ .  One can observe that the approximate solutions are in very good agreement with the numerical integration solutions Also, the results obtained by OAFM are remarkably good when compared to the numerical integration results in the phase plane. All these prove the excellent accuracy of the analytical results, as shown in Figures 2-8.  It is important to mention that the value of the frequency obtained through numerical integration is ωNUM=0.662478, and should be emphasized that the OAFM lead to high precision results. Figure 2 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach. Figure 3 shows the residual (21) corresponding to the first approximate solution θ .  One can observe that the approximate solutions are in very good agreement with the numerical integration solutions Also, the results obtained by OAFM are remarkably good when compared to the numerical integration results in the phase plane. All these prove the excellent accuracy of the analytical results, as shown in Figures 2-8. One can observe that the approximate solutions are in very good agreement with the numerical integration solutions Also, the results obtained by OAFM are remarkably good when compared to the numerical integration results in the phase plane. All these prove the excellent accuracy of the analytical results, as shown in Figures 2-8.

Case 2
In the second case it is considered A=0.06; α=0.07; β=0.08; ω0=1.2; ωt=0.1; L=20; E=0.05. In this case, following the same procedure, one obtains the following optimal values of the convergence-control parameters and frequency In this case, the value obtained through numerical integration for the frequency is ωNUM=0.676796, which means that a very good approximation of the frequency is obtained through OAFM. Figure 9 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach while Figure 10 shows the residual corresponding to the first approximate solution θ in the second case.
It is observed in this case as well that we obtain an excellent approximation of the solution ) ( τ X and also, the results plotted in the phase plane emphasize the excellent accuracy as shown in Figures  9-15.

Case 2
In the second case it is considered A = 0.06; α = 0.07; β = 0.08; ω0 = 1.2; ωt = 0.1; L = 20; E = 0.05. In this case, following the same procedure, one obtains the following optimal values of the convergence-control parameters and frequency C 1 = −1.0405329739831501; C 2 = −0.0023584330911787105; C 3 = −1.091819177088501; In this case, the value obtained through numerical integration for the frequency is ω NUM = 0.676796, which means that a very good approximation of the frequency is obtained through OAFM. Figure 9 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach while Figure 10 shows the residual corresponding to the first approximate solution θ in the second case.
In this case, the value obtained through numerical integration for the frequency is ωNUM=0.676796, which means that a very good approximation of the frequency is obtained through OAFM. Figure 9 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach while Figure 10 shows the residual corresponding to the first approximate solution θ in the second case.
It is observed in this case as well that we obtain an excellent approximation of the solution ) ( τ X and also, the results plotted in the phase plane emphasize the excellent accuracy as shown in Figures  9-15. In this case, the value obtained through numerical integration for the frequency is ωNUM=0.676796, which means that a very good approximation of the frequency is obtained through OAFM. Figure 9 presents a comparison between the solution (49) and the numerical solution obtained by means of a fourth order Runge-Kutta approach while Figure 10 shows the residual corresponding to the first approximate solution θ in the second case.
It is observed in this case as well that we obtain an excellent approximation of the solution ) ( τ X and also, the results plotted in the phase plane emphasize the excellent accuracy as shown in Figures  9-15.  It is observed in this case as well that we obtain an excellent approximation of the solution X(τ) and also, the results plotted in the phase plane emphasize the excellent accuracy as shown in

Case 3
In The numerical integration result in this case for the frequency is ωNUM=0.662157, and therefore the error of the approximate frequency obtained by means of OAFM in this case is 0.04%. In the same case, in [21], the approximate frequency obtained by using the method of multiple scales is ωMMS=0.664, having an error of 0.27%, which emphasizes the advantage of OAFM. The key innovation of OAFM, which ensures such a good precision, relies on the involvement of the auxiliary functions, which have a decisive role in simplifying computations, and the involvement of the convergence-control parameters, whose optimal values are rigorously identified using various mathematical procedures in order to ensure a fast convergence, even after the first iteration, avoiding simplifying hypotheses or small parameters.
Comparison between our approximate solutions and numerical integration results for the Cases 5.1-5.3 are presented in Figures 2-22. From these figures one can see that the results obtained using the present procedure are almost identical with those obtained by means of numerical integration using a fourth-order Runge-Kuttamethod implemented by means of Mathematica software, which  The numerical integration result in this case for the frequency is ω NUM = 0.662157, and therefore the error of the approximate frequency obtained by means of OAFM in this case is 0.04%. In the same case, in [21], the approximate frequency obtained by using the method of multiple scales is ω MMS = 0.664, having an error of 0.27%, which emphasizes the advantage of OAFM. The key innovation of OAFM, which ensures such a good precision, relies on the involvement of the auxiliary functions, which have a decisive role in simplifying computations, and the involvement of the convergence-control parameters, whose optimal values are rigorously identified using various mathematical procedures in order to ensure a fast convergence, even after the first iteration, avoiding simplifying hypotheses or small parameters.
Comparison between our approximate solutions and numerical integration results for the Cases 5.1-5.3 are presented in Figures 2-22. From these figures one can see that the results obtained using the present procedure are almost identical with those obtained by means of numerical integration using a fourth-order Runge-Kuttamethod implemented by means of Mathematica software, which validate the analytical approach. The errors between the analytical and numerical solutions are remarkable good. Also, it can be observed that the errors between the analytical and numerical results for the frequencies are very good. Some deviations are observed in the phase plane, especially for the last considered example, since this is not an exact solution, but an approximate explicit solution with a good error of approximation in the presence of strong nonlinearity. The precision of the approximate solution is also demonstrated by the small residual obtained by replacing the approximate solution into the original governing equation.
The numerical integration result in this case for the frequency is ωNUM=0.662157, and therefore the error of the approximate frequency obtained by means of OAFM in this case is 0.04%. In the same case, in [21], the approximate frequency obtained by using the method of multiple scales is ωMMS=0.664, having an error of 0.27%, which emphasizes the advantage of OAFM. The key innovation of OAFM, which ensures such a good precision, relies on the involvement of the auxiliary functions, which have a decisive role in simplifying computations, and the involvement of the convergence-control parameters, whose optimal values are rigorously identified using various mathematical procedures in order to ensure a fast convergence, even after the first iteration, avoiding simplifying hypotheses or small parameters.
Comparison between our approximate solutions and numerical integration results for the Cases 5.1-5.3 are presented in Figures 2-22. From these figures one can see that the results obtained using the present procedure are almost identical with those obtained by means of numerical integration using a fourth-order Runge-Kuttamethod implemented by means of Mathematica software, which validate the analytical approach. The errors between the analytical and numerical solutions are remarkable good. Also, it can be observed that the errors between the analytical and numerical results for the frequencies are very good. Some deviations are observed in the phase plane, especially for the last considered example, since this is not an exact solution, but an approximate explicit solution with a good error of approximation in the presence of strong nonlinearity. The precision of the approximate solution is also demonstrated by the small residual obtained by replacing the approximate solution into the original governing equation.

Conclusions
The present paper proposes a new analytical approximate procedure to solve the strongly nonlinear differential equations with variable coefficients governing the motion of the multirotor system with a flexible coupling on spring supports. We have considered the effects of dynamic angular misalignment and mass eccentricity. The system of nonlinear equations is analytically

Conclusions
The present paper proposes a new analytical approximate procedure to solve the strongly nonlinear differential equations with variable coefficients governing the motion of the multirotor system with a flexible coupling on spring supports. We have considered the effects of dynamic angular misalignment and mass eccentricity. The system of nonlinear equations is analytically

Conclusions
The present paper proposes a new analytical approximate procedure to solve the strongly nonlinear differential equations with variable coefficients governing the motion of the multirotor system with a flexible coupling on spring supports. We have considered the effects of dynamic angular misalignment and mass eccentricity. The system of nonlinear equations is analytically solved by means of a new method, namely the optimal auxiliary functions method (OAFM) and extensive numerical simulations were performed in order to validate the proposed approach.
We remark that the angular frequency ω is increased by increasing the non-dimensional angular frequency ωt (cases 2 and 3) but the frequency is insensitive to the increasing of the non-dimensional mass eccentricity (cases 2 and 3). It is observed that the variation of the angle θ between the two misaligned rotors is harmonic having the period T = 2π/ω.
The proposed technique, OAFM, accelerates the convergence of the approximate analytical solutions obtained for the multirotor system with a flexible coupling taking into consideration the effect of the dynamic angular misalignment. Also, this method leads to very good values for the approximate frequencies in all considered cases. A major advantage of the OAFM is that it is independent of any small or large parameters. The construction of the first iteration is different from any other approach mainly due to the optimal auxiliary functions depending on some initially unknown convergence-control parameters Ci, which are rigorously identified.
It should be emphasized that in our construction given by Equations (13) and (20), especially for very complicated equations such as (10), (11), and (12) it is not necessary the presence of the entire nonlinear functions N and therefore this construction lead to a considerable simplification of the calculus of the first approximation. On the other hand, the auxiliary functions Ai and Ai + 1, i = 1,2,3, have a decisive role in solving the problem.
The convergence-control parameters Ci involved in the auxiliary functions ensure a fast convergence of the approximate solutions using only the first iteration, their optimal values being identified using rigorous mathematical criterions. Consequently, very accurate explicit analytical solutions were obtained through the proposed procedure for complicated governing equations.
The most important advantage of the OAFM is represented by the feature related to optimal control of the convergence of solutions by means of the auxiliary functions Ai and Ai + 1. The approximate analytical solutions are in excellent agreement with the numerical integration results, which underlines the validity of the proposed procedure and on the other hand, it proves that this procedure is very efficient in practice.