Nonlinear Dynamic Response Analysis of Cable–Buoy Structure Under Marine Environment

Abstract

The nonlinear dynamics of the cable–buoy structure in marine engineering present significant analytical challenges due to the complex motion of the buoy, which impacts the system’s dynamic response. The drag force acting on the structure can be categorized into the absolute velocity and relative velocity models, distinguished by their reference frames. The absolute velocity model incorporates flow velocity coupling terms, offering higher accuracy but at the expense of increased computational complexity. In contrast, the relative velocity model is computationally simpler and therefore more widely adopted. Nevertheless, the accuracy and applicability of these simplified models remain open to further in-depth investigation. To address these challenges, this study derives coupled differential equations for the cable–buoy structure based on the two drag force models. Galerkin discretization is then employed to construct coupled systems that account for nonlinear buoy motion, as well as decoupled systems assuming linear buoy motion. The modulation equations for the system’s primary resonance response are derived using the method of multiple scales. Numerical results indicate that changes in cable parameters lead to complex modal coupling behaviors in the system. The flow velocity coupling terms in the absolute velocity drag force model enhance the system’s damping effect, and the relative velocity drag force model, which omits these coupling terms, results in increased system response amplitudes. Although neglecting nonlinear buoy motion has little impact on the cable’s dynamic response, it significantly reduces the amplitude of the buoy’s dynamic motion. The relative velocity drag force model and the decoupled system can serve as effective simplifications for analyzing the dynamic responses of cable–buoy systems, providing a balance between computational efficiency and result accuracy. Variations in system parameters cause both qualitative and quantitative changes in the system’s nonlinear stiffness characteristics.

1. Introduction

Cable–buoy structures play a crucial role in marine engineering, supporting ocean observation, offshore energy development, and underwater infrastructure deployment [1,2]. These systems operate in highly dynamic marine environments, and are subjected to complex forces such as ocean currents and nonlinear fluid–structure interactions. Accurate nonlinear vibration analysis is essential for understanding their dynamic behavior, optimizing design parameters, and improving system stability and reliability [3,4].

An important aspect of analyzing the dynamic behavior of cable–buoy systems is the accurate representation of drag forces in marine environments. The drag force acting on these systems is typically modeled using the Morison equation, which can be formulated in two ways based on the choice of reference frame: the absolute velocity model and the relative velocity model. The absolute velocity drag force model, which adopts a stationary reference frame, determines drag based on the relative velocity between the system and the environmental flow. This model provides a comprehensive representation of fluid–structure interactions [5], making it particularly suitable for dynamic analyses under complex marine conditions, such as ocean currents and waves [6]. However, the absolute velocity drag force model results in relatively complex equations of motion, which are typically analyzed through numerical simulations implemented in finite element software [7,8]. To enhance computational accuracy, correction factors, such as variations in cable length and diameter, are incorporated [9]. Consequently, the absolute velocity model is widely employed to investigate the effects of parameters like mooring tension, vertical displacement, and tilt angle on the system’s dynamic response [10,11,12,13,14], offering valuable theoretical support for structural optimization. However, the inherent complexity of this model often precludes direct theoretical solutions.

To address the computational challenges associated with the absolute velocity model, some researchers have modified the linear drag force model to improve its accuracy, making it comparable to the absolute velocity drag force model in specific scenarios [15,16,17]. This modification typically introduces a correction factor, which, however, requires the simulation of nonlinear drag forces to accurately calculate. In contrast, the relative velocity drag force model offers a more computationally efficient alternative. By adopting the fluid as the reference frame, this model calculates drag based solely on the system’s velocity, eliminating the need for coupling terms with flow velocities [18,19,20,21]. This simplification significantly reduces analytical complexity while maintaining reasonable accuracy for most applications. The relative velocity model is particularly advantageous in frequency-domain methods, where it facilitates the application of equivalent linearization and perturbation techniques to obtain analytical solutions [22]. This approach enhances computational efficiency, enabling more comprehensive nonlinear dynamic analyses of submerged cable systems.

The second challenge lies in the choice between coupled and decoupled system modeling approaches. Under the influence of different types of wave forces, the buoy exhibits distinct response characteristics, such as periodic harmonic response, symmetric bifurcation, and non-periodic response [23,24,25]. Meanwhile, the streamwise motion of the cable is more sensitive to variations in excitation amplitude. As excitation amplitude increases, the streamwise motion exhibits more intense and chaotic behavior [26,27]. Given these dynamic behaviors, the cable–buoy system can be modeled as a parametric system to better analyze the response under varying excitation conditions [28,29,30]. However, when excitation is applied only to the cable and not the buoy, the complexity of the governing equations depends on whether the buoy’s nonlinear motion is considered. In Galerkin discretization, considering the nonlinear motion of the buoy results in a coupled system, while assuming linear motion leads to a decoupled system [31]. In the coupled system, there are both linear and nonlinear coupling terms between the buoy and the cable, which limit the perturbation accuracy in nonlinear vibration analysis [32]. In contrast, the decoupled system adopts the analytical approach of linear vibration analysis [33], achieving significantly improved perturbation accuracy in the multiscale method by neglecting only the buoy’s nonlinear motion. However, neglecting the buoy’s nonlinear motion may impact the accuracy of the system’s dynamic response. Therefore, a comprehensive comparison of the dynamic responses of the coupled and decoupled systems is necessary to assess the trade-offs between computational efficiency and model accuracy.

In light of the above discussion, this paper presents a novel simplified model that combines the relative velocity drag force model with a decoupled system to analyze the nonlinear dynamic response of the cable–buoy system. The relative velocity drag force model is particularly well suited for obtaining approximate solutions through perturbation methods, while the decoupled system improves perturbation accuracy during the approximation process. A comparison of the response errors between the relative velocity and absolute velocity drag force models is made, alongside an evaluation of the applicability of the decoupled system relative to the coupled system. The findings offer a new approach for the simplified analysis of cable–buoy systems, characterized by the coupling of rigid and flexible structures, and provide valuable insights for improving the accuracy and efficiency of dynamic analyses in marine engineering applications.

2. Problem Formulation

2.1. Dynamic Equation

In this section, the coupled cable–buoy model is established. It consists of cables anchored to the seabed at one end and attached to the buoy at the other, as illustrated in Figure 1. Here, an O-xy coordinate system is introduced, where the displacement functions u(s,t) and v(s,t) represent the displacements along the tangential direction t and the normal direction n of the cable in the local coordinate system, respectively. Let V and ρ_f denote the uniform fluid velocity and fluid density; E, L, D and ρ_c represent the Young’s modulus, length, diameter and density of the cable; M and D_M are the mass and diameter of the buoy.

Assuming the cable is a uniform, linearly elastic one-dimensional continuum, its torsional stiffness, bending stiffness, and shear stiffness are neglected. Based on Hamilton’s principle, the equations of motion for the system depicted in Figure 1 are written as follows [33,34]:

$$\ddot{u}\left(1,t\right)+{\displaystyle \frac{c_u}{M_a}}\dot{u}\left(1,t\right)+{\omega }_1^{2}u\left(1,t\right)={\omega }_1^2{\displaystyle {\int }_0^1\left(\kappa v-\frac{1}{2}{v^{\prime }}^2\right)\mathrm{d}s}+{\displaystyle \frac{H_{DT}-H_{DT}^i}{M_{a}L}}$$

