Quasi-static and dynamic characteristics analysis of a marine gear system under coupling misalignment

Diaphragm couplings are critical in marine gear trains but inevitably suffer from both parallel and angular misalignment, which degrade tooth-mesh conditions and excite low-frequency vibrations. This paper presents a quasi-static contact analysis model for gear systems that accounts for coupling misalignment using a coupled iterative calculation framework that integrates the misalignment mechanics model, the sliced static model, and the loaded tooth contact analysis (LTCA) model. The study analyzes the impact of coupling misalignment on the contact characteristics of the gear pair under non-ideal mesh conditions and establishes the differential equation of motion for the marine gear system. Findings reveal that coupling misalignment significantly affects mesh stiffness, which fluctuates at the misalignment frequency (twice the shaft frequency). As the imposed misalignment increases, the differences in meshing misalignment and load distribution between the left- and right-hand teeth of the double-helical gear pair increase progressively. At light loads, the quasi-static contact characteristics are more susceptible to coupling misalignment, resulting in regions of contact loss on the tooth surfaces. Coupling-misalignment excitation also induces system resonance. As misalignment and load torque increase, system vibration intensifies, with load torque primarily influencing vibration responses to mesh frequency excitation. These insights deepen the understanding of marine transmission systems and inform strategies for addressing misalignment-induced vibrations, which can enhance the performance, reliability, and longevity of marine gear systems in practical applications.

1 Introduction

Available failure data associate rotor misalignment with over 60 % of faults in rotating equipment (Moghe and Mahure, 2017). Tests have shown that misalignment increases the vibration amplitude of the gearbox by 20 % to 40 %, leading to bearing damage (Kumar and Hirani, 2021). Diaphragm couplings are widely used in marine gear systems because they can accommodate a certain degree of misalignment. However, due to deviations in base height, coupling assembly errors, hull deformation, and other factors, couplings inevitably experience parallel and angular misalignment. Forces and moments induced by coupling misalignment compromise gear system stability. These forces and moments cause shaft deformation, gear pair disengagement, system vibration, and bearing wear. Accordingly, the effects of coupling misalignment on the quasi-static contact behavior and dynamic characteristics of marine gear systems should be systematically examined.

Most marine gear systems comprise multiple series-connected gearboxes, resulting in multiple sources of mesh excitation and coupling misalignment between the gearboxes. Dynamic interactions among the gearboxes induce system-wide coupled vibration. Changes in transmission chain flexibility under different operating conditions induce variations in mesh excitation, thereby giving rise to highly complex vibration characteristics of the system. Gong et al. (2021) proposed a hierarchical dynamic modeling approach based on a whole-system-first, substructure-second strategy. Meng et al. (2014) established a transient dynamic model of the mode-switching process for the fuel-combined-power cross-connect gearbox (CCG) using ADAMS and Simulink. Chen et al. (2022) developed a computational method for determining the clutch impact load during rapid CCG paralleling. Research to date has centered on the enhancement of operational smoothness and energy efficiency and the formulation of control schemes for transitions between operating modes.

Coupling misalignment is classified as parallel, angular, or combined misalignment. Coupling misalignment produces angular deflection and deformation. As the shaft system rotates, coupling misalignment generates time-varying forces. Coupling-misalignment modeling approaches include the mechanical method (Wang and Jiang, 2018), the equivalent shaft segment method (Hujare and Karnik, 2018), and the Lagrangian energy method (Al-Hussain and Redmond, 2002). The mechanical method uses the coupling geometry to establish the force–deformation relationship and to derive the forces induced by coupling misalignment. Gibbons (1976) derived expressions for coupling forces and moments induced by parallel misalignment under static loading conditions. Sekhar and Prabhu (1995) subsequently refined Gibbons' model by incorporating harmonic terms representing the forces and moments induced by coupling misalignment into the dynamic formulation. Xu and Marangoni (1994) derived an expression for the moment induced by angular misalignment using the kinematic theory of a Hooke joint. Xia et al. (2019) derived expressions for the forces and moments induced by parallel and angular misalignment of hexagonal flexible couplings based on deformation at the bolt locations. In the equivalent shaft segment method, also termed the finite-element method, the coupling is replaced by a shaft segment of equivalent bending stiffness. Lee and Lee (1999) and Pennacchi et al. (2012) modeled coupling misalignments in axially misaligned beam elements, with parallel and angular misalignments expressed in terms of rotational displacement vectors at the coupling nodes. The Lagrangian energy method uses Lagrange's equations and coordinate transformations to establish the system energy equation and differential equation of motion directly. This method avoids the complicated force analysis of the coupling. Al-Hussain (2003) derived an expression for the potential energy associated with angular misalignment in the flexible coupling. Zhen et al. (2020) established a dynamic model for parallel and angular misalignment subject to rotor integrity constraints. Li et al. (2017) considered the displacement constraints of adjacent rotors and the geometric relationships governing their motion to derive the kinetic and potential energy equations of a misaligned coupling-rotor system.

Research on coupling misalignment has mainly focused on rotor-bearing-coupling systems. The stability (Li and He, 2010; Huang et al., 2011), nonlinear vibration characteristics (Wu et al., 2020; Fu et al., 2022), and friction problems (Jin et al., 2021) of rotor systems under misalignment of couplings have been extensively studied in rotor dynamics. Patel and Darpe (2009a, b), Prabhakar et al. (2002), and Lu et al. (2021) found, by theoretical and experimental methods, that coupling misalignment causes a strong component at twice the shaft frequency in the rotor system vibration response. Lee and Lee (1999) confirmed through simulations and experiments that the shaft center orbit becomes elliptical under coupling misalignment, which is a typical characteristic of coupling misalignment. Wang et al. (2022b) indicated that misalignment raises the critical speed of the system and intensifies its axial vibration. Wu et al. (2020) found through simulation that the dynamic misalignment of the coupling is affected by external excitations. When other external excitations weaken the dynamic misalignment force/moment, the rotor system vibration response is linear. Wan et al. (2012) and Zhen et al. (2020) examined the stability of a multi-span rotor system supported by sliding bearings under diaphragm-coupling misalignment. Their results indicated that coupling misalignment progressively destabilizes the sliding bearings as the rotor speed increases. A limited number of studies have examined coupling misalignment in gear systems. Kumar and Hirani (2021) conducted experiments to study the impact of misalignment on gearboxes and found that misalignment can cause severe damage to bearings and uneven loading on the tooth surfaces. Wang et al. (2022a) and Xu et al. (2015) incorporated coupling-misalignment excitation into gear system dynamic models. They studied the vibration response under misalignment faults but did not account for the effects of forces and moments induced by misalignment on the meshing excitation of the gear pair. Coupling misalignment causes shaft misalignment and gear mesh misalignment. Yang et al. (2022) studied the changes in gear mesh stiffness under shaft misalignment based on an improved loaded tooth contact analysis (LTCA) model. Wang et al. (2021) proposed a gear pair mesh stiffness model incorporating inter-slice coupling and examined the influence of mesh misalignment on the resulting stiffness.

Further investigations have examined the effects of misalignment on bearing support characteristics, nonlinear contact behavior, load distribution, and more complex rotor dynamic responses. Bearing and coupling misalignment can alter vibration amplitudes and frequency components, system stability, shaft center trajectories, nonlinear contact behavior, and load distribution (Ma et al., 2015; Xu et al., 2022, 2024; Guan et al., 2025; Zhou et al., 2025; Nan et al., 2024; Ma et al., 2025; Jiang et al., 2025). These studies have provided valuable insights into the dynamic mechanisms associated with misalignment excitation and its interaction with bearing support, rub impact, and coupling contact. In gear transmission systems, installation errors have also been found to interact with structural faults and further amplify their effects on the system's dynamic response (Huang et al., 2025). However, most existing studies have focused primarily on the overall vibration response of rotor-bearing or coupling systems, while the mechanism by which coupling-misalignment excitation affects gear mesh misalignment, tooth surface contact states, and the system's dynamic response has received comparatively limited attention. In particular, the coupled effects of coupling misalignment on the left- and right-hand tooth surfaces of herringbone gears and the resulting redistribution of mesh misalignment and contact loads remain insufficiently investigated.

