A High-Precision Random Forest-Based Maximum Lyapunov Exponent Prediction Model for Spherical Porous Gas Bearing Systems

Spherical porous air bearing (SPAB) systems have been extensively used in various mechanical engineering applications. SPABs are promising materials in high-rotational speed, high-precision, and high-stiffness instruments. In SPAB systems, a rotor is supported by gas bearings, which provides higher rotational speed and lower heat generation environment than oil bearings do. Furthermore, SPAB does not cause deformation. Although, the supporting force of gas bearings is less, their stability is better than that of oil films. However, because the pressure distribution in the gas films is nonlinear, they are prone to failure at specific critical speeds, rotor imbalances, or inappropriate operations, which results in nonperiodic or chaotic motion and causes structural fatigue to the system. To understand and control the operating conditions of the SPAB systems during the nonperiodic motion, first, the governing equations of the SPAB system were solved to obtain the dynamic behavior of the rotor center. Then, the performance of the SPAB system were examined under different operating conditions by generating the maximum Lyapunov exponents (MLEs). However, the calculation process of MLE is extremely time consuming and complex. To solve this problem efficiently, a high-precision machine learning (ML)–based MLE prediction model was proposed in this study. The results show that the training process can be finished within few minutes, and the prediction process is able to be completed within few seconds. Meanwhile, the results demonstrate the merit of using the machine learning method for solving the MLE prediction problem and shorten the calculation time significantly. The proposed prediction model achieves excellent prediction outcome and it is more efficient and precise than traditional iteration scheme for the calculation of MLE. The feasibility of the proposed model is validated and the results also are the major contribution of this study.


I. INTRODUCTION
Spherical porous air bearing (SPAB) systems exhibit several attractive features such as low-noise rotation, high load capacity, improved damping properties, and zero friction. SPABs are simpler and cheaper than externally pressurized bearings. SPABs generate less heat and provide higher The associate editor coordinating the review of this manuscript and approving it for publication was Yunhua Li. accuracy than oil bearings do. SPABs are applied in highspeed spindles and machine tools.
In 1964, Sneck and Yen [1] considered a one-dimensional flow for radial direction in a porous medium and developed a perturbation solution for a finite journal bearing. They subsequently experimentally verified the theoretical solutions in 1965 [2].
In 2001, Sinha et al. [3] analyzed externally pressurized conical bearings with a constant porous gap when the slider VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ was rotated with a uniform angular velocity. The governing equations of the system included coupled momentum and energy equations, which were analyzed using finite difference method (FDM) and various bearing characteristics. The results revealed that the inlet pressure decreased remarkably for highly porous surfaces, thus reducing the load capacity of the bearing. However, the torque of the system with respect to variation in permeability remained unaffected.
Rashidi et al. [4], [5] studied the nonlinear dynamic behavior of rigid rotors supported by noncircular air journal bearing systems with two-, three-, and four-lobe bearings In their study, the Reynolds equation for air pressure distribution was analyzed using the finite element method and the rotor dynamic equations were calculated using the Runge-Kutta method. The numerical results revealed that the bifurcation phenomenon with different bearing numbers and rotor masses including periodic and nonperiodic motions and chaos at different operational situations. The bearing number and rotor mass were the major parameters for the bearing system.
The mathematical model of SPAB includes determining the Reynolds equation to obtain the gas pressure distribution and subsequently combining it with the rotor dynamics of the rotor to solve the displacement and vibrations of the rotor center. However, it is not easy to distinguish whether the motion of the rotor center is periodic. Dynamic orbits, power spectra, bifurcation diagram, and Poincaré maps are used to confirm the types of behavior and the changes of rotor motions as the rotor mass or bearing number is increased [6]. The maximum Lyapunov exponent (MLE) is applied to guarantee the occurrence of chaos and should be avoided the chaotic area as the bearing system is operated. The MLE is complex and influenced by coupled parameters such as the rotor mass and bearing number. Therefore, its calculation is time consuming and prediction of the correct chaotic range is unreliable.
The major contributions of this paper are as follows: 1) providing the derivations of the MLE calculation process; 2) proposing the MLE prediction model to reduce the computation time; 3) validating the feasibility of the proposed framework; 4) comparing the performances of the machine models in this topic.
The system architecture and derivations of the MLE calculation process are introduced in Section II. The proposed method is presented in Section III. Section IV presents the experimental results and suggestions. Finally, Section V presents the conclusions.