(1)

$$\begin{array}{c}\ddot{v}+{\displaystyle \frac{c_v}{m_a}}\dot{v}={\displaystyle \frac{P}{m_a{L}^2}}{v}^{″}+{\displaystyle \frac{EA}{m_a{L}^2}}u\left(1,t\right)v^{″}-\left[{\displaystyle \frac{EA}{m_a{L}^2}}{\displaystyle {\int }_0^1\left(\kappa v-\frac{1}{2}{v^{\prime }}^2\right)\mathrm{d}s}\right]v^{″}\\ +{\displaystyle \frac{\kappa EA}{m_a{L}^2}}\left[u\left(1,t\right)-{\displaystyle {\int }_0^1\left(\kappa v-\frac{1}{2}{v^{\prime }}^2\right)\mathrm{d}s}\right]+{\displaystyle \frac{F_{DN}-F_{DN}^i}{m_{a}L}}+{\displaystyle \frac{F_v\left(s,t\right)}{m_{a}L}}\end{array}$$

(2)

with the boundary condition at s = 0 is v(0,t) = 0, and at s = 1, the boundary condition is

$${\displaystyle \frac{M_a}{m_{a}L}}\ddot{v}\left(1,t\right)={\displaystyle \frac{H_{DN}-H_{DN}^i}{m_a{L}^2}}-{\displaystyle \frac{P}{m_a{L}^2}}{v}^{\prime }\left(1,t\right)-{\displaystyle \frac{EA}{m_a{L}^2}}\left[u\left(1,t\right)-{\displaystyle {\int }_0^1\left(\kappa v-\frac{1}{2}{v^{\prime }}^2\right)\mathrm{d}s}\right]v^{\prime }\left(1,t\right)$$

(3)

where the equation above uses the dimensionless variables s* = s/L, κ* = κL, u* = u/L and v* = v/L, with the asterisks omitted for simplicity in the expressions. M_a = M + C_aBρ_fV_B represents the buoy mass, including the added fluid mass, where C_aB is the added-mass coefficient of the buoy, and V_B is the buoy’s volume. m_a = (ρ_c + C_acρ_f)A represents the mass per unit length of the cable, including the added fluid mass, where C_ac is the added-mass coefficient of the cable and A is the cross-sectional area of the cable. ${\mathsf{\omega }}_1^2$ = EA/(M_aL) is the vibration frequency of the buoy modeled as a particle. c_u and c_v represent the damping coefficients of the system. F_v(s,t) = f_vcos(Ωt) represents the excitation transmitted to the cable by ships, waves, or internal waves. P and κ represent the tension and curvature of the cable, respectively, with the assumption that P(s) ≌ P(L) and κ(s) ≌ κ(L) in the subsequent analysis [33]. For the sake of continuity in the subsequent analysis, it is assumed that the external excitation does not affect the buoy. Finally, F_DN, H_DT, and H_DN represent the drag forces, where the superscript ( )ⁱ denotes the fluid forces at static equilibrium, and their detailed expressions are as follows.

The drag force acting on the cable element in the underwater direction along the cable normal is denoted as F_DN, while H_DT and H_DN represent the drag forces acting on the buoy in the underwater tangential and normal directions at the cable end, respectively. For the absolute velocity drag force model, the dimensionless expression for the drag force is given by