To address this limitation, the influence of coupling misalignment on both the system vibration response and the actual meshing state of the gear pair should be considered simultaneously. For marine gear systems, the contact characteristics of the gear pair under coupling misalignment can be analyzed, and the resulting gear pair meshing excitation, together with the excitation induced by coupling misalignment, can then be incorporated into the dynamic model to improve the accuracy of system vibration–response predictions. The method proposed in this study is primarily applicable to multistage parallel gear transmission systems. The modeling procedure is relatively straightforward and is suitable for various types of cylindrical gears, including spur, helical, and herringbone gears.

This paper investigates the transmission chain of a marine multi-engine, multi-shaft gear system. A mechanical model for parallel and angular misalignment of the diaphragm coupling is established based on geometric deformation relationships. A sliced static model incorporating misalignment forces and moments is coupled with the loaded tooth contact analysis method to analyze the quasi-static contact characteristics of gear pairs under coupling misalignment. The dynamic characteristics of the system under mixed excitation are investigated by considering coupling-misalignment excitation and gear mesh excitation simultaneously. The effects of coupling misalignment and load torque on the quasi-static and dynamic characteristics of the gear system are analyzed.

2 Research objects and parameters

Figure 1 shows the structural schematic of the marine multi-engine multi-shaft gear system, in which the engines drive the propeller shafts individually or jointly through clutch engagement and disengagement. To facilitate the study, the transmission chain comprising the M1 engine, clutch gearbox, and parallel gearbox is selected for analysis. The effect of coupling misalignment between the two gearboxes on the gear mesh state and system vibration characteristics is analyzed. The gear parameters for each gearbox are listed in Table 1. Each pair of double-helical gears has a face width on each side greater than 10 times the modulus and is therefore classified as a wide-faced gear pair. Tables 2 to 4 show the shaft parameters, bearing support stiffness, and diaphragm coupling stiffness, respectively. The gear shafts are supported by sliding bearings arranged symmetrically about the gear pair. The stiffness and damping coefficients of the journal bearings are calculated based on classical hydrodynamic lubrication theory. The theoretical formulations of the direct and cross-coupled stiffness and damping coefficients have been systematically derived by Gong et al. (2025). Therefore, only the resulting bearing dynamic coefficients are employed in the present model, and the general derivation is not repeated here.

Figure 1 Configuration of the research object.

Table 1 Gear parameters.

Table 2 Shaft parameters.

Table 3 Bearing stiffness.

Table 4 Diaphragm coupling stiffness.

3 Quasi-static contact analysis model for a gear system considering coupling misalignment

3.1 Mechanical model of coupling misalignment

Following the method presented by Xia et al. (2019), a misalignment model of the diaphragm coupling is established, with the left and right coupling halves used to define the coordinate system for parallel misalignment, as shown in Fig. 2. Coupling misalignment is defined as the relative offset between the left and right coupling halves. The right half of the coupling is used as the reference in the coordinate system, and it is assumed that the left half of the coupling undergoes parallel misalignment. Trigonometric functions are used to calculate the deformation at bolt holes A, B, and C in the left coupling half along each direction. These deformations are multiplied by the corresponding directional diaphragm stiffnesses to obtain the coupling force induced by parallel misalignment. The centers of the left- and right-hand coupling sections are located at points O and O^′, respectively, and the angle between the y axis and the OO^′ line is β. For the left coupling half, ϵ_x and ϵ_y represent the misalignments along the x and y axes, respectively. R denotes the distance from the rotor axis to the center of the bolt hole, and θ represents the angular position of bolt hole A about the coupling center O^′. The relationship between β and the parallel offsets ϵ_x and ϵ_y is $\tan \mathit{\beta }={\mathit{ϵ}}_x/{\mathit{ϵ}}_y$, and the eccentricity is OO${}^{\prime }=\sqrt{\left({\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}\right)}$.

Figure 2 Schematic of parallel coupling misalignment.

According to the law of cosines, the deformations at bolt holes A, B, and C can be obtained as

$$\begin{array}{}\text{(1)}& {\displaystyle }\begin{array}{rl}{\mathit{\mu }}_{\mathrm{1}}& =\text{OA}-{\text{O}}^{\prime }\text{A}\\ & =\sqrt{(R^{\mathrm{2}}+{\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}+\mathrm{2}R\sqrt{\left({\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}\right)}\sin (\mathit{\beta }+\mathit{\theta }\left)\right)}-R,\end{array}\text{(2)}& {\displaystyle }\begin{array}{rl}{\mathit{\mu }}_{\mathrm{2}}& =\text{OB}-{\text{O}}^{\prime }\text{B}\\ & =\sqrt{(R^{\mathrm{2}}+{\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}+\mathrm{2}R\sqrt{\left({\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}\right)}\cos (\mathit{\pi }/\mathrm{6}+\mathit{\beta }+\mathit{\theta }\left)\right)}-R,\end{array}\text{(3)}& {\displaystyle }\begin{array}{rl}{\mathit{\mu }}_{\mathrm{3}}& ={\text{O}}^{\prime }\text{C}-\text{OC}\\ & =R-\sqrt{(R^{\mathrm{2}}+{\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}-\mathrm{2}R\sqrt{\left({\mathit{ϵ}}_x^{\mathrm{2}}+{\mathit{ϵ}}_y^{\mathrm{2}}\right)}\sin (\mathit{\pi }/\mathrm{3}+\mathit{\beta }+\mathit{\theta }\left)\right)}.\end{array}\end{array}$$

The deformation components along the x and y axes are

$$\begin{array}{}\text{(4)}& {\displaystyle }{\mathit{\mu }}_x={\mathit{\mu }}_{\mathrm{1}}\cos \left(\mathit{\theta }\right)+{\mathit{\mu }}_{\mathrm{2}}\cos (\mathit{\theta }+\mathrm{2}\mathit{\pi }/\mathrm{3})+{\mathit{\mu }}_{\mathrm{3}}\cos (\mathit{\theta }+\mathrm{4}\mathit{\pi }/\mathrm{3}),\text{(5)}& {\displaystyle }{\mathit{\mu }}_y={\mathit{\mu }}_{\mathrm{1}}\sin \left(\mathit{\theta }\right)+{\mathit{\mu }}_{\mathrm{2}}\sin (\mathit{\theta }+\mathrm{2}\mathit{\pi }/\mathrm{3})+{\mathit{\mu }}_{\mathrm{3}}\sin (\mathit{\theta }+\mathrm{4}\mathit{\pi }/\mathrm{3}).\end{array}$$

The x- and y-direction components of the coupling force induced by misalignment are then given by Eqs. (6) and (7), respectively, where k_D is the diaphragm coupling stiffness, and μ_x and μ_y are the deformation components along the x and y axes.

$$\begin{array}{}\text{(6)}& {\displaystyle }{F}_{px}=k_{\mathrm{D}}{\mathit{\mu }}_x\text{(7)}& {\displaystyle }{F}_{py}=k_{\mathrm{D}}{\mathit{\mu }}_y\end{array}$$

