Numerical Study for Comparison of Pseudo Modal and Direct Method in Predicting Critical Speed of Coaxial Dual Rotor System

Machines with high speed rotation, high power, and lighter weight are the challenge in design of modern rotation machinery. In analysis of coaxial dual rotor system, the problem leads to a multi-degree of freedom system (1multi-DOF). The resulting equation of motion cannot be solved analytically. For multi-DOF system, numerical method has to be implemented. In this study, the equations of motion of co-axial dual rotor system was solved numericallyusing MATLAB programming language. This study compares the result of two methods used: the pseudo-modal method and the direct method, with the aim of getting the most effective method in predicting critical rotational speedby plotting Campbell diagram. It was foundthat the difference in the value of critical speeds between Pseudo Modal and the Direct Method is 0%-9,56% for the first 6 lowest frequencies and 0.09%-3.58% for the first 8 lowest frequencies.


Introduction
Machines with high speed rotational, high power, and lighter weight are the tendencies in design of rotational machinery in modern engineering. This develops in line with the increasing ability in prediction of rotor dynamic characteristics. The use of dual rotor systems (coaxial rotors) can be found in many applications such as aviation (helicopter), turbine and other fields. To get a machine with high speed, large power and light weight, the designers of dynamic rotor systems require certain techniques in predicting the dynamic characteristics of the rotor. All the dynamic characteristics of the rotor are closely related to the rotation of the rotor, so the direction of the mode is very important in the dynamic rotor. When the rotor system starts to spin, two different modes will appear, known as the forward and backward modes. The modal characteristics of a dynamic rotor system that are related to the forward and backward modes, differ from the conventional characteristics of the non-rotational structure [1]. The application of the classical modal testing method has developed well to identify the modal parameters of various structures, except for rotation machines. In analysing the dynamic characteristics of rotor, there are several methods that might be used such as the finite element method, the transfer matrix method, the time series numerical integration method and the Hubolt method for modelling the rotor shaft system and obtain faster computational time. Budiwantoro [2] and also Huang [3] investigated the dynamic characteristics of coaxial dual rotor systems using finite element method. The results obtained was in form of critical sonar velocity, transient displacement, and trajectory diagrams. Yanto and Hidayat [4] carried out research on the dynamic behaviour of 3D rotor shaft system (1 level) using finite element analysis and obtained results of a critical rotation variation of 192 frequencies. In an analysis of rotor dynamics, Hutahaen [5] implemented the substructure method to analyse the dynamic behaviour of the rotor (turbine 1 level). His study focused on analysing the dynamic characteristics of the rotor in form of natural frequencies plotted in Campbell diagrams. Determination of the direction of whirling orbit of turbines having high speed rotation always has its own challenges. Simulations using finite element method with the help of ANSYS software have been carried out by Jahromi [6] to determine the direction of whirling orbit and calculate the gyroscopic effect. The numerical results obtained was then verified experimentally to confirm the phenomenon of whirling orbit direction.
The pseudo modal method in a computational program was used by Zulhendri et al., to analyse the dynamic characteristics of a coaxial dual rotor system. The method was also used by Mahartana [8] for modeling and simulation of rotor balancing to see the increase in work cycle of the dynamic response.
This paper presents results of numerical study on the dynamic characteristics of a coaxial dual rotor system in the form of natural frequency (the Campbell diagram) and response to mass unbalance. This study compared the results obtained using the pseudo modal and the direct method. The computer program that has been compiled using MATLAB would then be validated using experimental studies published elsewhere.

Methods
By using the differential equation of motion, the dynamic characteristics of a coaxial dual rotor system can be solved. The discretization method can be used for the equation of the coaxial dual rotor system. By applying the Lagrange equation, a linear differential equation which is useful can be obtained.

Equation of Motion
The equation of motion of the rotor system can be written as follows: By completing Equation (1) a natural frequency can be obtained as a function of rotational speed of the rotor, Campbell's diagram, critical rotation, and excitation response, such as unbalanced mass response. For tiered rotor shaft system as shown in Figure 1, the equation of motion can be written as: The three matrices , , and are obtained by arranging matrices 1 , 2 , 1 (Ω 1 ), C 2 (Ω 2 ), and on bearings and seals In a coaxial dual rotor, analysis generally leads to a multi-degree of freedom system (multi-DOF), and to solve the equation of motion, a numerical method must be used. There are two methods here used for solving the rotor system equation, namely the pseudo-modal method and the direct method.

Pseudo Modal Method
The pseudo-modal method can be used for real cases in industry. In the pseudo-modal method, the damping value C in Equation (1) is assumed to be small, so the mechanism is unknown and for that reason the damping is not installed in the rotor. Equation (1) can be defined as a modal base equation [9] as given below:

Figure 1. Coaxial dual rotor design system
The value of is the mass matrix and K* is the stiffness matrix obtained from the value of , where the values of and can be ignored as indicated by the bearing being removed. The lowest number of first vibration modes (N), ( ), 1 ,….. from Equation (3) can be obtained by iteration technique, and in matrix form can be written as in Equation (4): From Equation (3) the following new equation can be obtained.