$$\begin{array}{c}{F}_{DN}={\displaystyle \frac{1}{2}}{\rho }_f{C}_{D}D\left(-V\sin \theta -L\dot{v}\right)\left|-V\sin \theta -L\dot{v}\right|\\ H_{DT}={\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B\left[V\cos {\theta }_1-L\dot{u}\left(1,t\right)\right]\left|V\cos {\theta }_1-L\dot{u}\left(1,t\right)\right|\\ H_{DN}={\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B\left[-V\sin {\theta }_1-L\dot{v}\left(1,t\right)\right]\left|V\sin {\theta }_1-L\dot{v}\left(1,t\right)\right|\end{array}$$

(4)

where θ₁ represents the angle value at s = 1.

In the relative velocity drag force model, the total drag force consists of two components: the fluid forces under static equilibrium and the drag force due to the relative velocity. The dimensionless expression for the drag force is given by

$$\begin{array}{c}{F}_{DN}=F_{DN}^i-{\displaystyle \frac{1}{2}}{\rho }_f{C}_{D}D\left(L\dot{v}\left|L\dot{v}\right|\right)\\ H_{DT}=H_{DT}^i-{\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B\left[L\dot{u}\left(1,t\right)\left|L\dot{u}\left(1,t\right)\right|\right]\\ H_{DN}=H_{DN}^i-{\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B\left[L\dot{v}\left(1,t\right)\left|L\dot{v}\left(1,t\right)\right|\right]\end{array}$$

(5)

where $F_{\mathit{DN}}^i$ represents the fluid force acting on the cable along the normal direction in static equilibrium, while $H_{\mathit{DT}}^i$ and $H_{\mathit{DN}}^i$ represent the fluid forces acting on the buoy along the tangential and normal directions of the cable in static equilibrium, respectively. Their expressions are as follows:

$$F_{DN}^i=-{\displaystyle \frac{1}{2}}{\rho }_f{C}_{D}D{V}^2{\sin }^2\theta ,H_{DT}^i={\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B{V}^2{\cos }^2{\theta }_1,H_{DN}^i=-{\displaystyle \frac{1}{2}}{\rho }_f{C}_{DB}{A}_B{V}^2{\sin }^2{\theta }_1$$

(6)

2.2. Modal Approach

A modal approach is then employed in order to discretize the system dynamics depicted in Equation (2). For the coupled system, the displacement field of the cable is expanded as follows:

$$v\left(s,t\right)={\displaystyle \sum _{n=1}^N{\varphi }_n\left(s\right)q_n\left(t\right)}$$

(7)

where ϕ_n(s) is the eigenshape function of the n-th linear mode associated with modal coordinate q_n(t), and N is the total number of modes that kept in the expansion. By taking into account the complex boundary condition, the eigenshape functions together with the eigenfrequencies are computed via

$${\varphi }_n\left(s\right)=C_1\left(\cos {\eta }_{n}s+{\displaystyle \frac{z_1{\eta }_n\cos {\eta }_n+\sin {\eta }_n-z_1{\eta }_n}{\cos {\eta }_n-z_1{\eta }_n\sin {\eta }_n}}\sin {\eta }_{n}s-1\right)$$

(8)

$$\left[-{\lambda }^2-\left(z_1{\lambda }^2+{\displaystyle \frac{EA}{P}}\right){\eta }_n^2+z_1{\eta }_n^4\right]\sin {\eta }_n+\left[\left({\lambda }^2-2z_1{\lambda }^2+{\displaystyle \frac{EA}{z_{1}P}}\right){\eta }_n-{\eta }_n^3\right]\cos {\eta }_n+2z_1{\lambda }^2{\eta }_n=0$$

(9)

where z₁ = M_a/(m_aL) represents the ratio of the buoy mass with added fluid mass to the cable mass with added fluid mass; λ² = κ²EA/P represents the Irvine parameter (cable parameter); ${\eta }_n^2={\omega }_n^2{m}_a{L}^2/P$ is the eigenvalue parameter used for solving the system’s frequencies, where ω_n is the n-th natural frequency of the system; C₁ is a constant that can be determined through modal normalization.

Moreover, owing to the system’s complex boundary conditions, the modal functions are required to satisfy the following orthogonality relations during the Galerkin discretization process:

$${\displaystyle {\int }_0^1{\varphi }_n{\varphi }_k\mathrm{d}s}+{\displaystyle \frac{M_a}{m_{a}L}}{\varphi }_n\left(1\right){\varphi }_k\left(1\right)=\left\{\begin{array}{l}1,n=k\\ 0,n\ne k\end{array}\right.$$

(10)

Therefore, the system of ordinary differential equations for the coupled system can be obtained as follows:

$$\ddot{U}+2{\mu }_1\dot{U}+{\omega }_1^{2}U=\kappa {\omega }_1^2{\displaystyle \sum _{n=1}^N{A}_n}{q}_n-{\displaystyle \frac{{\omega }_1^2}{2}}{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^N{P}_{ij}}}{q}_i{q}_j+{\displaystyle \frac{H_{DT}-H_{DT}^i}{M_{a}L}}$$

(11)

$$\begin{array}{c}{\ddot{q}}_k+2{\mu }_k{\dot{q}}_k+{\omega }_{vk}^2{q}_k={\displaystyle \frac{\kappa EA}{m_a{L}^2}}{A}_{k}U-{\displaystyle \frac{EA}{m_a{L}^2}}{\displaystyle \sum _{n=1}^N{P}_{nk}}U{q}_n+{\displaystyle \frac{\kappa EA}{2m_a{L}^2}}{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^N\left(2P_{ik}{A}_j+P_{ij}{A}_k\right)}}{q}_i{q}_j\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-{\displaystyle \frac{EA}{2m_a{L}^2}}{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^N{\displaystyle \sum _{n=1}^N{P}_{ij}{P}_{nk}}}}{q}_i{q}_j{q}_n+{\displaystyle {\int }_0^1\frac{F_{DN}-F_{DN}^i}{m_{a}L}{\varphi }_k\mathrm{d}s}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\displaystyle \frac{H_{DN}-H_{DN}^i}{m_a{L}^2}}{\varphi }_k\left(1\right)+{\displaystyle \frac{A_k}{m_{a}L}}{f}_v\cos \left(\Omega t\right)\end{array}$$

(12)

where U(t) = u(1,t), and

$$A_n={\displaystyle {\int }_0^1{\varphi }_n\mathrm{d}s},P_{ij}={\displaystyle {\int }_0^1{{\varphi }^{\prime }}_i{{\varphi }^{\prime }}_j\mathrm{d}s},{\mu }_1={\displaystyle \frac{c_u}{2M_a}},{\mu }_k={\displaystyle \frac{1}{2m_a}}{\displaystyle {\int }_0^1{c}_v{\varphi }_k^2\mathrm{d}s},{\omega }_{vk}^2=\left[{\displaystyle \frac{P}{m_a{L}^2}}{P}_{kk}+{\displaystyle \frac{{\kappa }^{2}EA}{m_a{L}^2}}{A}_k^2\right]$$

(13)

It can be observed that the coupled cable–buoy system comprises both nonlinear and linear coupling terms. The presence of linear coupling terms means that ω₁ and ω_vk cannot directly represent the actual vibration frequencies of the system [32]. Consequently, employing high-order perturbation methods for approximate solutions will significantly increase computational complexity.

To enhance the analytical accuracy, the nonlinear terms associated with the buoy are neglected based on the theory of linear vibration. Building upon the coupled system framework, the displacement field of the buoy is expressed as follows:

$$U\left(t\right)={\displaystyle \sum _{n=1}^N{\xi }_n{q}_n\left(t\right)}$$

(14)

where ξ_n represents the amplitude constant of the eigenshape function for the n-th linear mode. Its relationship with the constant C₁ is expressed as ξ_n = C₁m_aL/(κM_a).

Thus, substituting Equation (14) into Equation (12), the ordinary differential equations for the decoupled system can be derived as follows:

$$\begin{array}{c}{\ddot{q}}_k+2{\mu }_k{\dot{q}}_k+{\omega }_k^2{q}_k={\displaystyle \frac{\kappa EA}{2m_a{L}^2}}{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^N\left(2P_{ik}{A}_j+P_{ij}{A}_k-\frac{2}{\kappa }{P}_{ik}{\xi }_j\right)}}{q}_i{q}_j\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-{\displaystyle \frac{EA}{2m_a{L}^2}}{\displaystyle \sum _{i=1}^N{\displaystyle \sum _{j=1}^N{\displaystyle \sum _{n=1}^N{P}_{ij}{P}_{nk}}}}{q}_i{q}_j{q}_n+{\displaystyle {\int }_0^1\frac{F_{DN}-F_{DN}^i}{m_{a}L}{\varphi }_k\mathrm{d}s}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+{\displaystyle \frac{H_{DN}-H_{DN}^i}{m_a{L}^2}}{\varphi }_k\left(1\right)+{\displaystyle \frac{A_k}{m_{a}L}}{f}_v\cos \left(\Omega t\right)\end{array}$$

(15)

where

$${\omega }_k^2=\left[{\displaystyle \frac{P}{m_a{L}^2}}{P}_{kk}+{\displaystyle \frac{{\kappa }^{2}EA}{m_a{L}^2}}{A}_k^2-{\displaystyle \frac{\kappa EA}{m_a{L}^2}}{\xi }_k{A}_k\right]$$

(16)

The coupled and decoupled systems obtained above will next be solved using perturbation methods.

3. Multiple-Time-Scale Analysis

3.1. Coupled System

In the approximate analysis, only the k-th mode is considered, and the relative velocity drag force model is adopted. Under this condition, Equations (11) and (12) can be simplified as follows:

$$\ddot{U}+2{\mu }_1\dot{U}+{\omega }_1^{2}U=c_1{q}_k+c_2{q}_k^2-c_3{Q}_{1,k}\dot{U}$$

(17)

$${\ddot{q}}_k+2{\mu }_k{\dot{q}}_k+{\omega }_{vk}^2{q}_k=d_{1}U+d_{2}U{q}_k+d_3{q}_k^2+d_4{q}_k^3-d_5{Q}_{2,k}{\dot{q}}_k+p_v\cos \left(\Omega t\right)$$

(18)

where the nonlinear drag force is linearized, and

$$\begin{array}{l}{c}_1=\kappa {\omega }_1^2{A}_k,c_2=-{\displaystyle \frac{{\omega }_1^2}{2}}{P}_{kk},c_3={\displaystyle \frac{{\rho }_f{C}_{DB}{A}_{B}L}{2M_a}},Q_{1,k}={\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{U}}^2\left|\dot{U}\right|\mathrm{d}t}/{\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{U}}^2\mathrm{d}t}\\ d_1={\displaystyle \frac{\kappa EA}{m_a{L}^2}}{A}_k,d_2=-{\displaystyle \frac{EA}{m_a{L}^2}}{P}_{kk},d_3={\displaystyle \frac{3\kappa EA}{2m_a{L}^2}}{P}_{kk}{A}_k,d_4=-{\displaystyle \frac{EA}{2m_a{L}^2}}{P}_{kk}^2,p_v={\displaystyle \frac{A_k}{m_{a}L}}{f}_v\\ d_5=\left[{\displaystyle \frac{{\rho }_f{C}_{D}DL}{2m_a}}{\displaystyle {\int }_0^1{\varphi }_k^2\left|{\varphi }_k\right|\mathrm{d}s}+{\displaystyle \frac{{\rho }_f{C}_{DB}{A}_B}{2m_a}}{\varphi }_k^2\left(1\right)\left|{\varphi }_k\left(1\right)\right|\right],Q_{2,k}={\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{q}}_k^2\left|{\dot{q}}_k\right|\mathrm{d}t}/{\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{q}}_k^2\mathrm{d}t}\end{array}$$

