Analysis of the influence of piston–cylinder friction on the torsional vibration characteristics of compressor crankshaft system

With the development of compressors toward high speed and multiple columns, the torsional vibration phenomenon has become a major factor affecting the service life and reliability of the shaft system. Therefore, this paper considers the influence of friction between piston and cylinder on the instantaneous inertia of the crank connecting rod mechanism, establishes a nonlinear torsional vibration mechanics model of shale gas compressor shaft system, solves the natural frequency of the shaft system under undamped and damped conditions using the eigenvector method, and investigates the influence of the friction coefficient between piston and cylinder and the operating speed on the torsional vibration response of the shaft system under the self-excitation of the shaft system and the action of and excitation moment by the Runge–Kutta methods. The results show that after considering the friction between the piston and the cylinder, the second-order natural frequency of the shaft system shows a "high-low–high" fluctuation pattern. As the friction coefficient increases, the amplitude of the shaft system and the peak vibration speed show a rising trend. Meanwhile, when the speed increases, the vibration of the shaft system changes from chaotic \ period to the proposed periodic state, but the amplitude shows a decreasing trend. The research in this paper aims to improve the theory of nonlinear dynamics of compressor shaft systems, and the determined nonlinear parameters can be used to guide the operation and maintenance of compressors in engineering.


Introduction
Due to the process of shale gas development, the pressure decay is fast, resulting in variable compressor boosting conditions, coupled with the current development of the compressor in the direction of highspeed multiple columns [1], making the compressor shaft system torsional vibration phenomenon is increasingly prominent, which seriously affects the delivery efficiency and service life of the compressor. For this reason, many scholars have carried out research on the torsional vibration response law of the shaft system. Initially, the torsional vibration mechanics model only considered the influence of linear factors of the shaft system. For example, Feng et al. [2] solved the natural frequency of the shaft system based on the eigenvector method and avoided the torsional vibration of the shaft system by adding additional masses. Liu et al. [3] solved the torsional vibration response of the shaft system based on the MSM (Modal Superposition Method) and modified the structural parameters of the shaft system based on the SA (Sensitivity Analysis) method. Karimaei et al. [4] considered the effect of linear damping on the torsional vibration response of the system and proposed the selection of suitable dampers to reduce the torsional vibration amplitude of the shaft system. However, the above studies neglected the effects of variable inertia of the crank connecting rod, damping nonlinearity, and nonlinearity of friction between the piston and the cylinder on the torsional vibration of the shaft system when modeling the lumped mass of the shaft system, resulting in the inability to accurately predict the torsional vibration response of the shaft system [5,6].
For this reason, some scholars take into account the nonlinear factors of the shaft system when establishing the lumped mass model, such as Zhao [7] carried out the law of the influence of the nonlinear stiffness of the coupling on the torsional vibration response of the shaft system, pointing out that the larger the hard nonlinear stiffness, the larger the peak angular displacement. Babagiray et al. [8] carried out a study on the effect of system dry friction on crankshaft speed fluctuation. Karabulutd et al. [9] explored the effect of combined crank and flywheel moments of inertia on crankshaft speed fluctuation pattern. Zhang et al. [10] considered the distributed torsional flexibility of the crankshaft and pointed out that the flexibility of the crankshaft has a significant effect on the mechanical behavior of the system. Some scholars have also carried out research on the influence law of flywheel rotational inertia [10], main bearing structural parameters [11] and oil supply pressure [12] on the torsion angle of the shaft system.
Meanwhile, many scholars pointed out that neglecting the variable inertia characteristics of the crank-link mechanism would significantly reduce the calculation accuracy of the torsional vibration response of the shaft system [13,14]. Therefore, Pasricha [15], Wang [16] and Chen et al. [17] established a single-cylinder torsional vibration mechanics model considering the effect of variable inertia and discussed the stable torsional vibration states corresponding to different parameters. And Han [18] further improved the single-cylinder model by considering the effect of the connecting rod ratio on the variable inertia on this basis. However, there are still large errors in the nonlinear equations between the single-cylinder model and the actual multi-cylinder model. Therefore, Metallidis et al. [19] developed a multi-cylinder torsional vibration mechanics model based on a single-cylinder model and discussed the effects of shaft system stiffness and damping on velocity fluctuations. Xiang et al. [20] solved the natural frequency of the system through a multi-cylinder dynamics model to determine whether there is a risk of torsional vibration in the shaft system. Zhu et al. [21] then refined the multi-cylinder dynamics model by introducing factors such as variable rotational inertia and nonlinear friction between the piston and the cylinder. Guzzomi et al. [22,23] pointed out that the friction before the piston and cylinder causes the variable inertia characteristic of the shaft system, which leads to the change of the natural frequency of the shaft system, and experimentally verified this theory in a subsequent study, but Guzzomi did not carry out a follow-up study because the relationship between the two was too complicated.
At present, only the influence of the friction coefficient between piston and cylinder on the natural frequency of the shaft system has been explored, but the study of the response law of the shaft system to torsional vibration has not been carried out. Based on the shortcomings of the above study, this paper considers the influence of piston and cylinder friction on the characteristic of variable inertia of the shaft system and reasonably simplifies the relationship between them to establish a multi-cylinder torsional vibration mechanics model. The free vibration analysis of the shaft system under damped/undamped conditions was carried out, and the model results were compared with the finite element results to verify the reasonableness of the model. Finally, using the Runge-Kutta method, the variation law of shaft system torsional vibration response with friction coefficient and operating speed under the action of self-excitation and external torque is analyzed, which provides some theoretical reference to improve the establishment of shaft system multi-cylinder torsional vibration model. The rest of this paper is organized as follows. In Sect. 2, the kinetic equations of the crankshaft system are established, which take into account the effect of the friction coefficient between the piston and the cylinder on the crankshaft rotational inertia. In Sect. 3, the effect of the friction coefficient on the free vibration of the shaft system is explored. In Sect. 4, the effect of nonlinear factors on the vibration law when the shaft system is in two states of self-excitation and external excitation is explored. In Sect. 5, the findings of this paper and the implications for engineering guidance are summarized.
2 Dynamics modeling of crankshaft system