Parallel coupling misalignment is a planar mechanical problem, whereas angular coupling misalignment is a three-dimensional mechanical problem. Angular misalignment generates forces in the xy-plane, an axial force along the shaft axis, and bending moments about the x and y axes. As illustrated in the angular-misalignment schematic in Fig. 3, the left coupling half is fastened to the diaphragm by bolts passing through holes A, B, and C. γ is the offset angle. When one coupling half is angularly offset, the centers of the left and right coupling halves are separated by an axial distance ϵ_z, and the planes of the two halves intersect along a line l. The distance ΔR between this intersection line and the coupling center is calculated from ϵ_z and the offset angle γ using Eq. (8).

$$\begin{array}{}\text{(8)}& \mathrm{\Delta }R={\displaystyle \frac{{\mathit{ϵ}}_z}{\tan \mathit{\gamma }}}\end{array}$$

The deformations of the diaphragm at bolt holes A, B, and C are

$$\begin{array}{}\text{(9)}& {\displaystyle }\mathrm{AA}{}^{\prime }=\left(\mathrm{\Delta }R+R\sin \mathit{\theta }\right)\sqrt{\left(\mathrm{2}-\mathrm{2}\cos \mathit{\gamma }\right)},\text{(10)}& {\displaystyle }\mathrm{BB}{}^{\prime }=\left[\mathrm{\Delta }R+R\sin (\mathit{\theta }+{\displaystyle \frac{\mathrm{2}\mathit{\pi }}{\mathrm{3}}})\right]\sqrt{\left(\mathrm{2}-\mathrm{2}\cos \mathit{\gamma }\right)},\text{(11)}& {\displaystyle }\mathrm{CC}{}^{\prime }=\left[\mathrm{\Delta }R+R\sin (\mathit{\theta }+{\displaystyle \frac{\mathrm{4}\mathit{\pi }}{\mathrm{3}}})\right]\sqrt{\left(\mathrm{2}-\mathrm{2}\cos \mathit{\gamma }\right)}.\end{array}$$

The misalignment forces at each bolt hole are given by $F_{\mathrm{A}}=k_{\mathrm{s}}{\mathrm{AA}}^{\prime }$, $F_{\mathrm{B}}=k_{\mathrm{s}}{\mathrm{BB}}^{\prime }$, and $F_{\mathrm{C}}=k_{\mathrm{s}}{\mathrm{CC}}^{\prime }$, where k_s is the diaphragm bending stiffness. The misalignment forces at each bolt hole are resolved along the axes of the three-dimensional coordinate system to obtain the components F_Ax, F_Bx, F_Cx, F_Ay, F_By, F_Cy F_Az, F_Bz, and F_Cz. Summing the forces along each axis gives the resultant force induced by angular misalignment, as shown in Eq. (12). Since the torque T on the left coupling half acts along its axis, it must be resolved into components for the angularly offset right coupling half. The x- and y-direction moment components, T_x and T_y, are given in Eq. (13).