(19)

Following the multiple scales method, the non-dimensional parameter ε is introduced for bookkeeping purposes, and the solution of Equations (17) and (18) is assumed as the following:

$$U\left(t,\epsilon \right)={\displaystyle \sum _{m=0}^2{\epsilon }^m}{U}_m\left(T_0,T_1,T_2\right),q_k\left(t,\epsilon \right)={\displaystyle \sum _{m=0}^2{\epsilon }^m}{q}_{km}\left(T_0,T_1,T_2\right)$$

(20)

where T_m = ε^mt (m = 0, 1, 2) represents the multiple time scales.

In order to balance the damping, excitation, and nonlinear terms, μ₁, Q_1,k, μ_vk, Q_2,k and p_v are replaced with ε²μ₁, ε²Q_1,k, ε²μ_vk, ε²Q_1,k and ε³p_v, respectively. Substituting Equation (20) into Equations (17) and (18), equating the terms of like order in ε, yields the following:

ε⁰:

$$\begin{array}{c}{D}_0^2{U}_0+{\omega }_1^2{U}_0=c_1{q}_{k0}\\ D_0^2{q}_{k0}+{\omega }_{vk}^2{q}_{k0}=d_1{U}_0\end{array}$$

(21)

ε¹:

$$\begin{array}{c}{D}_0^2{U}_1+{\omega }_1^2{U}_1=-2D_0{D}_1{U}_1+c_1{q}_{k1}+c_2{q}_{k0}^2\\ D_0^2{q}_{k1}+{\omega }_{vk}^2{q}_{k1}=-2D_0{D}_1{q}_{k0}+d_1{U}_1+d_2{U}_0{q}_{k0}+d_3{q}_{k0}^2\end{array}$$

(22)

ε²:

$$\begin{array}{c}{D}_0^2{U}_2+{\omega }_1^2{U}_2=-2D_0{D}_1{U}_1-2D_0{D}_2{U}_0-D_1^2{U}_0\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-2{\mu }_1{D}_0{U}_0-c_3{Q}_{1,k}{D}_0{U}_0+c_1{q}_{k2}+2c_2{q}_{k0}{q}_{k1}\hfill \\ D_0^2{q}_{k2}+{\omega }_{vk}^2{q}_{k2}=-2D_0{D}_1{q}_{k1}-2D_0{D}_2{q}_{k0}-D_1^2{q}_{k0}-2{\mu }_k{D}_0{q}_{k0}-d_5{Q}_{2,k}{D}_0{q}_{k0}+d_1{U}_2\hfill \\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+d_2{U}_1{q}_{k0}+d_2{U}_0{q}_{k1}+2d_3{q}_{k0}{q}_{k1}+d_4{q}_{k0}^3+p_v\cos \left(\Omega t\right)\hfill \end{array}$$

(23)

where D_k = ∂/∂T_k and ${\mathrm{D}}_{\mathrm{k}}^2$ = ∂²/∂${\mathrm{T}}_{\mathrm{k}}^2$ are differential operators.

The solutions to Equation (21) can be expressed as follows:

$$U_0=B_k\left(T_1,T_2\right){\mathrm{e}}^{\mathrm{i}{\omega }_k{T}_0}+cc,q_{k0}=A_k\left(T_1,T_2\right){\mathrm{e}}^{\mathrm{i}{\omega }_k{T}_0}+cc$$

(24)

where B_k and A_k are complex-valued amplitudes of the modes, which will be determined by imposing the solvability conditions. cc denotes the complex conjugate of the preceding terms on the right-hand side. By substituting Equation (24) into Equation (22), the relationship between B_k and A_k can be obtained:

$$B_k={\displaystyle \frac{c_1}{{\omega }_1^2-{\omega }_k^2}}{A}_k\left(\mathrm{or}{B}_k={\displaystyle \frac{{\omega }_{vk}^2-{\omega }_k^2}{d_1}}{A}_k\right)$$

(25)

As time goes on, the homogeneous solutions of Equation (22) will disappear, leaving only the particular solutions. Substituting Equation (24) into Equation (22) yields

$$\begin{array}{c}{U}_1={\displaystyle \frac{c_2}{{\omega }_1^2}}{A}_k{\overline{A}}_k-{\displaystyle \frac{c_2}{4{\omega }_k^2-{\omega }_1^2}}{A}_k^2{\mathrm{e}}^{2\mathrm{i}{T}_0{\omega }_k}+cc\\ q_{k1}={\displaystyle \frac{d_3}{{\omega }_{vk}^2}}{A}_k{\overline{A}}_k+{\displaystyle \frac{d_2}{{\omega }_{vk}^2}}{A}_k{\overline{B}}_k-{\displaystyle \frac{d_2}{4{\omega }_k^2-{\omega }_{vk}^2}}{A}_k{B}_k{\mathrm{e}}^{2\mathrm{i}{T}_0{\omega }_k}-{\displaystyle \frac{d_3}{4{\omega }_k^2-{\omega }_{vk}^2}}{A}_k^2{\mathrm{e}}^{2\mathrm{i}{T}_0{\omega }_k}+cc\end{array}$$

(26)

and eliminating the secular terms yields

$$D_1{A}_k=D_1{\overline{A}}_k=D_1{B}_k=D_1{\overline{B}}_k=0$$

(27)

The detuning parameter σ₁ is introduced to quantify the degree of proximity between the external excitation frequency and the system’s natural frequency, facilitating the investigation of the detuning effect on the system’s dynamic behavior. It is defined as Ω = ω_k + ε²σ₁. Substituting Equations (22) and (26) into Equation (23), and introducing Equation (27), by eliminating the secular terms, one obtains

$$\begin{array}{c}\mathrm{i}{\omega }_k\left({\displaystyle \frac{1}{2}}{c}_3{Q}_{1,k}{B}_k+{\mu }_1{B}_k+D_2{B}_k\right)={\displaystyle \frac{{\Delta }_4}{2}}{A}_k^2{\overline{A}}_k+{\displaystyle \frac{{\Delta }_5}{2}}{A}_k{B}_k{\overline{A}}_k+{\displaystyle \frac{c_2{d}_2}{{\omega }_{vk}^2}}{A}_k^2{\overline{B}}_k\\ \mathrm{i}{\omega }_k\left({\displaystyle \frac{1}{2}}{d}_5{Q}_{2,k}{A}_k+{\mu }_k{A}_k+D_2{A}_k\right)={\Delta }_1{A}_k^2{\overline{A}}_k+{\Delta }_2{A}_k{B}_k{\overline{A}}_k+{\Delta }_3{A}_k^2{\overline{B}}_k+{\displaystyle \frac{d_2^2}{{\omega }_{vk}^2}}{B}_k^2{\overline{A}}_k\\ +{\displaystyle \frac{d_2^2}{{\omega }_{vk}^2}}{A}_k{B}_k{\overline{B}}_k-{\displaystyle \frac{d_2^2}{4{\omega }_k^2-{\omega }_{vk}^2}}{A}_k{B}_k{\overline{B}}_k+{\displaystyle \frac{1}{2}}{p}_v{\mathrm{e}}^{\mathrm{i}{\sigma }_1{T}_2}\end{array}$$