Torsional vibration concentrated mass of crankshaft system Modeling
In this paper, the crankshaft system of DTY500 shale gas compressor is taken as the research object, and its three-dimensional model is shown in Fig. 1. The crankshaft system dynamics model is established based on the lumped mass method, which mainly consists of crankshaft, coupling and motor and considers the piston-cylinder friction and internal damping between the shaft segments. The model has 15 degrees of freedom, as shown in Fig. 2, where I 0 -I 3 denotes the motor, I 4 -I 5 denotes the coupling, and I 6 -I 14 denotes the crankshaft. In Fig. 2, J i (i = 0, 1, …14) is the rotational inertia of each crankshaft segment. K j, j?1 (j = 0, 1, …13) is the torsional stiffness between the crankshaft segments. C m, m?1 (m = 0, 1, …13) is the internal damping factor. D n (n = 0, 1, …14) is the external damping factor. The damping factor can be determined by Eq. (1) [24].
where x is the angular frequency of crankshaft.

Derivation of single-cylinder instantaneous inertia expression
The crankshaft revolves along its axis, and the piston reciprocates motion through the connecting rod. Here, it is assumed that the piston center is co-linear with the crankshaft rotation center, so the single-cylinder crank connecting rod can be simplified to a crank-slider mechanism, as shown in Fig. 3. O, G, A, and B denote the rotation center of the crank, the mass center of the connecting rod, center of big end of connecting rod and the center of the crosshead pin, respectively. k denotes the ratio of AG to AB. r denotes the radius of the crank, and l denotes the length of the connecting rod.h,c, x denote the crank rotation angle, connecting rod swing angle and displacement of the piston, respectively. Without considering the effect of piston-cylinder friction, the instantaneous kinetic energy of a single crank-connecting rod mechanism is expressed as follows [25]: where m 1 is connecting rod large head mass. m 2 is connecting rod small head mass. m p is piston assembly mass. J c is crank rotational inertia. _ h is the rotation angular speed of the crank. _ x is the instantaneous velocity of the piston.
. Based on the above relationship, we obtain.
The total mass of the connecting rod remains the same before and after the transformation, and the position of the center of mass remains the same.
Based on the instantaneous kinetic energy equivalence method, the equivalent inertia J E h ð Þ of the single-cylinder crank-connecting rod mechanism is found as follows: Fig. 1 DTY500 compressor crankshaft system three-dimensional model Table 1 presents the structural parameters of the crank-connecting rod mechanism.
Guzzomi points out that the friction between the piston and the cylinder affects the effective inertia of the crank-connecting rod, so this paper introduces the G l ð Þ function, which is related to the speed of the piston and the direction of the lateral force F N , as shown in Eq. (5) [22]: However, the calculation time of F N is very long [23], which makes the modified equivalent inertia very complicated, so Guzzomi only carried out the effect of the friction coefficient on the natural frequency of the single-cylinder system. Therefore, in this paper, the G l ð Þ function is simplified so that the G l ð Þ function is not affected by _ x and F N , which is convenient for the subsequent establishment of the torsional vibration mechanics model of the crankshaft system.
MATLAB was used to plot the graph before and after the simplification of G l ð Þ, as shown in Fig. 4. The smaller the value of l, the smaller the difference between the G 3 l ð Þ ¼ cos hÀl 2 tan c sin h 1Àl 2 tan 2 c function and the original function, the piston ring material studied in this paper is PTFE and the cylinder liner material is cast iron, when under lubricated conditions l ¼ 0:35 $ 0:5, so the simplified G 3 l ð Þ function can describe the change law of the original function.
The equivalent inertia after correction is: where cos c ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi Since k ¼ 1=5 À 1=3, h h ð Þ = k cos h and ignoring the high-level multinomial, Eq. (7) is transformed as follows: Let kG 3 l ð Þ cos c þ 1 ¼ k Ã l, the effect of different values of k Ã on J E is explored, as shown in Fig. 5.
As can be seen from Fig. 5, when k* is taken as appropriate, the correlation coefficient between the simplified function and the original function can be made greater than 0.98, and the maximum error is less than 2.0%, and it can be analyzed that there is a nonlinear relationship between k* and l, that is, kÃ ¼ 14l 2 À 17:94l þ 7:542 ð Þ , so Eq. (8) can be simplified as follows: The above equation is simply given by: According to Ref. [24], under the condition that the friction between piston ring and cylinder liner is not considered, the traditional equivalent inertia of a single cylinder is given by: Taking cylinder #1 as an example, the average rotational inertia obtained using Eq. (12) is 0.0777 kg.m 2 . When l = 0.35, 0.4, 0.45, 0.5, the average rotational inertia of the cylinder obtained using Eq. (11) is 0.0785, 0.0785, 0.0783, 0.0783, respectively, and the maximum error of the calculated result with Eq. (12) is 0.8%. It is further verified that the simplification of G l ð Þ function and cylinder rotational inertia equation is reasonable. Lagrange equation corresponding to the four-cylinder crankshaft system is as follows: where the kinetic energy is and the potential energy is and the dissipative energy is where J i , K i , D i represent the inertia, torsional stiffness and damping of the mass./ i is the torsion angle. i is the mass number.
taking the compression cylinder as an example, the compressed gas produces an excitation torque, expressed as Q i ¼ M ri , and its dynamics equation is as follows: Compared with the constant inertia equation, _ J E h i ð Þ _ h 2 i represents the effect of introducing crank linkage variable inertia on the torsional vibration response of the shaft system [25].
Introducing the parameter s¼x t, and for a small torsion angle /, we obtain cos / % 1; sin / % /. Substituting Eq. (17), the single-cylinder vibration differential equation can be obtained as follows: In the literature [16,17], the linearized approach to Eq. (18) does not accurately reflect the effect of variable inertia on the shaft system torsional vibration, so this paper retains the / as well as _ / i quadratic terms.
n i represents the phase difference between the cylinders.
Here, it is specified that n 7 ¼ n 8 ¼ 0,n 12 ¼ n 13 ¼ p/2. The non-dimensional conversions are expressed as follows: where