$$\begin{array}{}\text{(12)}& {\displaystyle }\left\{\begin{array}{l}{F}_x=F_{\mathrm{A}x}+F_{\mathrm{B}x}+F_{\mathrm{C}x}\\ F_y=F_{\mathrm{A}y}+F_{\mathrm{B}y}+F_{\mathrm{C}y}\\ F_z=F_{\mathrm{A}z}+F_{\mathrm{B}z}+F_{\mathrm{C}z}\end{array}\right.\text{(13)}& {\displaystyle }\left\{\begin{array}{l}{T}_x=T\sin \mathit{\gamma }\cos \mathit{\theta }\\ T_y=T\sin \mathit{\gamma }\sin \mathit{\theta }\end{array}\right.\end{array}$$

Figure 3 Schematic of angular coupling misalignment.

3.2 Sliced static model and mesh misalignment

Within the framework of generalized finite-element theory, the system is represented by sliced mesh, bearing, shaft, and coupling elements, thereby forming the sliced static model illustrated in Fig. 4. At each meshing position, the mesh stiffness K_m is obtained by iteratively solving the static equilibrium equations and the loaded tooth contact analysis model, with mesh misalignment calculated from the gear nodal displacements. The sliced coupling stiffness K_slice is simulated using Timoshenko beam stiffness. The method proposed by Gong et al. (2022b) is utilized to obtain the bearing element stiffness matrix K_b and the shaft and coupling element stiffness matrices K_s and K_c, respectively, followed by assembly of the element matrices. Section 3.1 shows that the coupling-misalignment forces and moments are equal in magnitude and opposite in direction for the left and right coupling halves. Therefore, the misalignment force/moment is applied to the respective nodes on either side of the coupling element in the static formulation. The equilibrium relation governing the sliced model is formulated in matrix form in Eq. (14), where P₀ is the load column vector. The right-hand-side term ${\mathit{F}}_{\mathrm{c}}({\mathit{ϵ}}_x,{\mathit{ϵ}}_y,{\mathit{ϵ}}_z,\mathit{\gamma },\mathit{\theta }$) represents the misalignment-induced force/moment column vector of the coupling, which is related to the offset, offset angle, and rotation angle. The explicit expression for the misalignment force/moment column vector is shown in Eq. (15). For each paired nodal degree of freedom, the two coupling nodes carry oppositely directed forces and moments.

$$\begin{array}{}\text{(14)}& {\displaystyle }\begin{array}{rl}& \left[{\mathbf{K}}_{\mathrm{m}}\left(\mathit{\theta }\right)+{\mathbf{K}}_{\mathrm{slice}}+{\mathbf{K}}_{\mathrm{b}}+{\mathbf{K}}_{\mathrm{s}}+{\mathbf{K}}_{\mathrm{c}}\right]q\left(\mathit{\theta }\right)\\ & ={\mathit{P}}_{\mathrm{0}}+{\mathit{F}}_{\mathrm{c}}\left({\mathit{ϵ}}_x,{\mathit{ϵ}}_y,{\mathit{ϵ}}_z,\mathit{\gamma },\mathit{\theta }\right)\end{array}\text{(15)}& {\displaystyle }\begin{array}{rl}& {\mathit{F}}_{\mathrm{C}}({\mathit{ϵ}}_x,{\mathit{ϵ}}_y,{\mathit{ϵ}}_z,\mathit{\gamma },\mathit{\theta })=\left\{\mathrm{0},\mathrm{\dots },F_x,F_y,F_z,T_x,T_y,\mathrm{0},\mathrm{\dots },\right.\\ & \left.-F_x,-F_y,-F_z,-T_x,-T_y,\mathrm{0},\mathrm{\dots },\mathrm{0}\right\}{}^{\prime }\end{array}\end{array}$$

Each gear sub-slice is represented by a two-node element with 12 degrees of freedom. The generalized coordinates of the jth slice are defined as $q_{\mathrm{m}}^j=\left\{x_{\mathrm{p}}^j,y_{\mathrm{p}}^j,z_{\mathrm{p}}^j,{\mathit{\theta }}_{x\mathrm{p}}^j,{\mathit{\theta }}_{y\mathrm{p}}^j,{\mathit{\theta }}_{z\mathrm{p}}^j,x_{\mathrm{g}}^j,y_{\mathrm{g}}^j,z_{\mathrm{g}}^j,{\mathit{\theta }}_{x\mathrm{g}}^j,{\mathit{\theta }}_{y\mathrm{g}}^j,{\mathit{\theta }}_{z\mathrm{g}}^j\right\}$. As shown in Fig. 5, O_p and O_g are the rotation centers of the driving and driven gears, β_b is the base helix angle, $\mathit{\gamma }=\mathit{\alpha }\mp {\mathit{\theta }}_{\mathrm{p}}$ denotes the coordinate system projection parameters, the symbol “∓” indicates the direction of rotation of the gear, α is the pressure angle, θ_p is the installation phase angle, and r_p and r_g are the base circle radii of the driving and driven gears. The driving and driven gears are connected through the slice mesh stiffness $k_{\mathrm{m}}^j$, and j is the number of slices. The slices are connected through Timoshenko beam elements. The static equations of the jth gear pair slice are given in Eq. (16), where ${\mathit{\sigma }}_{\mathrm{m}}^j$ and $e_{\mathrm{m}}^j$ are the displacement and error along the line of action of the jth sliced mesh element. $k_{\mathrm{p}\mathit{\_}i}(i=x,y,z,{\mathit{\theta }}_{\mathrm{x}},{\mathit{\theta }}_y,{\mathit{\theta }}_z)$ and $k_{\mathrm{g}\mathit{\_}i}(i=x,y,z,{\mathit{\theta }}_x,{\mathit{\theta }}_y,{\mathit{\theta }}_z)$ denote the values of isotropic coupling stiffness between the slices. The matrix form of the static equilibrium equation of the sliced mesh element is shown in Eq. (17). The stiffness and displacement in the equations are related to the gear mesh position. When the relative displacement of a gear slice along the line of action is negative, the slice is out of contact. In the static equilibrium equation, the mesh stiffness of the corresponding slice should be set to zero. Based on the definition and calculation principles of mesh misalignment specified in clause 7.5.3.4 of ISO 6336-1:2019, the mesh misalignment formulation is further derived and extended according to the structural characteristics of the proposed model, as given in Eq. (18).

Figure 4 Sliced static model incorporating misalignment forces and moments.

Figure 5 Sliced mesh element model.

3.3 Loaded tooth contact analysis model (LTCA model)

The deformation compatibility relationship before and after loading of the double-helical gear pair is shown in Fig. 6. The contact equations under load are established using the deformation compatibility relationship, as shown in Eq. (19), where δ_{L_i} and δ_{R_i} are the deformations at the ith contact point after loading, including macroscopic deformation and contact deformation. x_{L_i} and x_{R_i} are the rigid-body approaches, and d_{L_i} and d_{R_i} are the residual clearances after the contact point is loaded and deformed. ϵ_{L_i} and ϵ_{R_i} represent the initial gaps associated with the contact points, and the initial clearance is caused by the mesh misalignment of the double-helical gear pair due to the system deformation.

$$\begin{array}{}\text{(16)}& {\displaystyle }\left\{\begin{array}{l}\begin{array}{rl}& k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\cos {\mathit{\beta }}_b\sin \mathit{\gamma }-k_{\mathrm{p}\mathit{\_}x}(x_{\mathrm{p}}^{j-\mathrm{1}}-x_{\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}x}(x_{\mathrm{p}}^j-x_{\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\cos {\mathit{\beta }}_b\cos \mathit{\gamma }-k_{\mathrm{p}\mathit{\_}y}(y_{\mathrm{p}}^{j-\mathrm{1}}-y_{\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}y}(y_{\mathrm{p}}^j-y_{\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \mp k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\sin {\mathit{\beta }}_b-k_{\mathrm{p}\mathit{\_}z}(z_{\mathrm{p}}^{j-\mathrm{1}}-z_{\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}z}(z_{\mathrm{p}}^j-z_{\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\sin {\mathit{\beta }}_b\sin \mathit{\gamma }-k_{\mathrm{p}\mathit{\_}\mathit{\theta }x}({\mathit{\theta }}_{x\mathrm{p}}^{j-\mathrm{1}}-{\mathit{\theta }}_{x\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}\mathit{\theta }x}({\mathit{\theta }}_{x\mathrm{p}}^j-{\mathit{\theta }}_{x\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\sin {\mathit{\beta }}_b\cos \mathit{\gamma }-k_{\mathrm{p}\mathit{\_}\mathit{\theta }y}({\mathit{\theta }}_{y\mathrm{p}}^{j-\mathrm{1}}-{\mathit{\theta }}_{y\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}\mathit{\theta }y}({\mathit{\theta }}_{y\mathrm{p}}^j-{\mathit{\theta }}_{y\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\cos {\mathit{\beta }}_b-k_{\mathrm{p}\mathit{\_}\mathit{\theta }z}({\mathit{\theta }}_{z\mathrm{p}}^{j-\mathrm{1}}-{\mathit{\theta }}_{z\mathrm{p}}^j)\\ & +k_{\mathrm{p}\mathit{\_}\mathit{\theta }z}({\mathit{\theta }}_{z\mathrm{p}}^j-{\mathit{\theta }}_{z\mathrm{p}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& -k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\cos {\mathit{\beta }}_b\sin \mathit{\gamma }-k_{\mathrm{g}\mathit{\_}x}(x_{\mathrm{g}}^{j-\mathrm{1}}-x_{\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}x}(x_{\mathrm{g}}^j-x_{\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \mp k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\cos {\mathit{\beta }}_b\cos \mathit{\gamma }-k_{\mathrm{g}\mathit{\_}y}(y_{\mathrm{g}}^{j-\mathrm{1}}-y_{\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}y}(y_{\mathrm{g}}^j-y_{\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)\sin {\mathit{\beta }}_b-k_{\mathrm{g}\mathit{\_}z}(z_{\mathrm{g}}^{j-\mathrm{1}}-z_{\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}z}(z_{\mathrm{g}}^j-z_{\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\sin {\mathit{\beta }}_b\sin \mathit{\gamma }-k_{\mathrm{g}\mathit{\_}\mathit{\theta }x}({\mathit{\theta }}_{x\mathrm{g}}^{j-\mathrm{1}}-{\mathit{\theta }}_{x\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}\mathit{\theta }x}({\mathit{\theta }}_{x\mathrm{g}}^j-{\mathit{\theta }}_{x\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\sin {\mathit{\beta }}_b\cos \mathit{\gamma }-k_{\mathrm{g}\mathit{\_}\mathit{\theta }y}({\mathit{\theta }}_{y\mathrm{g}}^{j-\mathrm{1}}-{\mathit{\theta }}_{y\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}\mathit{\theta }y}({\mathit{\theta }}_{y\mathrm{g}}^j-{\mathit{\theta }}_{y\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\\ \begin{array}{rl}& \pm k_{\mathrm{m}}^j({\mathit{\sigma }}_{\mathrm{m}}^j-e_{\mathrm{m}}^j)r_{\mathrm{p}}\cos {\mathit{\beta }}_b-k_{\mathrm{g}\mathit{\_}\mathit{\theta }z}({\mathit{\theta }}_{z\mathrm{g}}^{j-\mathrm{1}}-{\mathit{\theta }}_{z\mathrm{g}}^j)\\ & +k_{\mathrm{g}\mathit{\_}\mathit{\theta }z}({\mathit{\theta }}_{z\mathrm{g}}^j-{\mathit{\theta }}_{z\mathrm{g}}^{j+\mathrm{1}})=\mathrm{0}\end{array}\end{array}\right.\text{(17)}& {\displaystyle }\left\{\begin{array}{l}\left({\mathbf{K}}_{\mathrm{m}}^j\left(t\right)+{\mathbf{K}}_{\mathrm{slice}}\right)q_{\mathrm{m}}^j\left(t\right)=\mathrm{0}\\ \mathrm{If}{\mathit{\sigma }}_{\mathrm{m}}^j\le \mathrm{0},{\mathbf{K}}_{\mathrm{m}}^j=\mathrm{0}\end{array}\right.\text{(18)}& {\displaystyle }{\mathit{\xi }}_{\mathrm{m}}^j=max\left({\mathit{\sigma }}_{\mathrm{m}}^j\right)-{\mathit{\sigma }}_{\mathrm{m}}^j\text{(19)}& {\displaystyle }\left\{\begin{array}{l}{\mathit{\delta }}_{\mathrm{L}\mathit{\_}i}^{\mathrm{P}}+{\mathit{\delta }}_{\mathrm{L}\mathit{\_}i}^{\mathrm{G}}+{\mathit{ϵ}}_{\mathrm{L}\mathit{\_}i}-x_{\mathrm{L}\mathit{\_}i}-d_{\mathrm{L}\mathit{\_}i}=\mathrm{0}\\ {\mathit{\delta }}_{\mathrm{R}\mathit{\_}i}^{\mathrm{P}}+{\mathit{\delta }}_{\mathrm{R}\mathit{\_}i}^{\mathrm{G}}+{\mathit{ϵ}}_{\mathrm{R}\mathit{\_}i}-x_{\mathrm{R}\mathit{\_}i}-d_{\mathrm{R}\mathit{\_}i}=\mathrm{0}\\ \text{When point}i\text{is in contact,}{d}_{\mathrm{L}\mathit{\_}i}=\mathrm{0},d_{\mathrm{R}\mathit{\_}i}=\mathrm{0}\\ \text{When point}i\text{is not in contact,}{d}_{\mathrm{L}\mathit{\_}i}>\mathrm{0},d_{\mathrm{R}\mathit{\_}i}>\mathrm{0}\end{array}\right.\end{array}$$

The total deformation is partitioned into macroscopic and contact components. Equation (20) gives the loaded contact relationship at each meshing position. The macroscopic deformation compliance matrix of the contact points along the meshing line is denoted by λ. u_contact denotes the contact deformation vector for points along the meshing line. x_STE represents the rigid-body approach vector, equivalent to the static transmission error. d is defined as the residual clearance vector at the contact points. ϵ refers to the initial clearance vector at the contact point. F represents the column vector of loads acting at the contact points, and P is the total load on the contact area. The slice mesh stiffness is obtained by summing the stiffnesses of the contact points within the corresponding gear slice.

$$\begin{array}{}\text{(20)}& \left\{\begin{array}{l}-\mathit{\lambda }\mathit{F}-{\mathit{u}}_{\mathrm{contact}}+{\mathit{x}}_{\mathrm{STE}}+\mathit{d}=\mathit{ϵ}\\ \sum _{i=\mathrm{1}}^n{\mathit{F}}_i=P\\ \mathrm{If}{\mathit{F}}_i>\mathrm{0},d_i=\mathrm{0}\\ \mathrm{If}{\mathit{F}}_i=\mathrm{0},d_i>\mathrm{0}\end{array}\right.\end{array}$$

In the loaded tooth contact analysis of a gear pair, the initial clearance at the contact point originates from the geometric errors of the gears and their assembly. Minor profile and lead deviations, installation tolerances of bearings and housings, and even bearing internal clearance can cause a discrepancy between the actual and theoretical meshing locations of the two gears, thereby forming an initial backlash when unloaded or lightly loaded. During operation, thermal expansion and shaft elastic deformation can further change the magnitude and distribution of this clearance. Because this study examines how coupling misalignment affects the contact and dynamic behavior of the gear system, the contact point clearance is attributed to the gear–mesh misalignment caused by shaft deformation induced by coupling misalignment. In Eq. (20), the contact point clearance is determined by the gear meshing misalignment through Eq. (21). Here, j denotes the number of contact points.

$$\begin{array}{}\text{(21)}& \mathit{ϵ}=\left[{\mathit{\xi }}^{\mathrm{1}}\mathrm{\cdots }{\mathit{\xi }}^j\right]\end{array}$$

Figure 6 Deformation compatibility relationship of the double-helical gear pair.

3.4 Iterative calculation process

To analyze the quasi-static contact characteristics of gear pairs considering coupling misalignment, an iterative procedure is established to solve the static equilibrium equation and LTCA model and to update the mesh misalignment, as shown in Fig. 8. The corresponding steps are outlined below:

1. Discretize the gear system. Obtain the misalignment force vector from the coupling-misalignment model. Use the LTCA model to calculate the initial mesh stiffness of the sliced mesh elements.

2. Establish the static equilibrium equations for the gear system and solve them to obtain the initial deformation values of the gear nodes.

3. Determine the gear pair mesh misalignment from the nodal deformations.

4. Assign the misalignment-induced initial clearances to the contact points and solve the LTCA model to obtain the stiffness of each corresponding sliced mesh element.

5. Update the static equilibrium equation and solve for the nodal deformations.

6. Determine whether the relative error between the nodal deformations obtained in successive iterations is less than the convergence tolerance. Failure to satisfy the prescribed criterion redirects the solution procedure to step 3 for another update cycle; once the iterative sequence converges, the resulting mesh stiffness, static transmission error, tooth surface load distribution, and mesh misalignment are retained as the final outputs.

The model was experimentally verified with a test rig specifically designed to investigate the dynamic characteristics of a double-helical gearbox under various loading conditions. The setup, as detailed by Gong et al. (2022a), incorporates a purpose-built, isolation-mounted, power-circulating double-helical gear test rig, which allows for precise measurement of vibration acceleration and torque under controlled conditions. The experimental results, which focus on vibration response, torque effects, and resonance frequencies, closely match the model's predicted behavior. This agreement between the experimental and numerical results supports the capability of the adopted model to reproduce the overall dynamic behavior of the gear transmission system, providing a solid foundation for future gear system analysis and design. It should be noted that the previous experimental validation focused mainly on the overall dynamic response of the basic gear transmission model, while the coupling-misalignment-related effects introduced in the present study are investigated on the basis of this validated modeling framework, with corresponding experimental investigations to be conducted in subsequent work.

4 Quasi-static characteristics of a gear system considering coupling misalignment

4.1 Effect of coupling misalignment

This section considers six coupling-misalignment cases. ϵ_x is 0.05, ϵ_y is 0.05, ϵ_z is 0.25 mm, and γ is 0.25° for case 1. The amount of misalignment for case 2 to case 6 is 2 to 6 times that of case 1, respectively. First, the effect of each misalignment parameter is analyzed using case 1 as the reference condition. While keeping all other parameters constant, one misalignment parameter at a time is increased by specified factors to investigate the independent effects of each parameter on the quasi-static contact characteristics of the gear pair. As shown in Fig. 7, the horizontal axis represents different coupling-misalignment parameters, and the vertical axis represents the magnification factor relative to the case-1 reference condition, where 1× denotes the reference condition, and 2× to 5× represent values that are 2 to 5 times the reference value, respectively. This single-factor control method effectively eliminates the coupled effects caused by simultaneous changes in different misalignment parameters, thereby allowing for a comparison of the independent effects of each parameter on meshing stiffness. The results show that the extent of the influence of different coupling-misalignment parameters on meshing stiffness varies significantly, with ϵ_x having the most significant effect, followed by ϵ_y, while γ has a relatively minor effect, and ϵ_z has no significant effect on meshing stiffness. As ϵ_x and ϵ_y increase, the meshing stiffness exhibits a generally pronounced upward trend, indicating a strong positive correlation between these parameters and meshing stiffness; in contrast, the effect of γ on meshing stiffness is weaker, and its variation does not follow a strictly monotonic pattern. The primary reason for this is that the meshing misalignment of a gear pair is more sensitive to torsional deformation of the shaft system, whereas deformation in the axial direction has a relatively limited effect on meshing misalignment. Furthermore, different misalignment parameters affect the form and extent of shaft system deformation in varying ways, which in turn leads to differing degrees of influence on meshing misalignment and meshing stiffness. It can thus be seen that the quasi-static contact characteristics of the gear pair are primarily sensitive to ϵ_x and ϵ_y, while the effects of γ and ϵ_z are relatively limited. It should be noted that this section primarily analyzes the impact of changes in coupling-misalignment parameters on the quasi-static contact characteristics of the gear pair, with a focus on the variation patterns of meshing stiffness. Since meshing stiffness is a critical internal excitation parameter in subsequent dynamic analyses, its variation will further influence the system's dynamic response. Therefore, the patterns of parameter sensitivity revealed by this single-factor quasi-static analysis can provide a basis for subsequently assessing the relative effects of different coupling-misalignment parameters on the system's dynamic response.

Figure 7 Changes in mesh stiffness under different misalignment parameters.

Figure 8 Iterative calculation process.

As shown in Fig. 9, the g2 gear pair mesh stiffness curve exhibits not only short-period fluctuations, that is, fluctuations at the mesh frequency, but also obvious long-period fluctuations. The fluctuation frequency is the same as the coupling-misalignment force frequency, which is twice the shaft frequency. The results indicate that the contact condition of the g2 gear pair is governed by the force/moment excitation induced by coupling misalignment. With increasing misalignment, the shaft undergoes greater bending deformation, which further intensifies the fluctuation of mesh stiffness over one rotational period. For the left helical gear, since complete tooth surface contact is established during meshing cycle intervals 42 to 86 and 142 to 186, the actual contact line length does not change, and so the mesh stiffness is the same. Outside these intervals, the reduction in mesh stiffness results from partial disengagement of the tooth surface, and the curve shape changes significantly. For the right helical gear, during mesh cycles 42 to 86 and 142 to 186, the actual contact line length changes due to partial disengagement of the tooth surface, which markedly alters both the level of mesh stiffness and its fluctuation characteristics.

Figure 9 Mesh stiffness of g2 gear pair under different coupling misalignments.

Since the misalignment force/moment of the coupling changes periodically with variations in the rotational angle, the load-induced deformation of the system also varies over the period corresponding to twice the shaft frequency. Figure 10 shows the mesh misalignment of the g2 gear pair in the 13th and 66th mesh cycles. In the 13th mesh cycle, for the g2 gear pair, the left tooth exhibits markedly greater mesh misalignment than the right, and this left-side misalignment increases progressively with increasing coupling misalignment. Misalignment-induced forces and moments have an increasing influence on the bending deformation of the shaft carrying the g2 gear pair, causing the disengaged area on the left tooth surface to increase and the actual contact line length to decrease. Therefore, the mesh stiffness of the left tooth decreases significantly with the increase in misalignment in the 13th mesh cycle. Because the peak mesh misalignment on the right tooth is only 11.14 µm across the six prescribed misalignment levels and because the tooth surface remains in complete contact during the 13th meshing cycle, the right mesh stiffness remains almost unchanged as the misalignment level increases.

Figure 10 Mesh misalignment of the g2 gear pair during the 13th and 66th mesh cycles.

With increasing coupling misalignment, the force/moment excitation induced by misalignment exerts a progressively stronger effect on the bending deformation of the shaft carrying the g2 gear pair. In the 66th mesh cycle, as the coupling-misalignment level rises, the g2 left tooth exhibits a gradual reduction in mesh misalignment magnitude; in case 6, the peak misalignment is displaced toward the left-hand region of the tooth surface. The mesh misalignment of the right tooth surface gradually increases with the increase in coupling misalignment. Because the left tooth mesh misalignment is smaller in cases 2 to 6, the mesh stiffness of the left tooth of the g2 gear pair is the same in this mesh cycle; the curve shape is also the same, and the tooth surfaces are in complete contact. Figure 11 presents the tooth surface load pattern calculated for the g2 gear pair during the 66th mesh cycle under different coupling-misalignment amounts. In case 1, disengagement is observed only on the left tooth surface, and the mesh stiffness is significantly lower than that in the other misalignment cases. The shaft carrying the g2 gear pair bends under the combined action of the meshing force, and the force/moment excitation induced by misalignment. The bending deformation caused by its meshing force is symmetrically distributed along the axis. As the misalignment force/moment increases, the influence on shaft bending deformation gradually increases, and the differences in mesh misalignment and tooth surface load distribution between the left and right increase accordingly.

Figure 11 Tooth surface load distribution of the g2 gear pair during the 66th mesh cycle.

Figure 12 Mesh stiffness of gear pairs under different load torques.

By comparing the tooth surface contact states shown in Fig. 11 with the face width mesh misalignment at the 66th meshing cycle in Fig. 10, it can be observed that the central relief groove of the herringbone gear separates the left- and right-hand helical tooth surfaces. Therefore, its meshing behavior can be approximately regarded to be the combined action of two single-helical gear pairs, allowing the mesh misalignment and contact states of the left and right tooth surfaces to be analyzed separately. It should be noted that the mesh misalignment of the gear pair is not governed by a single factor but results from the combined effects of the tooth meshing force under the actual loading condition and the additional forces and moments induced by coupling misalignment. Consequently, the contact states of the left and right tooth surfaces do not exhibit a simple positive correlation with the magnitude of coupling misalignment. In case 1, the coupling eccentricity is relatively small, and the additional forces and moments induced by the misalignment have only a limited influence on the shaft system attitude. Therefore, the mesh misalignment of the gear pair is mainly governed by the tooth meshing force under the actual loading condition. Although the mesh misalignment values of the left and right tooth surfaces are relatively close, the misalignment on the left tooth surface is slightly larger, resulting in local disengagement, whereas the right tooth surface maintains a relatively favorable contact state. This indicates that tooth disengagement in a herringbone gear cannot be determined solely by the magnitude of coupling eccentricity but depends primarily on the actual distribution of mesh misalignment between the left and right tooth surfaces under the combined effects of the tooth meshing force and coupling-misalignment excitation. Further examination of the results at the 66th meshing cycle in Fig. 10 shows that, as the coupling eccentricity increases, the additional forces and moments induced by the misalignment gradually increase, and their influence on the displacement and orientation of the shaft system becomes more pronounced. As a result, the original mesh misalignment distribution dominated by the tooth meshing force is progressively altered. During this process, the mesh misalignment of the left and right tooth surfaces exhibits opposite trends: the mesh misalignment on the left tooth surface gradually decreases with increasing misalignment excitation, whereas that on the right tooth surface gradually increases. Correspondingly, the tooth-surface contact states shown in Fig. 11 also change significantly, evolving from local disengagement on the left tooth surface with good contact on the right tooth surface in case 1 to good contact on the left tooth surface with disengagement on the right tooth surface in case 6. These results demonstrate that, for a herringbone gear transmission system, the magnitude of coupling misalignment itself does not directly determine the tooth surface contact state. Instead, coupling misalignment modifies the additional forces and moments acting on the shaft system, which interact with the tooth meshing force to redistribute the mesh misalignment between the left and right tooth surfaces, ultimately leading to corresponding changes in the tooth contact region and local load-carrying state. The above analysis shows that the system deformation and tooth surface contact state vary under coupling-misalignment forces and moments over the period corresponding to twice the shaft frequency. The mesh stiffness exhibits long-period fluctuations. With increasing misalignment magnitude, the mesh stiffness fluctuation of the gear pair, which is greatly affected by the misalignment force/moment, will also increase. Larger misalignment forces and moments lead to greater differences in mesh misalignment between the left and right tooth surfaces and a more asymmetric tooth surface load distribution. Under the combined effect of meshing force and misalignment force/moment, the mesh misalignment of the gear pair within certain mesh cycle ranges is small. When the tooth surface reaches a complete contact state, the actual contact line length is the same, and the mesh stiffness within this range is the same.

4.2 Effect of load torque

The load torque determines the meshing force of the gear pair. Shaft deformation at the gear pair location is governed jointly by meshing force excitation and the force/moment components generated by coupling misalignment. In this section, the coupling misalignment ϵ_x is 0.15, ϵ_y is 0.15, ϵ_z is 0.75 mm, and γ is 0.75°. The quasi-static characteristics of the gear system are analyzed at load torques of 22 705, 15 137, 11 353, 7568, 5676, and 3784 N m⁻¹. Figure 12 shows the mesh stiffness of the gear pair under different load torques. As the load torque gradually increases, the g2 and g3 gear pairs exhibit progressively reduced mesh stiffness fluctuations. Owing to the shorter distance between the g2 left helical gear and the coupling, it is more affected by coupling misalignment. Across the load torque levels considered, mesh stiffness fluctuation is more pronounced in the left helical gear than in the right. The tooth surfaces of the g2 left helical gear are in complete contact during certain mesh cycles, and so the mesh stiffness value does not change significantly. The g2 right helical gear has a disengaged area on the tooth surface during certain mesh cycles, and the mesh stiffness value changes significantly. Since the g3 gear pair has partial disengagement during the rotation period, the long-term fluctuation trend of the mesh stiffness curve under different load torques is sinusoidal. The greater distance between the g3 gear pair and the coupling weakens the influence of the misalignment-induced force. Therefore, the g3 gear pair has smaller long-term fluctuations compared with the g2 gear pair.

Figure 13 Gear mesh misalignment under different load torques during the 66th mesh cycle.

Figure 14 The rms dynamic meshing force of the g2 left helical gear.

Figure 15 Spectral map of the dynamic meshing force of the g2 gear pair over 0 to 4000 r min⁻¹.

In general, higher load torque produces larger system deformation under load and correspondingly greater mesh misalignment. The formation of a disengagement region on the tooth surface reduces the actual contact line length; the mesh stiffness value will decrease, and the fluctuation amount will increase. However, when the coupling is misaligned, the misalignment force/moment and the gear pair meshing force jointly affect the system deformation. The meshing force is very small under a light load. The system deformation is mainly caused by the coupling-misalignment force/moment. Consequently, an extensive region of contact loss develops on the tooth surface. The actual contact line length changes with the mesh cycle, and the mesh stiffness value also changes accordingly, and so the long-period fluctuation increases under light loads. As the load torque gradually increases, the gear pair meshing force also increases. As the meshing force increases, the disengaged area on the tooth surface decreases until full contact is reached. Only limited variation occurs in the actual tooth surface contact line length, and so the mesh stiffness value changes little, and the long-period fluctuation is reduced.

Figure 13 shows the mesh misalignment of the gear pair under different load torques during the 66th mesh cycle. As the load torque increases, the influence of the g2 gear pair meshing force on the bending deformation of the shaft increases, and the mesh misalignment discrepancy between the left and right teeth is gradually reduced. For the g3 gear pair, the mesh misalignment is symmetrically distributed under the six load torques and is less affected by the misalignment force. As the torque increases, the meshing force increases, increasing the mesh misalignment. The proximity of the right tooth surface to the power output end causes the mesh misalignment associated with torsional deformation of the shaft to be slightly larger than that on the left tooth surface.

Figure 16 The rms of bearing reaction force.

This section shows that coupling misalignment causes the mesh stiffness of the gear pair to fluctuate at twice the shaft frequency. Since mesh stiffness is an important source of excitation in the gear system, it can be inferred that the gear system dynamic model incorporating coupling misalignment presented below contains an excitation component at twice the shaft frequency.

5 Dynamic characteristics of a gear system considering coupling-misalignment excitation

5.1 System differential equation of motion

According to the study by Gong et al. (2022b), the differential equation of motion of the gear system considering coupling misalignment is established, and its matrix form is given in Eq. (22). t is the time required for the slowest shaft in the system to rotate one revolution. The coupling-misalignment forces and moments are included as excitation terms in the gear system dynamics model, where F_c(t) represents the coupling-misalignment excitation vector. According to the above misalignment excitation equation, it is a low-frequency excitation associated with the rotational angle. M_global, C_global, and K_global are the system mass, damping, and stiffness matrices, respectively. x, $\dot{x}$, and $\ddot{x}$ are displacement, velocity, and acceleration, respectively. e(t) is the meshing error. Since the research object of this paper includes multiple meshing excitations and coupling-misalignment excitations, to improve calculation efficiency, the parameters of Eq. (22) are normalized as shown in Eq. (23). Here, ${\stackrel{\mathrm{‾}}{\mathbf{K}}}_{\mathrm{global}}$ represents the mean stiffness, and Δ denotes its fluctuation. The system excitation force is calculated using the frequency domain solution method. The system vibration responses to time-varying mesh stiffness excitation, error excitation, and coupling-misalignment excitation in each gearbox are calculated separately using Fourier series expansions. The mesh stiffness, transmission error, and coupling-misalignment excitations are all expanded in Fourier series, with the gear meshing frequency taken as the fundamental frequency. The detailed formulation of the dynamic model and the corresponding solution procedure are provided in Gong et al. (2025). The system vibration response is then obtained by superposing the vibration responses associated with identical nodes and corresponding degrees of freedom. The dynamic meshing force of the gear pair is calculated using Eq. (24), where k_m and c_m denote the mesh stiffness and mesh damping coefficient, respectively; e_m is the composite mesh error; and DTE is the dynamic transmission error. The vector q_m represents the nodal vibration displacements of the gears, and V is the projection vector along the line of action.

$$\begin{array}{}\text{(22)}& {\displaystyle }\begin{array}{rl}& {\mathbf{M}}_{\mathrm{global}}\stackrel{\mathrm{..}}{x}\left(t\right)+{\mathbf{C}}_{\mathrm{global}}\stackrel{.}{x}\left(t\right)+{\mathbf{K}}_{\mathrm{global}}\left(t\right)\left[x\left(t\right)-e\left(t\right)\right]\\ & ={\mathit{P}}_{\mathrm{0}}+{\mathit{F}}_{\mathrm{C}}\left(t\right)\end{array}\text{(23)}& {\displaystyle }\begin{array}{rl}& {\mathbf{M}}_{\mathrm{global}}\mathrm{\Delta }\stackrel{\mathrm{..}}{x}+{\mathbf{C}}_{\mathrm{global}}\mathrm{\Delta }\stackrel{.}{x}+{\stackrel{\mathrm{‾}}{\mathbf{K}}}_{\mathrm{global}}\mathrm{\Delta }x\\ & ={\stackrel{\mathrm{‾}}{\mathbf{K}}}_{\mathrm{global}}\mathrm{\Delta }e-\mathbf{\Delta }\mathbf{K}\left[x_{\mathrm{STE}}\left(t\right)-e\left(t\right)\right]+{\mathit{F}}_{\mathbf{C}}\left(t\right)+{\mathit{P}}_{\mathrm{0}}\end{array}\text{(24)}& {\displaystyle }{F}_{\mathrm{d}}=k_{\mathrm{m}}\left(\mathrm{DTE}-e_{\mathrm{m}}\right)+c_{\mathrm{m}}\left(\mathrm{D}\dot{\mathrm{T}}\mathrm{E}-{\dot{e}}_{\mathrm{m}}\right)\text{(25)}& {\displaystyle }\mathrm{DTE}=\mathit{V}\cdot {\mathit{q}}_{\mathrm{m}}\end{array}$$

5.2 Effect of combined high- and low-frequency excitation

Figure 14 shows the root mean square (rms) of the dynamic meshing force of the g2 left helical gear. Coupling-misalignment excitation produces a new resonance peak for the g2 gear pair at 2350 r min⁻¹. In the resonance region induced by misalignment excitation, the linear model predicts a marked increase in dynamic meshing force fluctuations. The dynamic meshing force spectral map of the g2 gear pair is shown in Fig. 15. A distinct frequency band at twice the shaft frequency (2f_n) is observed in the spectrum. f_n denotes the rotational frequency of the shaft carrying the coupling. The maximum amplitude in the dynamic meshing force spectral map of the g2 gear pair appears at 2350 r min⁻¹ in the band at twice the shaft frequency. The mesh frequency band, f_m1, of the adjacent gearbox GB1 appears in the dynamic meshing force spectrum of the g2 gear pair, which shows a coupling resonance between the GB1 and GB2 gearboxes at 2300 r min⁻¹.

Figure 17 The rms of the dynamic meshing force under different coupling misalignments.

The bearing reaction force under coupling misalignment is characterized in Fig. 16 using its rms value. The misalignment of the coupling significantly increases the root mean square reaction force of the b3 and b4 bearings in the GB1 gearbox that are closer to the coupling, and a larger resonance peak appears at 3100 r min⁻¹. The misalignment of the coupling significantly increases the root mean square reaction force of the b5 to b10 bearings in the GB2 gearbox, and large resonance peaks appear at 2350 and 3100 r min⁻¹. As the distance from the coupling gradually increases, the root mean square of the bearing reaction force gradually decreases, and the vibration weakens. The vibration response of the b9 and b10 bearings farthest from the coupling is also affected by coupling misalignment.

5.3 Effect of coupling misalignment

Figure 17 shows the root mean square of the dynamic meshing force under different coupling misalignments. For the g1 right helical gear, the RMS dynamic meshing force curves coincide under the aligned condition and cases 1 to 3. This is because the GB1 gearbox bearings have high stiffness and are unaffected by coupling misalignment. The mesh frequency causes the resonance peak at 1350 r min⁻¹. The dynamic meshing force of the g2 left helical gear fluctuates periodically at twice the shaft frequency. The amount of fluctuation increases with the amount of misalignment. Therefore, the RMS dynamic meshing force of the g2 left helical gear increases as the coupling misalignment becomes larger. Notably, when the coupling-misalignment excitation produces resonance at 2350 r min⁻¹, the resonance peaks for the aligned condition and cases 1 to 3 are 909.06, 30 579.46, 60 550.19, and 91 213.82 N, respectively. These results indicate a progressive enhancement of the vibration response of the g2 left helical gear. By contrast, the g3 left helical gear shows little variation in RMS dynamic meshing force within 0 to 850 r min⁻¹ as the coupling misalignment increases. This is because the vibration in this speed range is mainly caused by mesh frequency excitation and is less affected by coupling-misalignment excitation. Between 850 and 4000 r min⁻¹, resonance is caused by coupling-misalignment excitation. As misalignment increases, the resonance peak amplitude increases, and gear pair vibration intensifies. The change in coupling misalignment does not change the resonance speed of the gear pair.

Figure 18 Dynamic meshing force spectrum of the g2 left helical gear under different coupling misalignments at 900 r min⁻¹.

Figure 19 The rms dynamic meshing force under different load torques.

Figure 18 shows the meshing force spectrum of the g2 left helical gear under different coupling misalignments. The corresponding amplitudes of the component at twice the shaft frequency, 2f_n, in cases 1 to 3 are 1936.64, 3337.74, and 3735.16 N, respectively. This shows that, as the amount of misalignment increases, the impact of coupling-misalignment excitation on the vibration response of the gear pair increases. In the dynamic meshing force spectrum under cases 1∼3, the amplitude of the adjacent gearbox mesh frequency f_m1 and the amplitude of its mesh frequency f_m2 do not change significantly with the change in the misalignment amount. The component at the mesh frequency f_m2 is modulated by the component at twice the shaft frequency, 2f_n, producing a sideband near f_m2.

5.4 Effect of load torque

The coupling misalignment ϵ_x is 0.15, ϵ_y is 0.15, ϵ_z is 0.75 mm, and γ is 0.75°. On this basis, the vibration response characteristics of the gear system under different load torques are analyzed in Fig. 19. Three load torques of 5676, 11 353, and 22 705 N m⁻¹ are selected. For the g1 right helical gear, the rms dynamic meshing force increases with load torque, and the increase is more obvious at the resonance speed. This is because the resonance peak of the g1 right helical gear in the range of 0 to 4000 r min⁻¹ is excited by its mesh frequency. Higher torque produces a larger excitation force in the gear pair, a higher rms dynamic meshing force, and more severe vibration. The resonance peak amplitude caused by the mesh frequency of the g2 and g3 gear pairs below 1000 r min⁻¹ increases with the increase in load torque. Since the misalignment excitation of the coupling is the same under different load torques, the resonance peak value of the g2 and g3 gear pairs caused by the misalignment excitation does not change much under the three load torques. The load torque has little effect on the resonance speed of the gear pair. It mainly affects the vibration response under mesh frequency excitation. The greater the torque, the more severe the vibration.

6 Discussion and conclusions

This paper proposes an efficient quasi-static contact analysis model for marine gear systems considering coupling misalignment. The coupling misalignment excitation is introduced into the gear system differential equation of motion. The quasi-static contact behavior and dynamic characteristics of the marine gear system are analyzed under different coupling-misalignment and load torque conditions. The main research conclusions are as follows:

1. Owing to periodic variations in coupling-misalignment forces and moments, the mesh stiffness fluctuates over long periods. As the amount of misalignment increases, long-period fluctuations in mesh stiffness will also increase. A higher level of misalignment-induced force/moment excitation produces more pronounced shaft bending deformation, accompanied by larger differences between the left and right tooth surfaces in mesh misalignment and tooth surface load distribution of the double-helical gear pair.

2. The gear pair meshing force and misalignment force/moment jointly affect the system deformation. Under light loads, the quasi-static contact behavior of the gear pair is more susceptible to coupling misalignment, causing the tooth surface to disengage and the long-period fluctuations in mesh stiffness to increase. As the torque gradually increases, system deformation is primarily governed by the meshing force. For the symmetrically supported double-helical gear pair, the mesh misalignment is symmetrically distributed between the left and right tooth surfaces, and the disengaged area decreases until it disappears.

3. As the rotational speed gradually increases, the coupling misalignment excitation causes the system to resonate. In the resonance region induced by misalignment excitation, the linear model predicts a marked increase in dynamic meshing force fluctuations. Under the action of combined high- and low-frequency excitations, the modulation phenomenon is evident in the vibration response spectrum of the gear pair. When the coupling is misaligned, the bearings on both sides exhibit strong vibration dominated by a component at twice the shaft frequency.

4. As the coupling misalignment increases, the system vibration gradually intensifies. The change in coupling misalignment does not change the resonance speed of the gear pair. The load torque has little effect on the system resonance speed. It mainly affects the system vibration response under mesh-frequency excitation.

In future research, we will focus on enhancing experimental validation, exploring additional types of misalignments, and evaluating long-term system performance. Furthermore, nonlinear effects will be addressed by incorporating gear backlash and two-way iterative coupling between tooth contact states and vibration responses, with surrogate models used to improve computational efficiency. Extending the analysis to varying lubrication regimes and multiple operating scenarios will further enhance its practical applicability, supporting reliability assessment and maintenance planning for marine gear systems. These future directions build upon the findings of this work and aim to refine the model for real-world implementation.