(28)

where

$$\begin{array}{l}{\Delta }_1=-{\displaystyle \frac{c_2{d}_2}{4{\omega }_k^2-{\omega }_1^2}}+{\displaystyle \frac{2c_2{d}_2}{{\omega }_1^2}}-{\displaystyle \frac{2d_3^2}{4{\omega }_k^2-{\omega }_{vk}^2}}+{\displaystyle \frac{4d_3^2}{{\omega }_{vk}^2}}+3d_4,{\Delta }_2={\displaystyle \frac{4d_2{d}_3}{{\omega }_{vk}^2}}-{\displaystyle \frac{2d_2{d}_3}{4{\omega }_k^2-{\omega }_{vk}^2}}\\ {\Delta }_3={\displaystyle \frac{2d_2{d}_3}{{\omega }_{vk}^2}}-{\displaystyle \frac{d_2{d}_3}{4{\omega }_k^2-{\omega }_{vk}^2}},{\Delta }_4={\displaystyle \frac{4c_2{d}_3}{{\omega }_{vk}^2}}-{\displaystyle \frac{2c_2{d}_3}{4{\omega }_k^2-{\omega }_{vk}^2}},{\Delta }_5={\displaystyle \frac{2c_2{d}_2}{{\omega }_{vk}^2}}-{\displaystyle \frac{2c_2{d}_2}{4{\omega }_k^2-{\omega }_{vk}^2}}\end{array}$$

(29)

The polar form of A_k(T₁,T₂) can be expressed as follows:

$$A_k\left(T_1,T_2\right)={\displaystyle \frac{1}{2}}a\left(T_1,T_2\right){\mathrm{e}}^{\mathrm{i}\beta \left(T_1,T_2\right)}$$

(30)

where a_k and β_k are the amplitude and phase angle of A_k, respectively. By substituting Equations (25) and (30) into Equation (28) and separating the real and imaginary parts, the autonomous modulation equations in the polar form can be obtained as follows:

$$\dot{a}=-{\displaystyle \frac{1}{2}}{d}_5{Q}_{2,k}a-{\mu }_{k}a+{\displaystyle \frac{p_v}{2{\omega }_k}}\sin \gamma ,a\dot{\gamma }={\Delta }_6{a}^3+{\sigma }_{1}a+{\displaystyle \frac{p_v}{2{\omega }_k}}\cos \gamma$$

(31)

where γ = σ₁t−β and

$${\Delta }_6=-{\displaystyle \frac{{\Delta }_1}{8{\omega }_k}}-{\displaystyle \frac{c_1\left({\Delta }_2+{\Delta }_3\right)}{8{\omega }_k\left({\omega }_1^2-{\omega }_k^2\right)}}-{\displaystyle \frac{c_1^2{d}_2^2}{4{\omega }_k{\left({\omega }_1^2-{\omega }_k^2\right)}^2{\omega }_{vk}^2}}+{\displaystyle \frac{c_1^2{d}_2^2}{8{\omega }_k{\left({\omega }_1^2-{\omega }_k^2\right)}^2\left(4{\omega }_k^2-{\omega }_{vk}^2\right)}}$$

(32)

The steady-state solutions can be found by sets $\dot{\mathrm{a}}$ = 0 and $\dot{\mathsf{\gamma }}$ = 0. Solving the Equation (31), the frequency response equation is obtained:

$${\left({\displaystyle \frac{1}{2}}{d}_5{Q}_{2,k}a+{\mu }_{k}a\right)}^2+{\left({\Delta }_6{a}^3+{\sigma }_{1}a\right)}^2={\left({\displaystyle \frac{p_v}{2{\omega }_k}}\right)}^2$$

(33)

Furthermore, small perturbations are introduced as a = a₀ + δa and γ = γ₀ + δγ to analyze the stability of the steady-state solution. By substituting these expressions into Equation (31), the stability is determined by the eigenvalues of the resulting coefficient matrix. The steady-state solution is considered stable if and only if all eigenvalues have real parts that are less than or equal to zero.

3.2. Decoupled System

Similar to the treatment of the coupled system, when applying the relative velocity drag force model and considering the k-th mode, the nonlinear drag force is linearized. Under these conditions, Equation (15) can be simplified as follows:

$${\ddot{q}}_k+2{\mu }_k{\dot{q}}_k+{\omega }_k^2{q}_k=c_1{q}_k^2+c_2{q}_k^3-c_3{Q}_k{\dot{q}}_k+p_v\cos \left(\Omega t\right)$$

(34)

where

$$\begin{array}{c}{c}_1={\displaystyle \frac{\kappa EA}{m_a{L}^2}}\left(3A_k{P}_{kk}-{\displaystyle \frac{2}{\kappa }}{\xi }_k{P}_{kk}\right),c_2=-{\displaystyle \frac{EA}{2m_a{L}^2}}{P}_{kk}^2,Q_k={\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{q}}_k^2\left|{\dot{q}}_k\right|\mathrm{d}t}/{\displaystyle {\int }_0^{\frac{2\pi }{{\omega }_k}}{\dot{q}}_k^2\mathrm{d}t}\\ c_3=\left[{\displaystyle \frac{{\rho }_f{C}_{D}DL}{2m_a}}{\displaystyle {\int }_0^1{\varphi }_k^2\left|{\varphi }_k\right|\mathrm{d}s}+{\displaystyle \frac{{\rho }_f{C}_{DB}{A}_B}{2m_a}}{\varphi }_k^2\left(1\right)\left|{\varphi }_k\left(1\right)\right|\right]\end{array}$$

(35)

Since there is no coupling in this system, the accuracy of the approximate analysis can be improved. The solution to Equation (34) is given by

$$q_k\left(t,\epsilon \right)={\displaystyle \sum _{m=0}^4{\epsilon }^m}{q}_{km}\left(T_0,T_1,T_2,T_3,T_4\right)$$

(36)

To balance the damping, excitation, and nonlinear terms, μ_k, Q_k, and p_v are replaced by ε²μ_k, ε²Q_k, and ε³p_v, respectively. Substituting Equation (36) into Equation (34) and equating the terms of like order in ε yields

ε⁰:

$$D_0^2{q}_{k0}+{\omega }_k^2{q}_{k0}=0$$

(37)

ε¹:

$$D_0^2{q}_{k1}+{\omega }_k^2{q}_{k1}=-2D_0{D}_1{q}_{k0}+c_1{q}_{k0}^2$$

(38)

ε²:

$$\begin{array}{c}{D}_0^2{q}_{k2}+{\omega }_k^2{q}_{k2}=-2D_0{D}_1{q}_{k1}-2D_0{D}_2{q}_{k0}-D_1^2{q}_{k0}+2c_1{q}_{k0}{q}_{k1}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+c_2{q}_{k0}^3-2{\mu }_v{D}_0{q}_{k0}-c_3{Q}_k{D}_0{q}_{k0}+p_v\cos (\Omega t)\end{array}$$

(39)

ε³:

$$\begin{array}{c}{D}_0^2{q}_{k3}+{\omega }_k^2{q}_{k3}=-2D_0{D}_1{q}_{k2}-2D_0{D}_2{q}_{k1}-D_1^2{q}_{k1}-2D_0{D}_3{q}_{k0}-2D_1{D}_2{q}_{k0}+3c_2{q}_{k0}^2{q}_{k1}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}+c_1\left(q_{k1}^2+2q_{k0}{q}_{k2}\right)-2{\mu }_k\left(D_0{q}_{k1}+D_1{q}_{k0}\right)-c_3{Q}_k\left(D_0{q}_{k1}+D_1{q}_{k0}\right)\end{array}$$

(40)

ε⁴:

$$\begin{array}{c}{D}_0^2{q}_{k4}+{\omega }_k^2{q}_{k4}=-2D_0{D}_1{q}_{k3}-2D_0{D}_2{q}_{k2}-2D_0{D}_3{q}_{k1}-2D_0{D}_4{q}_{k0}-D_1^2{q}_{k2}-2D_1{D}_2{q}_{k1}\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-2D_1{D}_3{q}_{k0}-D_2^2{q}_{k0}+2c_1\left(q_{k0}{q}_{k3}+q_{k1}{q}_{k2}\right)+3c_2\left(q_{k0}^2{q}_{k2}+q_{k0}{q}_{k1}^2\right)\\ \hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}\hspace{1em}-\left(2{\mu }_k+c_3{Q}_k\right)\left(D_0{q}_{k2}+D_1{q}_{k1}+D_2{q}_{k0}\right)\end{array}$$

(41)

The solutions to Equation (37) can be expressed as follows:

$$q_{k0}=A\left(T_1,T_2,T_3,T_4\right){\mathrm{e}}^{\mathrm{i}{\omega }_k{T}_0}+cc$$

(42)

The subsequent solution process follows a similar approach to that presented in Section 3.1. The detuning parameters σ₁ and σ₂ are defined such that Ω = ω_k + ε²σ₁ + ε⁴σ₂. Equations (38)–(41) are solved sequentially, and the secular terms are eliminated, yielding

$$D_{1}A=D_1\overline{A}=D_{3}A=D_3\overline{A}=0$$

(43)

$$\begin{array}{c}{D}_{2}A=-{\displaystyle \frac{1}{2}}{c}_3{Q}_{k}A-{\mu }_{k}A-\left({\displaystyle \frac{5c_1^2}{3{\omega }_k^3}}+{\displaystyle \frac{3c_2}{2{\omega }_k}}\right)A^2\overline{A}\mathrm{i}-{\displaystyle \frac{p_v}{4{\omega }_k}}{\mathrm{ie}}^{\mathrm{i}{T}_2{\sigma }_1+\mathrm{i}{T}_4{\sigma }_2}\\ D_{4}A={\displaystyle \frac{2c_3{c}_1^2{Q}_k}{9{\omega }_k^4}}{A}^2\overline{A}+{\displaystyle \frac{4c_1^2{\mu }_k}{9{\omega }_k^4}}{A}^2\overline{A}+\left({\displaystyle \frac{3c_2^2}{16{\omega }_k^3}}-{\displaystyle \frac{335c_1^4}{108{\omega }_k^7}}-{\displaystyle \frac{143c_2{c}_1^2}{12{\omega }_k^5}}\right)A^3{\overline{A}}^2\mathrm{i}\\ \hspace{1em}\hspace{1em}+{\displaystyle \frac{c_3{Q}_k}{2{\omega }_k}}\mathrm{i}{D}_{2}A+{\displaystyle \frac{{\mu }_k}{{\omega }_k}}\mathrm{i}{D}_{2}A+{\displaystyle \frac{1}{2{\omega }_k}}\mathrm{i}{D}_2^{2}A+{\displaystyle \frac{2c_1^2{p}_v}{9{\omega }_k^5}}\mathrm{i}A\overline{A}{\mathrm{e}}^{\mathrm{i}{\sigma }_1{T}_2+\mathrm{i}{\sigma }_2{T}_4}\end{array}$$

(44)

The polar form of A(T₁,T₂,T₃,T₄) can be expressed as follows:

$$A\left(T_1,T_2,T_3,T_4\right)={\displaystyle \frac{1}{2}}a\left(T_1,T_2\right){\mathrm{e}}^{\mathrm{i}\beta \left(T_1,T_2\right)}$$

(45)

where a and β are the amplitude and phase angle of A, respectively. We substitute Equation (45) into Equation (44), and apply the method of reorganization dA(t,ε)/dt = ε²D₂A + ε⁴D₄A. By separating the real and imaginary parts, the modulation equation is derived, from which the corresponding frequency response equation can be obtained:

$$\begin{array}{l}{\left({\displaystyle \frac{1}{2}}{c}_3{Q}_{k}a+{\mu }_{k}a\right)}^2+{\left({\sigma }_{1}a+{\Delta }_1{a}^3\right)}^2={\left({\displaystyle \frac{p_v}{2{\omega }_k}}\right)}^2\\ {\left({\displaystyle \frac{p_v{c}_3}{4{\omega }_k^2}}{Q}_k+{\displaystyle \frac{p_v{\mu }_k}{2{\omega }_k^2}}\right)}^2+{\left({\displaystyle \frac{p_v{c}_1^2}{9{\omega }_k^5}}{a}^2\right)}^2={\left({\Delta }_2{a}^3+{\Delta }_3{Q}_k{a}^3\right)}^2+{\left({\Delta }_4{a}^5+{\sigma }_{2}a+{\displaystyle \frac{{\sigma }_1^2}{2{\omega }_k}}a+{\displaystyle \frac{{\mu }_k^2}{{\omega }_k}}a+{\displaystyle \frac{c_3{\mu }_k}{{\omega }_k}}{Q}_{k}a+{\displaystyle \frac{c_3^2}{4{\omega }_k}}{Q}_k^{2}a\right)}^2\end{array}$$

(46)

where

$${\Delta }_1={\displaystyle \frac{5c_1^2}{12{\omega }_k^3}}+{\displaystyle \frac{3c_2}{8{\omega }_k}},{\Delta }_2={\displaystyle \frac{c_1^2{\mu }_{vk}}{9{\omega }_k^4}}+{\displaystyle \frac{{\mu }_{vk}{\Delta }_1}{{\omega }_k}},{\Delta }_3={\displaystyle \frac{c_1^2{c}_3}{18{\omega }_k^4}}+{\displaystyle \frac{c_3{\Delta }_1}{2{\omega }_k}},{\Delta }_4={\displaystyle \frac{335c_1^4}{1728{\omega }_k^7}}+{\displaystyle \frac{143c_1^2{c}_2}{192{\omega }_k^5}}-{\displaystyle \frac{3c_2^2}{256{\omega }_k^3}}$$

(47)

Finally, the stability analysis of the solution and the treatment of the coupled system are approached in the same manner. The next section presents the numerical study analysis.

4. Numerical Results and Discussion

In this section, the dynamic behaviors of the cable–buoy system are analyzed and discussed, utilizing the governing equations of different models and their theoretical solutions. A summary of the key parameters is provided in

Table 1. Key parameters of the cable–buoy system.

Object	Parameters	Value
Fluid	Velocity V	1 m/s
	Density ρf	1025 kg/m3
Cable	Length L	40 m
	Diameter D	0.0155 m
	Density ρc	4104.52 kg/m3
	Added-mass coefficient Cac	1.0
	Drag coefficient CD	1.05