Analysis of crankshaft system excitation torque
The externally excited torque of the crankshaft system consists of the torque M c generated by the compressed gas, the torque M j generated by the reciprocating inertia force, and the internally excited torque M i caused by the non-constant mass. The expression for the internal excitation torque is given by: The external excitation moment is expressed as follows: where M r is the external excitation torque. D is the piston diameter; d is the piston rod diameter. P C and P H are the cover-side pressure and shaft-side pressure, respectively, which can be obtained by professional thermodynamic software. m j is the reciprocating motion mass.
Applying the Fourier series expansion to M r yields the following equation: where M 0 is the average and harmonic component of the excitation torque. n denote the harmonic orders. Taking cylinder 1 as an example, the Fourier series expansion of M r is carried out at different harmonic orders, here n = 5,10,15, and the errors corresponding to different harmonic orders are compared, as shown in Fig. 6. From Fig. 6, it can be seen that there is a certain deviation between the fitted function and the actual excitation moment M r when n = 5, while the RMSE (Root mean square error) value of the system is less than 25 for n = 10 or 15, which indicates that the fitted function can better reflect the change trend of M r . Meanwhile, since the contribution of higher-order harmonics to the torsional amplitude of the shaft system is not significant, the unfolding harmonic order n is set to 10 in this paper.