= (5)
Substitute Equation (5) into the equation of motion of the rotor shaft system and multiply it by .
The diagonal matrix is added to the modal damping , which is obtained from the analogy of the single degree of freedom (SDOF) damped spring mass system: with the value of the modal damping factor obtained from the designer's experience. So that the solution of the equation of motion (1) without external force can be written as: The value of p can be assumed in the form: The values of m, c and k can be written as: The stiffness * * is expressed as: Equation (9) can be written in the form of a matrix equation, as in Equation (14).
The solution to the problem of eigenvalue or natural frequency and vibration mode of equation (14) can be obtained in complex quantities.
The value of is a natural frequency and is a viscous damping factor. If equation (15) > 0, the system is not stable. From the results obtained, Campbell's diagram can be drawn showing critical areas so that the critical velocity can be known. Unbalanced mass response may be obtained from the following equation: with the values 2 and 3 can be expressed as follows: Assuming an answer for under-damped vibration as follows: By substituting equation (19) into equation (16) the following matrix form will be obtained: The solution of equation (20) and from equation (19), the following displacement vectors are obtained:

Direct Method
The difference between the pseudo modal method and the direct method in solving the coupled equation of motion is that there is no equation reduction in the direct method. In the direct method, the answers are considered in the form of equation below: = Δ (22) When Equation (22) and its derivative are then substituted into Equation (2), the following equation is obtained. ( Equation (23) can be written in the form of state space.
Equation (24) is eigenvalue problem. The first two eigenvalues and corresponding eigenvector can be obtained in complex quantities as given in Equation (25) and (26) respectively: The eigenvector values can be seen in the following equation.
Campbell's diagram of the natural frequency of the system can be obtained by solving those equation.

Results and Discussion
Equation of motion of coaxial dual rotor system can be expressed as in Equation (2). The matrix values , , and in Equation (2) are obtained by arranging the matrix (disk mass matrix and inner axis), and (disk mass matrix and outer shaft). While the value of matrix is arranged based on the gyroscopic effect of the disk, the shaft (inside and outside), and the effect of bearings and seals attenuation. The matrix is obtained from the axial force matrix (inside and outside), as well as bearings and seals thickness. The form of Equation (2) can be arranged into the global equation matrix as follows: The coaxial dual rotor is a system that has two shafts that coincide on both axes and experience rotation simultaneously with different or equal rotational speeds. The two shaft systems are connected by intershaft bearings as shown in Figure 1. The model of the graded rotor shaft system can be seen in Figure1. Data and dimensions of the coaxial dual rotor system are given in Tables 1 and 2, characteristics of bearings and seals are given in Table 3, and material properties of the shaft are given in Table 4:   Table 4. Material properties of the shaft By using the finite element method to model the system, the rotor shaft system is divided into 7 elements (inside) and 4 elements (outside), so that there are 11 elements in modelling coaxial dual rotor system as shown in Figure 2. The location of nodes on the inner and outer shaft are given in Tables 5 and 6 respectively. Data for finite element is given in Table 7.  Numerical studies were carried out to compile a program used in predicting the dynamic characteristics of coaxial dual rotor system. In this research, the program was written using MATLAB software. Two solution methods were written for solving the equation of motion: the the pseudo modal and the direct method The speed rotation was in the range of 0 to 14,000 rpm. The first six lowest frequencies computed using the pseudo modal and the direct method is shown in Figure 3 and 4 respectively. Figure 5 and 6 show the first eight lowest frequencies obtained using the pseudo modal and the direct method respectively.    The values of Critical speed obtained using the Pseudo Modal method and the Direct Method are given in Table 8 for the first six lowest frequency and in Table 9 for the first eight frequency.  It can be observed from Table 8 and 9 that there are difference in values of critical rotation speed between the first six lowest frequency and the first eight lowest frequency if the Pseudo Modal method is used. This is due to using the Pseudo Modal method, the value of the matrix damping factor C is not included in the global equation system as shown in Equation (2). This causes the critical rotation speed value obtained by the Pseudo Modal method to be higher when compared to the critical rotation speed value obtained using the Direct method. The percentage change in the value of critical speed using the Pseudo Modal method between the first six lowest frequencies and the first eight lowest frequencies are given in Table 10. If the results obtained using the Pseudo Modal method are compared with the results obtained using the Direct method, the results can be seen in Table 11 for the first six lowest frequencies, and Table 12 for the first eight lowest frequencies. It should be noted that the Pseudo Modal method in the process of solving the equations of the rotor shaft system uses the number of modes (NoM) vibrations that want to be analyzed. Meanwhile, the Direct method uses the maximum number of elements (NMAX), on the rotor shaft. These differences can be seen in the following MATLAB program script.

Conclusion
The Pseudo Modal method and the Direct Method can be used to predict the critical rotation speed of a coaxial dual rotor shaft system. By using the Pseudo Modal method, the value of the damping matrix C is considered small so that it was not taken into account in the completion of the motion system. Assuming small value of C matrix, the value of critical speed obtained using the Pseudo Modal method is greater than those obtained using the Direct method. There is a difference in the value of critical speed using the Pseudo Modal method between the number of vibrational mode (NoM) between the first six lowest frequencies and the first eight lowest frequencies. Further research is needed on experimental study of coaxial dual rotor shafts to validate the best results between the Pseudo Modal method and the Direct Method.