Buoy	Mass M	170 kg
	Diameter DM	0.8 m
	Added-mass coefficient CaB	0.5
	Drag coefficient CDB	0.4

[33].

Under the above parameter conditions, the buoy mass parameter is calculated as z₁ ≈ 8. The relationship curve between the modal frequencies ω₂/ω_s and the Irvine parameter λ/π for the first six modes of the cable–buoy system in the streamwise direction is shown in Figure 2. Here, the modal frequencies ω₂ are non-dimensionalized by the first-order frequency ω_s = 2.7164 rad/s, which corresponds to the fixed-end condition of the anchor cable. The points labeled A, B, and C indicate the occurrence of veering phenomena between modal frequencies for specific Irvine parameters.

The veering point is a characteristic point where the modal frequency curve changes direction, with two modal frequencies approaching each other before rapidly diverging. At this point, there may be strong modal coupling and complex modal interactions between the modes. Therefore, in the subsequent analysis, the Irvine parameter λ = 1.5π is selected, corresponding to point B in the figure, while the Young’s modulus of the cable is E = 7.6 × 10⁸ Pa. Meanwhile, the numerical solutions in subsequent steps primarily adopt the fourth-order Runge–Kutta method due to its simplicity, efficiency, and suitability for systems with moderate to low stiffness and small time steps.

4.1. Comparison of Drag Force Models

To compare the effects of drag forces on the dynamic response of the cable–buoy system under absolute and relative velocity models, the fourth-order Runge–Kutta method is used to directly solve Equations (11) and (12) numerically. The displacement–time history, velocity–time history, and phase plots for the single-mode cable–buoy coupled system under different drag force models are illustrated in Figure 3. The black solid line represents the absolute velocity drag force model, while the red solid line represents the relative velocity drag force model.

It can be observed that, under the relative velocity drag force model, both the longitudinal displacement of the buoy and the streamwise displacement of the cable approach equilibrium near zero. In contrast, under the absolute velocity drag force model, the equilibrium positions of both the buoy and the cable exhibit significant offsets. This discrepancy arises because, in the relative velocity drag force model, the motion velocity is referenced to the flow velocity, ensuring that the equilibrium position remains consistent under steady flow conditions. In the absolute velocity drag force model, however, the system’s motion velocity is referenced to a fixed ground, with the equilibrium position thus referenced to a stationary flow velocity (i.e., zero flow velocity), leading to the observed displacement shifts.

In addition, the absolute velocity drag force model includes a coupling term related to the flow velocity, whereas the relative velocity model does not. As a result, the system under the absolute velocity model experiences greater nonlinear damping forces, leading to a slightly lower amplitude of dynamic response compared to the relative velocity model. However, the omission of the coupling term has a limited impact on the overall system’s dynamic response, and the general trend remains consistent. From the phase plot, it can be observed that the phase trajectories in both plots form closed elliptical shapes, indicating that the motion of both the buoy and the cable is periodic.

As shown in the frequency spectrum in Figure 4, the overall dynamic response of the system exhibits some differences under the two drag force models. However, the relative velocity drag force, as a simplified model, amplifies the dynamic response trend. Therefore, this error is acceptable in the single-mode analysis.

To further analyze the system, the dynamic response comparison of the two drag force models for the fourth-order mode is shown in Figure 5, with excitation applied only to the higher-order modes. It can be observed that in the comparison of the two drag force models, only the displacement response of the lower-order modes of the cable shows a shift in the equilibrium position, while the higher-order modes do not exhibit such a shift. This is due to the strong coupling effects at the veering point. Additionally, the lower-order modes and the buoy exhibit periodic oscillations, while the higher-order modes do not. This further confirms that under multi-modal coupling effects, energy transfer occurs between coupled modes, leading to the complexity of the overall system’s dynamic characteristics.

Furthermore, the damping effects caused by the two drag force models, amplified by the modal coupling effects, result in a greater damping effect in the absolute velocity drag force model. This leads to significant differences in the dynamic response of the cable’s lower-order modes between the two models. However, the dynamic responses of the buoy and cable’s lower-order modes are one to two orders of magnitude smaller than those of the higher-order modes, meaning that the errors in the dynamic response of the higher-order modes can be considered as errors in the overall system response. From the perspective of the cable’s higher-order modal dynamic response, these errors are within a reasonable range, indicating that the relative velocity drag force model provides acceptable accuracy as a simplified approach.

In summary, in the relative velocity drag force model, the system’s velocity is calculated relative to the fluid velocity, so the equilibrium position can be directly treated as the reference position. In the absolute velocity drag force model, the system’s velocity is calculated relative to the ground, which results in the equilibrium position being referenced at zero fluid velocity, causing the displacement response to shift. The multi-modal analysis reveals that the differences in the system’s dynamic response under the two drag force models mainly stem from the modal coupling effects at the veering point. As shown in the frequency spectrum in Figure 6, the relative velocity drag force model similarly amplifies the dynamic response of the cable in the system, while the dynamic response of the buoy exhibits the opposite trend. Although the relative velocity drag force model generally amplifies the cable’s dynamic response, the goal in studying dynamic characteristics is to optimize structural design and reduce response amplitudes. Therefore, under this simplified model, amplifying the system’s dynamic response is preferable to reducing it, as the introduced error is deemed acceptable for practical applications. Thus, when studying the dynamic behavior of the cable, the relative velocity drag force model serves as a reasonable simplification of the absolute velocity drag force model.

4.2. Comparison of Coupled and Decoupled Systems

In this section, to compare the dynamic responses of the coupled and decoupled systems, the fourth-order Runge–Kutta method is used to directly solve the coupled system Equations (11) and (12) as well as the decoupled system (15). Similar to the comparison with drag forces, Figure 7 presents the displacement–time history curve, velocity–time history curve, and phase plot for the single-mode coupled and decoupled cable–buoy systems. The black solid line represents the coupled system, and the red solid line represents the decoupled system.

The results show that the longitudinal displacement and velocity response amplitudes of the buoy in the decoupled system are smaller than those in the coupled system. This is because, in the decoupled system, the nonlinear motion of the buoy is neglected, and only its linear motion characteristics are considered, leading to a reduction in the dynamic response amplitude. However, the neglect of the buoy’s nonlinear motion only has a limited impact on the response amplitude, and the overall error is within an acceptable range.

Furthermore, the displacement and velocity responses of the cable in the decoupled system are nearly identical to those in the coupled system, indicating that the nonlinear motion of the buoy has a minimal impact on the cable’s dynamic response. In the phase plot, both the buoy and cable exhibit closed elliptical trajectories, suggesting that both display typical periodic motion. As shown in the frequency spectrum in Figure 8, the overall dynamic response of the coupled and decoupled systems shows little difference. This indicates that, when focusing on the dynamic response of the cable, appropriate decoupling of the system can effectively simplify the analysis process while maintaining the accuracy and reliability of the results.

Similarly, to further analyze the system, Figure 9 presents a comparison of the dynamic responses of the coupled and decoupled cable–buoy systems under multi-mode analysis, considering the fourth-order mode, with excitation applied only to the higher-order modes. As shown, the dynamic response amplitude of the buoy in the decoupled system is smaller than that in the coupled system, which is consistent with the single-mode analysis. However, the difference lies in the periodic oscillations generated by the buoy and cable in the lower-order modes, reflecting the complexity of the system’s dynamic behavior at veering points. In contrast, the dynamic response of the cable in higher-order modes is nearly identical in both systems, exhibiting periodic motion.