Shaft system with undamped natural frequency solution
We assume that the damping coefficient has little effect on the natural frequency of the crankshaft system, and thus, the mathematical model of free vibration is simplified to The eigenvalue method x 2 J À K ¼ 0 is used to solve the natural frequency of the crankshaft system, and Ref [26] points out that the lumped mass method can only accurately solve the anterior 1/3 order mode of the corresponding model, so the first 2 order modes of the crankshaft system are extracted in this paper. Table 2 presents the dynamic parameters of the lumped mass model, where the inertia parameters are provided by the manufacturer, the stiffness Fig. 6 The fitting function corresponding to different harmonic orders parameters are obtained by the finite element software, and the damping parameters are derived from Eq. (1). As can be seen from Fig. 7, the first-and secondorder natural frequencies of the shaft system show a sinusoidal trend. For the same order natural frequency, the value of l does not affect its change law, but the second-order natural frequency will have a ''highlow-high'' fluctuation region near 180°. Reference [27] shows that when l ¼ 0, the instantaneous inertia change of the crank-link mechanism only affects the shaft system fluctuation amplitude, and there will not be a secondary fluctuation region of the natural frequency. Therefore, after considering the influence of friction between the piston ring and cylinder liner, it will lead to the complicated variation law of the highorder natural frequency of the shaft system, and the torsional vibration behavior of the shaft system cannot be accurately calculated.

Shaft system with damped natural frequency solution
Since the shaft system damping matrix has small nondiagonal elements in canonical coordinates, the shaft system modal damping ratio is calculated as follows: where cp i is the damping matrix in canonical coordinates. x i is the angular frequency. Using the relationship between the undamped and damped natural frequencies, i.e., , the values of the damped natural frequencies of the shaft system are obtained, and Table 3 is used to show the natural frequencies of the shaft system for the three models.  Fig. 7 The variation law of the first two orders of the crankshaft system natural frequency with angle From Table 3, it can be seen that the l value has less influence on the average value of the first-and secondorder natural frequencies of the shaft system, and the maximum error between the calculation results of the lumped mass model and the finite element model is 10.5% after considering the influence of damping. Since the lumped mass takes into account the influence of variable inertia of the crank-link mechanism, the obtained results are smaller than the finite element model, but the maximum error is still in an acceptable range, which proves that the parameter values in Table 2 are reasonable.