II. SYSTEM ARCHITECTURE
In the proposed SPAB model, the gas flow is assumed to be isothermal and the gas viscosity is constant. The pressure distribution in the gas film between the shaft and the bushing is modeled using the following Reynolds equation: The corresponding boundary conditions are represented in the following expression: where . P is the dimensionless pressure corresponding to the atmospheric pressure, H is the dimensionless gap between the rotating shaft and the bushing corresponding to the radial clearance (C r ), µ is the gas viscosity, θ and φ are the coordinates, is the bearing number, σ is the squeeze number, ε z is the eccentricity in the axial direction, ε r is the eccentricity in the radial direction, and r is the bearing radius.
The proposed system comprises a perfectly balanced flexible rotor supported symmetrically on two identical SPABs mounted on rigid pedestals. Because the rotor is assumed to be perfectly balanced and the SPAB is symmetric about its central axes, the system is confined to a single bearing supporting a rotor with two degrees of translatory oscillation in the transverse plane. The equation of the motion of the rotor in the transient state can be expressed as follows: in which F elx and F ely are the external loading force in the horizontal and vertical directions, respectively. Here, F gfx and F gfy are the internal force in the horizontal and vertical directions, respectively, estimated by integrating the air pressure distribution.
For calculating the air pressure in (1), an iterative numerical process for each time step was applied using a hybrid numerical method [17] that integrates the differential transformation method and FDM. Differential transformation is one of the most widely used techniques for solving differential equations because of its rapid convergence rate and minimal calculation error. An advantage of this method over the integral transformation approach is that this method can be used to solve nonlinear differential equations.
For solving the Reynolds equation for the current air bearing system, the differential transformation method was used for performing transformation with respect to the time domain τ . Thus, we obtain the following expression: where The FDM was then used to discretize (4) with respect to the coordinate directions. Equation (4) was discretized using the second-order-accurate central-difference scheme for both the first and the second derivatives. The rotor center behavior was analyzed using the iterative process including determination of the acceleration, velocity, and displacement; calculation of the change of the air gap between the rotor and the bushing; and integrating the air pressure to form the internal force. The rotor displacements generated through the iterative procedure were then applied to analyze the bifurcation diagram and MLE. The dynamic responses of the SPAB system over the specific ranges of the rotor mass and bearing number were also examined.
The bifurcation phenomenon of the rotor motion under the operating ranges of the rotor mass and bearing number are displayed in Figs. 1, 2, 3, and 4. The calculation process of the overall iteration is shown in Fig. 5.
The MLEs were used for the interpretation and verification of the chaotic behavior. Figs. 6 and 7 depict that when the MLE tends to zero, the system is nonchaotic. Fig. 6 indicates that as the rotor mass increases to 22.66 kg, the MLE is greater than zero, which proves that the motion state of the system is chaotic. The results also revealed that the region of the red curve of chaos (index greater than zero) occurs in the interval of 22.66 m r < 23.52 kg, and the rest of the blue curve occurs in the stable region. The results in Fig. 7 reveal that the region of chaos (MLE > zero) occurs in the interval 1.0 < 2.4.

III. PROPOSED METHOD
In this section, the adopted ML model is introduced briefly.
To address the MLE prediction problem, this task can be    treated as a regression problem. The goal of the proposed framework is training the model with just few training data and using the trained model to predict the MLE values in the unknown situation of the rotor mass and bearing number. Therefore, several ML models can be adopted to solve this regression problem. In this study, the MLP, SVM, DT, and RF were adopted to obtain the prediction result, and the performance of the models was compared. MLP is a basic neural network structure. The working of the human brain is the inspiration for the MLP concept. Many neurons are included in this model. The training process of the MLP can be divided into two parts. The first part is the feed forward process, and the second part is the  backpropagation process. All the weights between the hidden layers are adjusted in the backpropagation process. The SVM is a famous ML model. By adopting the suitable kernel, the SVM can be treated as a classifier or a regression model. However, because the SVM is a supervised learning method, it requires sufficient training samples. The structure of the DT is depicted in Fig. 8(a). The DT is a tree-like model. First, the input data can be classified using the root node. The path from the root to the leaf of the DT can determine the prediction result. Although DT is simple, it provides excellent classification and regression performance. The structure of the RF is represented in Fig. 8(b). The RF model consist of several DT models. Therefore, the RF is an ensemble model. The RF considers all the prediction results of DTs and then provides the final decision. Because several DTs are combined in an RF model, generally, the performance of RF is better than that of DT.
The training flowchart of the proposed method is illustrated in Fig. 9. After the training data (MLE map) is obtained, the data preprocessing procedure is executed, and the number of training sample will also be defined. In this progress, all the values are normalized within [0,1]. This is a crucial step for the ML model because all the input dimension weights must be maintained at the same influence level. Then, the ML model is created and initialized. Next, the training process is initiated. In this process, if the predetermined convergence condition is satisfied, then the training process is complete. Otherwise, the training process continues to pursue a better result. Finally, the performance of the trained model is evaluated. For this part, after the MLE map prediction is completed, the input value and the prediction are denormalized to fit the original scale. The prediction error in terms of rootmean-square error (RMSE) and mean absolute error (MAE) is also calculated. Then, in order to give a suggestion for the users, the statistical analysis and the performance ranking calculation will also be done. Besides, RMSE and MAE are defined as follows: where M is the number of the samples, y i represents the ground truth, andŷ i depicts the prediction results.