Further analysis reveals that when excitation is applied to the higher-order modes of the cable, the nonlinear damping effect of the underwater drag force leads to gradual energy dissipation during the transfer process, particularly as energy is transmitted to the buoy and the lower-order modes of the cable, causing significant amplitude attenuation. This energy dissipation effectively suppresses the dynamic response of the lower-order modes. However, since the magnitude of the lower-order mode response of the cable is relatively small, the error caused by the nonlinear motion of the buoy is amplified in the decoupled system, resulting in discrepancies between the coupled and decoupled systems in the lower-order mode responses. As seen in the frequency spectrum in Figure 10, the impact of neglecting the buoy’s nonlinear terms on the overall response of the cable is limited, due to the relatively small magnitude of the lower-order mode responses, but it significantly affects the buoy’s response.

In conclusion, although the decoupled system neglects the nonlinear motion of the buoy, it still effectively captures the dynamic response of the cable and the overall periodic characteristics of the system, maintaining similar dynamic behavior to the coupled system. Therefore, in more complex nonlinear dynamic analyses, decoupling can effectively simplify computational complexity, demonstrating its applicability.

4.3. Primary Resonance Response

In this section, the nonlinear dynamic behavior of the cable in the cable–buoy system is investigated. The external excitation amplitude is assumed to be f_v = 2 N/m, and the third-order mode of the system is selected. The frequency response curve equations are solved numerically, with the coupled system described by Equation (33) and the decoupled system by Equation (46). The basic parameters used in the analysis are listed in Table 1. However, for the subsequent analysis, the buoy mass M and the Young’s modulus of the cable E need to be determined based on the buoy mass parameter z₁ and the Irvine parameter λ. Additionally, the fluid velocity of V = 0.1 m/s is assumed to ensure that the transverse motion induced by vortex-induced vibrations is neglected under low-flow conditions.

Figure 11 presents a comparison of the frequency response curves for the coupled and decoupled cable–buoy systems under six different parameter sets. Here, the black lines represent the coupled system, whereas the red lines denote the decoupled system. And in the following figures, stable solutions are shown with solid lines, and unstable solutions are represented by dashed lines. To validate the theoretical results, numerical solutions are obtained by directly integrating Equation (15) with the fourth-order Runge–Kutta method, and the results are shown in the figures as circles.

In Figure 11a–c, the frequency response curves of the coupled and decoupled systems exhibit negligible differences for a buoy mass parameter z₁ = 10, across the three Irvine parameters λ = 3π, 6π, and 9π. This observation indicates that, for these parameter combinations, low-order perturbation analysis accurately captures the nonlinear dynamic characteristics of the cable, yielding results consistent with numerical simulations. However, as the buoy mass parameter increases to z₁ = 100, as illustrated in Figure 11d–f, significant differences emerge between the frequency response curves of the coupled and decoupled systems, despite the Irvine parameters λ remaining unchanged. This observation highlights the limitations of low-order perturbation analysis under these conditions, indicating its reduced accuracy in capturing the nonlinear dynamic behavior of the system for higher buoy mass parameters.

The above results indicate that, under specific parameter conditions, analyzing the frequency response curves using the decoupled system achieves higher precision with controllable computational complexity. Moreover, the decoupled system not only simplifies the computational process but also maintains the accuracy of the frequency response analysis in most cases. Therefore, for subsequent analysis, the frequency response equations derived from the decoupled system will be utilized.

To further investigate the impact of different parameters on the nonlinear behavior of the cable–buoy system, Figure 12 illustrates the effect of variations in the Irvine parameter λ on the system’s frequency response curves for three different buoy mass parameters, z₁ = 1, 10, and 100. In Figure 12a,b, it can be observed that for z₁ = 1 and z₁ = 10, the frequency response curves for the six different Irvine parameters λ exhibit similar trends. As λ increases, the system’s nonlinear characteristics gradually intensify, with the frequency response curve shifting to the right and displaying pronounced hard spring behavior. In these cases, the system exhibits a maximum of three steady-state solutions, two of which are stable and one unstable.

As the buoy mass parameter increases to z₁ = 100, the system’s nonlinear behavior becomes more pronounced, as shown in Figure 12c,d. As λ increases from 1.5π to 6π, Figure 12c reveals that the system’s nonlinear characteristics become significantly stronger, with an increased rightward shift of the frequency response curve, continuing to display hard spring behavior. In this case, the system also has three steady-state solutions, two stable and one unstable. However, Figure 12d shows a more complex nonlinear behavior as λ increases from 7.5π to 9π. Initially, the response curve shifts to the left and then to the right, exhibiting soft spring behavior followed by hard spring behavior. In this scenario, the system can have up to five steady-state solutions, three stable and two unstable. As λ increases further to 10.5π, the frequency response curve shifts only to the left, displaying clear soft spring behavior.

Therefore, for a given buoy mass parameter z₁, variations in the cable parameter λ not only significantly affect the nonlinear vibrational behavior of the system but also cause both qualitative and quantitative changes in the system’s dynamic characteristics.

Finally, the effect of variations in the buoy mass parameters z₁ on the nonlinear behavior of the cable–buoy system is analyzed, as illustrated in Figure 13. The results reveal that as z₁ decreases, the system’s nonlinear behavior becomes less complex. Specifically, the system transitions from a softening spring response to a hardening spring response and ultimately exhibits exclusively hardening spring characteristics. This hardening phenomenon implies that the system’s stiffness increases with amplitude, resulting in three steady-state solutions: two stable and one unstable. Such multistability is a hallmark of nonlinear systems, highlighting the intricate dynamic response modes the system can exhibit under certain parameter conditions.

These results indicate that, for a given buoy mass parameter z₁, selecting an appropriate cable parameter λ can substantially influence the system’s nonlinear response characteristics. By tuning these parameters, the resonance response of the system can be effectively mitigated, thereby improving its overall stability. This insight offers valuable theoretical guidance and practical significance for the design and control of cable–buoy systems.

In summary, variations in the buoy mass parameter z₁ and the cable parameter λ not only significantly influence the nonlinear vibration behavior of the system but also lead to qualitative and quantitative changes in its dynamic characteristics. These complex nonlinear dynamic phenomena indicate that, under specific parameter conditions, the system response may transition between hard spring and soft spring behaviors, thereby impacting the system’s stability and response characteristics.

5. Conclusions

The dynamic response of a cable–buoy system in marine engineering was analyzed with a focus on reference frame selection and buoy nonlinearities. The results show that the relative velocity drag force model amplifies the system’s response but remains a valid approximation under typical conditions. Neglecting buoy nonlinearities in the decoupled system leads to minimal variation in cable dynamics but significantly reduces buoy response, indicating that buoy nonlinearities can be ignored when focusing on the cable. Theoretical analysis using the multiple scales method demonstrates higher accuracy for the decoupled system under specific conditions, meeting high-precision requirements. Nonlinear vibration characteristics are highly sensitive to parameter variations, providing insights for structural optimization. Future studies can validate the conclusions of the simplified model through experimental verification. Additionally, further theoretical investigations could explore the system’s bifurcation behavior and potential chaotic characteristics.