Forced vibration analysis of crankshaft system
The theoretical basis indicates that the torsional behavior of the shaft system is influenced by the piston-cylinder friction coefficient l, the crankshaft angular velocity x, and the excitation moment M. In this paper, the numerical results are used to discuss the effects of these three parameters on the torsional behavior of the shaft system in order to better control the amplitude of the shaft system torsional vibration. The numerical results include the torsional response of the shaft system in the time domains, the phase plane, and the Poincare surface.

Effect of piston-cylinder friction coefficient on torsional vibration
The study of the free vibration of the shaft system shows that with the increase of the piston-cylinder friction coefficient l, the variation law of the higherorder natural frequency of the shaft system is more complicated, so this paper firstly explores the variation law of the torsional vibration response of the shaft system with the value of l under the self-excitation state, as shown in Table 4 and Fig. 8. From Fig. 8, it can be seen that: (1) When the value of l increases, the time domain waveform of the shaft system always shows an exponentially decreasing trend, which indicates that the vibration of the shaft system is stable and there is no jump phenomenon. (2) When the value of l increases, partial time domain waveforms show a similar pattern of change, but the maximum amplitude increases slightly, indicating that the value of torsional amplitude of the shaft system can be reduced by improving the lubrication between the piston ring and the cylinder liner. (3) When the value of l increases, the phase diagram and the Poincare surface of section are basically similar, the trajectories are multiple circular superimposed but do not overlap, but as the value of l increases, the phase trajectory will show a right shift phenomenon, which indicates that the vibration of the shaft system is the proposed periodic motion when the value of l changes in a small range, the amplitude increased by about 4.3%.

Effect of rotational speed change on torsional vibration
When the friction coefficient of piston ring and cylinder liner is constant, the torsional vibration response of the shaft system is analyzed by changing the operating speed of the shaft system. This paper focuses on the effect of four commonly used operating speeds (500r/min, 1000r/min, 1500r/min and 2000r/ min) on the torsional vibration response of the shaft system. Figure 9 shows the variation pattern of the self-excitation response of the shaft system at different speeds: For the time domain waveforms corresponding to different rotational speeds, we can roughly divide them into three regions, i.e., region I: Fluctuation state; region II: Convergence state; and region (1) When the speed of the shaft system gradually increases, the time when the shaft system is in the fluctuation region will increase. This is because the lower the speed, the greater the internal damping of the shaft system, which makes the shaft system dissipate more energy, so it can transition to the stable stage faster, but when the shaft system is in the stable stage, the lower the speed, the greater the vibration amplitude.
(2) For the convergence stage, when the shaft system speed n = 500r/min, the phase diagram  Fig. 8 The variation law of torsional vibration response of shaft system corresponding to different l values of the shaft system shows up-and-down fluctuation and does not form a circular orbit, which indicates that the motion of the shaft system is in a chaotic state. When the shaft system speed n = 1000r/min, the phase diagram of the shaft system is a ''jagged'' circular orbit, which indicates that the motion of the shaft system is in an unstable cyclic state. When the shaft system speed n = 1500r/min and 2000r/min, the phase diagram of the shaft system is a superimposed and non-repeating circle, which indicates that the motion of the shaft system is in a cyclic state.
In summary, try to ensure that the speed of the shaft system is above 1500r/min, which can not only ensure the reduction of the vibration amplitude of the shaft system in stable operation, but also the vibration state of the shaft system in all three stages can be predicted, which is conducive to the control of shaft system torsional vibration.