IV. EXPERIMENTAL RESULTS
This section presents the experimental results of the MLE prediction models. In the first part, the training data of the ML model is introduced. Next, the details of the experimental result are presented in the next subsection.

A. TRAINING DATA
To determine whether the SPAB system under various parameters is stable, this study used various rotor masses and bearing numbers as analysis parameters and used the MLE to construct the stable and unstable regions of the SPAB system (Fig. 10). The red and light green blocks in the figure belong to the unstable region, that is, chaotic behavior, and the blue or light blue blocks belong to the stable region. The unstable regions of the system are mostly concentrated in the larger rotor mass region. When the number of bearings is low (red area, bottom left of Fig. 10), the system exhibited nonperiodic chaotic motion at lower rotor mass. This result can effectively identify the stability of the SPAB system during operation and can be applied as a crucial reference basis for subsequent bearing system design. The obtained MLE map is depicted in Fig. 10.

B. EXPERIMENTAL RESULTS
The obtained MLE map is illustrated in Fig. 10. The MLP, SVM, DT, and RF were adopted to model the MLE map. In aspect of the MLP regressor, the size of the hidden layer  is 10, and the Rectified Linear Unit (ReLU) activation function is selected for the hidden neurons. Besides, the maximum iteration number is set to 500 epochs, and the Adam optimizer is also used for its training process. For the support vector regression model, Radial Basis Function (RBF) is applied for its kernel, and the regularization parameter is set to 1.0. Moreover, for DT and RF Regressor, the maximum depth of the tree is set to 20. These parameters are tuned by heuristic method. In this experiment, only a small part of the total samples was used in the training process. In this experiment, n/2, n/3, n/4, n/5, n/6, n/7, and n/8 of total samples were adopted as the training data of the proposed ML model where n denotes the total sample points in Fig. 10. Furthermore, Figs. 11-22 depict the prediction results by several ML model by using n/2, n/4, and n/8 of the total samples for the training data.
The final results revealed that RF achieved the best prediction performance. The DT prediction also achieved an acceptable result. By using this proposed method, the MLE value can be estimated swiftly in real time. However, although MLP and SVM exhibit excellent results in other research topics, these two algorithms achieved no advantage in the MLE prediction problem of the SPAB system. As depicted in Figs    the prediction result of DT and RF are high. Furthermore, a small quantity of training data may cause higher prediction errors. Nevertheless, as depicted in Fig. 22, the prediction error of RF with only 1/8 of the total training samples was small. DT provides slightly worse performance in the same training condition. This experimental result indicates that the ensemble model type (RF like) is suitable for addressing the MLE prediction problem. Tables I and II present the RMSE and MAE results, which revealed that the RF exhibited the best RMSE value and the DT exhibited the best performance in terms of MAE. MLP and SVM achieved unacceptable prediction results in this experiment. Therefore, the DT and RF models are suitable for addressing the MLE prediction problem. However, to compare the DT and RF, the standard deviation and the distribution range were calculated. Table 1 indicates that the RF provided the minimum average RMSE value. The    standard deviation and the distribution range of RF were smaller than that of the DT. Although DT achieved the minimum average MAE value (Table 2), the standard deviation and the distribution range of RF are still smaller than the DT value, which demonstrated that the overall performance of     RF is still better than that of DT. Therefore, the RF model gives high precision MLE prediction results and it is also preferred in this study. This method reduces the calculation time and keeps the prediction error in a certain acceptable range. The experimental result also indicates the feasibility and robustness of the proposed framework.
In this paper, the traditional iterative algorithm we used is hybrid method. For the MLE of each operational condition (ex: bearing number = 5.0 & rotor mass = 5.0 kg), it can be calculated and completed in 1 ∼ 2 minutes by PC with intel core i7-10510U CPU(1.8GHz). The detail calculation procedures are shown in Fig. 5. From Fig. 10, the mesh size of bearing number × rotor mass is 170 × 250, and we calculated each MLE by hybrid method, and totally needed about one month to obtain all MLEs in Fig. 10.
In this study, for the machine learning methods, the training process will be finished within just few minutes, and the prediction process will also be completed within few seconds. This result also demonstrates the merit of using the machine learning method for solving the MLE prediction problem. So, it is very efficient and precisely for the calculation of MLE by machine learning methods and also the results are the major contribution of this paper.

V. CONCLUSION
The calculation of MLE is crucial in the design of the SPAB system. However, the computation process is extremely time consuming and complex. To solve this problem, a RF-based MLE prediction model was proposed in this study. The experimental results revealed that the performance of the RF model was the best, and the RF model is most suitable for the prediction task. The prediction error of the RF model was within a certain acceptable range. The feasibility of the proposed framework was proved in this study. Thus, the stability VOLUME 8, 2020 estimation of the SPAB system can be accomplished rapidly by using the proposed model.