Torsional vibration response analysis of shaft system for excitation moment
The effect of different friction coefficients and rotational speeds on the torsional vibration response of the shaft Fig. 9 The variation law of torsional vibration response of shaft system corresponding to different n values system under self-excitation is discussed above, while this section focuses on the torsional vibration response of the shaft system under the action of the excitation moment. Meanwhile, in order to discuss the influence of various parameter interactions on the shaft system torsional vibration, four parameter combinations (l ¼ 0:35; n ¼ 500r= min,l ¼ 0:5; n ¼ 500r= minn ¼ 2000r= min,l ¼ 0:5; n ¼ 2000r= min) were chosen as the boundary conditions for the numerical analysis, and the variation law of the shaft system torsional vibration response was obtained, as shown in Fig. 10. From Fig. 10, it can be seen that: (1) The vibration of the shaft system never appears chaotic under the action of the excitation moment. However, as the speed of the shaft system increases, the vibration state of the shaft system changes from a periodic motion to a proposed periodic motion, indicating that the higher the operating speed, the closer the vibration of the shaft system is to the chaotic state. Meanwhile, when the working speed is 500 r/min, the maximum vibration amplitude is 2.0°, and when the working speed is 2000 r/min, the maximum amplitude is 1.06°, which is reduced by 47%, so in the actual project, ensuring that the compressor working speed is medium and high speed. (2) When the shaft system is in low-speed motion (500r/min), the phase diagram is a closed circle, indicating that the shaft system vibration is periodic motion. However, the circular trajectory is not smooth and shows a trend of ''sawtooth decay,'' and the sawtooth density decreases with the increase in friction coefficient l. This indicates that small changes in the initial conditions can lead to changes in the vibration state of the shaft system. (3) When the shaft system is in medium and highspeed motion (2000 r/min), the phase diagram is a superposition of multiple closed circles and does not overlap, indicating that the shaft system vibration is a proposed periodic motion. With the increase in friction coefficient l, the phase diagram appeared upward phenomenon, which shows that in the high-speed operation, the greater the friction coefficient shaft system speed fluctuations, so should ensure that the piston ring and cylinder liner lubrication conditions between good.

Conclusion
(1) In this paper, the equations for calculating the rotational inertia of the crank-link mechanism are simplified, and the average error in calculating the rotational inertia between the simplified equations and the actual measurement is 1.4%. Therefore, the simplified equation can reasonably reflect the variation law of rotational inertia of the crank linkage mechanism with angle and facilitate the establishment of the subsequent torsional vibration mechanics model. (2) After considering the effect of the friction coefficient between the piston ring and the cylinder block, a ''high-low-high'' fluctuation in the second-order natural frequency of the shaft system was observed under undamped conditions, indicating that the introduction of the friction coefficient may lead to inaccurate prediction of the torsional vibration behavior of the shaft system. Under the condition of damping, the error of the natural frequency obtained from the torsional vibration mechanics model and the finite element model is 10.5%, which verifies the reasonableness of the parameters of the torsional vibration mechanics model. (3) In the self-excited state of the shaft system, the vibration state of the shaft system does not change with the change of the friction coefficient, all are the state of the proposed periodic motion, but the larger the friction coefficient, the vibration amplitude of the shaft system will show a tendency to increase. At the same time, as the working speed increases, the vibration state of the shaft system transforms from chaotic state to the proposed periodic state, and the amplitude gradually decreases. (4) Under the action of excitation torque, with the increase in working speed, the vibration state of the shaft system is transformed from periodic motion to the proposed periodic state, but the vibration amplitude gradually decreases. When the working speed is the same, with the increase in friction coefficient, the vibration speed of the shaft system shows an increasing trend. Therefore, in the actual project should ensure that the compressor is in high-speed operation and the piston ring with the cylinder liner between the lubrication state is good. (5) The stiffness and damping in the torsional vibration mechanics model of the shaft system established in this paper are linearly varying, while in actual engineering, the stiffness and damping are mostly nonlinear, especially the stiffness of the coupling will show the characteristics of the hard or soft stiffness [28]. Therefore, the effect of hard or soft stiffness of the coupling on the torsional vibration behavior of the shaft system should be considered in subsequent studies.