Vortex-Induced Vibrations of A Free-Spanning Pipe Based on A Nonlinear Hysteretic Soil Model at the Shoulders

The pipe-soil interactions at shoulders can significantly affect the vortex-induced vibrations (VIV) of free-spanning pipes in the subsea. In this paper, the seabed soil reacting force on the pipe is directly calculated with a nonlinear hysteretic soil model. For the VIV in the middle span, a classic van der Pol wake oscillator is adopted. Based on the Euler-Bernoulli beam theory, the vibration equations of the pipe are obtained which are different in the middle span and at the two end shoulders. The static configuration of the pipe is firstly calculated and then the VIV is simulated. The present model is validated with the comparisons of VIV experiment, pipe-soil interaction experiment and the simulation results of VIV of free-spanning pipes in which the seabed soil is modelled with spring-dashpots. With the present model, the influence of seabed soil on the VIV of a free-spanning pipe is analyzed. The parametric studies show that when the seabed soil has a larger suction area, the pipe vibrates with smaller bending stresses and is safer. While with the increase of the shear strength of the seabed soil, the bending stresses increase and the pipe faces more danger.


Introduction
The submarine pipes are widely used in the ocean engineering to transport oil and natural gas. Because of uneven seabed, the pipes may become free spanning at some places (Leckie et al., 2015). Subjected to external currents, a freespanning pipe will experience vortex-induced vibrations (VIV). The seabed soil at the shoulders can restrict the vibrations of the pipe and may have significant influence on it (Larsen et al., 2004;Vedeld et al., 2016). The influence from the seabed soil should be carefully estimated when predicting the vibrations of free-spanning pipes.
The VIV of free-spanning pipes involve the fluid−structure interactions (i.e., VIV) in the middle span and the pipe−soil interactions at the two end shoulders. With respect to the VIV, many pieces of research have been carried out with experiments, CFD simulations and empirical models (Chaplin et al., 2005). Williamson and Govardhan (2004), Wu et al. (2012) and Hong and Shah (2018) have comprehensively reviewed the background knowledge of VIV including the occurrence mechanism and the researching processes. In order to predict VIV, some empirical models such as VIVANA, Shear7 and wake oscillators have been proposed. Recently, Wang et al. (2018) performed experiments to investigate the VIV of flexible pipes subjected to oscillating flows. The vibration phenomena have also been numerically simulated by Ulveseter et al. (2018) and Lu et al. (2019).
As to the pipe−soil interactions at the shoulders, they are usually modelled with springs and dashpots . The industry guideline DNV-RP-F105 (Det Norske Veritas, 2006) recommends that the dynamic pipe-soil interactions can be simulated with linear spring and dashpot. For different types of soil (sand and clay), the spring stiffness is different. Sollund et al. (2015aSollund et al. ( , 2015b theoretically analyzed the eigenfrequencies and modal responses of freespanning pipes, in which the seabed soil was modelled with linear elastic spring. In their further study for multi-span pipes , the linear spring was also taken at the shoulders. dos Reis et al. (2018) simplified the shoulder pipe−soil interactions with generalized pinned and clamped boundary conditions. Larsen et al. (2004) combined VIVANA and RIFLEX two softwares to predict VIV of a free-spanning pipe, in which the seabed soil was described with linear or nonlinear spring. The nonlinear spring could not exert tension on the pipe. Slingsby (2015) considered the suction force from the seabed soil and simulated with a nonlinear spring that could give a small tension on the pipe. Recently, Ulveseter et al. (2016) took the static configuration of a free-spanning pipe into account and placed different nonlinear springs along the pipe to simulate the pipe−soil interactions at different places.
For those linear and nonlinear springs and dashpots, the seabed soil is simplified and some important characteristics are neglected such as plastic deformation and nonlinear hysteresis, especially for the soft clay. With respect to the touch down zone of steel catenary risers, the pipe−soil interactions are also involved, which has attracted much more attentions from the researchers . In general, two nonlinear hysteretic soil models have been proposed and widely used. They were respectively developed by Aubeny and Biscontin (2009) and Randolph and Quiggin (2009). Wang et al. (2013Wang et al. ( , 2015 linearized the nonlinear hysteretic soil model of Aubeny into spring-dashpot and then analyzed the fatigue damage of a steel catenary riser subjected to VIV. Bai et al. (2015) introduced the soil model of Aubeny into the CABLE 3D finite element code and predicted the dynamic responses of the riser at the touch down zone. In the studies of Elosta et al. (2013Elosta et al. ( , 2014 for catenary risers under random loads and waves, the vertical pipe−soil interactions were simulated with springs based on the model of Randolph and Quiggin (2009), while the horizontal pipe−soil interactions were described with bilinear and trilinear models. As for the VIV of free-spanning pipes, some researchers have also employed the nonlinear hysteretic soil models to simulate the pipe−soil interactions at the shoulders (Zhou et al., 2017;Xu et al., 2018). However, in these studies, the soil force was usually modelled with springs and dashpots. The main objective of this paper is to provide a modified approach to simulate the soil reacting force without using spring-dashpots and then analyze the influence of seabed soil on the VIV of free-spanning pipes.
The structure of this paper is as follows. In Section 2, a dynamic model for a free-spanning pipe subjected to VIV is established. In Section 3, the numerical solution methods and the validations of the present model are carried out. In Section 4, the influence of the seabed soil at the shoulders on the VIV of the pipe is analyzed in detail. Finally, the conclusions are drawn in Section 5.

Pipe vibration equation
As shown in Fig. 1, a long flexible pipe is placed on an uneven seabed. In the middle span, there is an incoming fluid flow passing the pipe with the velocity U. At the two end shoulders, the seabed soil supports the pipe. Since the pipe is slender and flexible, the rotary inertia and shear deforma-tion can be neglected. Based on the Euler−Bernoulli beam theory, the pipe's vibration equation can be expressed as: ∑ F ex in which z is the vertical displacement. The over dot and prime represent the derivations with respect to the time t and the axial coordinate x, respectively. m is the total mass of the pipe for per unit length, which includes the structural mass, internal fluid mass and external fluid added mass. c s is the structural damping coefficient, T is the axial tension, and EI is the bending stiffness. The external forces acting on the pipe are collected in , which are different in the middle span and at the shoulders.

External forces in the middle span
In the middle span, when the external fluid flow passes the pipe, some vortices will release from the pipe alternatively and periodically (Williamson and Govardhan, 2004). This will cause hydrodynamic forces acting on the pipe. With the strip theory, a cross section is taken from the pipe and the external forces are demonstrated, as shown in Fig. 2.
Under the action of the external incoming fluid flow, the pipe would vibrate in the vertical (z) and horizontal (y) directions at the same time. The vertical vibration is mainly considered in this paper because its vibration amplitudes are much larger than those in the horizontal direction. As a result, the pipe−soil interaction in the vertical direction is more severe. The external forces acting on the pipe in the vertical direction can be expressed as:  GAO Xi-feng et al. China Ocean Eng., 2020, Vol. 34, No. 3, P. 328-340 329 in which F f is the total hydrodynamic force and w e is the submerged weight of the pipe (the weight in air minus the fluid buoyancy). As shown in Fig. 2, the total hydrodynamic force F f is composed of the drag force and the lift force , which can be written as: where is the angle between the relative velocity V ( ) and the fluid flow velocity U. sin /V and =U/V. According to the Morison's equation, the drag force and the lift force can be expressed as: (4) ρ in which C D and C L are coefficients of the drag force and lift force, respectively; is the density of the external fluid flow; D is the outer diameter of the pipe. When the Reynolds number of the external fluid flow is in the subcritical regime, C D can be selected as 1.2 (Blevins, 1990). On the other hand, C L can be modelled with a classic van der Pol wake oscillator (Facchinetti et al., 2004): ε ω s ω s = 2πStU/D in which C L0 is the lift force coefficient of a stationary circular cylinder and it can be taken as 0.3; q is the dimensionless wake variable; and S are the empirical parameters which can be taken as 0.3 and 12 (Facchinetti et al., 2004), respectively; is the vortex shedding angular frequency, , and St is the strouhal number which can be taken as 0.17 (Chaplin et al., 2005).
With the above-mentioned equations, the external forces acting on the pipe in the middle span can be expressed as: in which the wake variable q is calculated with Eq. (7).

External forces at the shoulders
Similarly, a cross section is selected from the shoulders and taken into account, as shown in Fig. 3a. The external forces acting on the pipe include the soil supporting force F s , the hydrodynamic force F f and the submerged weight w e : (9) The external fluid flow at the shoulders is assumed as static. With the Morison's equation, the hydrodynamic force F f can be written as: As for the soil supporting force F s , it can be calculated with the nonlinear hysteretic soil model proposed by Aubeny and Biscontin (2009). This model describes the soil force with respect to displacement of the pipe, as shown in Fig. 3b, where z s is the displacement of the pipe penetrating into the seabed and F s is the soil reacting force. When the pipe is at different positions and vibrates with different directions (loading or unloading), the soil force is different. At first, when the pipe is penetrated into the virgin seabed, the soil force increases steadily following the path 0-1, which is the backbone curve. The soil force can be calculated with: in which a and b are the coefficients which are linked with the penetration depth, the width of the trench and the roughness of the pipe surface. For a smooth pipe and a narrow trench, the values of a and b can be taken as a=4.97, b=0.23 for z s /D<0.5 and a=4.88, b=0.21 for z s /D>0.5. S u0 is the soil shear strength at the mudline and S ug is the strength gradient. Table 1 shows the values of S u0 and S ug for soft clay in the Gulf of Mexico (Nakhaee and Zhang, 2010). It can be found that as S u0 and S ug increase, the seabed soil becomes harder and harder. When the pipe moves upwardly, the soil force decreases sharply to zero and then turns into suction caused by the adhesion between the soil and the pipe, as shown in Fig. 3b with the path 1-2. This path is the elastic rebound path. The corresponding soil force is expressed as: ω where F s,1 and z s,1 are the soil force and the penetration depth at Point 1, and k 0 is the initial slope of this path which is related to the soil shear strength through k 0 ≈660×S u0 (Nakhaee and Zhang, 2010). The parameter controls the path's asymptote, which can be taken as 0.433 from experiments. The maximum suction force at Point 2 F s,2 is: in which f suc is the suction factor. Based on experimental results, f suc can be taken as 0.203 (Aubeny and Biscontin, 2009). The penetration depth at Point 2 z s,2 is: When the pipe is further uplifted, the suction force will gradually decrease to zero, as shown in Fig. 3b in the partial separation path 2-3. The penetration depth at Point 3 z s,3 is: (15) in which f sep is the separation factor which can be taken as 0.661 (Aubeny and Biscontin, 2009). The soil force of this path can be calculated with: in which After this, the pipe is completely detached from the seabed soil in the trench. The soil force is always zero when the pipe is still uplifted as demonstrated in the full separation path 3-0 in Fig. 3b.
When the pipe is pushed down and re-contacted with the soil in the trench, the soil force follows the path 0-3-1. In the path 0-3, the soil force is zero; while in the re-contact path 3-1, the soil force increases, which is calculated with: in which If the re-penetration of the pipe is larger than that of Point 1, the soil force follows the backbone curve to a new Point 1 * , as shown in Fig. 3b. Then the soil reacting force (the dash line) is updated with respect to the new Point 1 * . These paths 0-1-2-3-0-3-1 describe the boundary loop of the soil force for full contact, full separation and full re-contact. However, when the pipe reverses arbitrarily from the boundary loop or within the boundary loop, the soil force will be different. When the pipe is reversed from the path 1-2 or 3-1 or within the boundary loop, the soil force F s can be calculated with: χ in which F s,r and z s,r are the soil force and the penetration depth at the reversal point. is a sign parameter which is 1 and −1 for loading and unloading, respectively. When the pipe is reversed from the partial separation path 2-3, the soil force is: With these equations, the external forces acting on the pipe at the shoulders can be directly calculated. Especially for the seabed soil force F s , there is no need to construct springs and dashpots.

Solution methods
Since the vibration equation of the pipe is a partial differential equation related to the time and the axial coordinate, it can be solved with the central difference formulas and the Runge−Kutta method . The pipe is equally divided into n elements and the elements are connected by nodes which are numbered as 0, 1, …, n. The second-order central difference formulas are: in which l is the length of the element (l=L/n). These formulas are substituted into the vibration equation. The partial differential equation is transformed into a series of ordinary differential equations: in which F i,s is the soil force at node i. At first, the static configuration of the pipe on the uneven seabed caused by the submerged weight should be calculated. The time related terms in Eq. (25) are neglected and these equations are reduced to: in which the soil force F i,s is calculated with the backbone curve of Eq. (11). At the two ends of the pipe, w e =F i,s (i=0 and n). For these nonlinear equations of z i (i=0, 1, …, n), the Newton−Raphson method can be employed. Once the static configuration is obtained, the dynamic response of the pipe can be solved with the fourth-order Runge−Kutta method in the time domain. Since the vibrations of the pipe in the vertical direction are mainly considered in this paper, the soil friction force in the horizontal direction is neglected. The vibration bending stress at node i can be calculated with , where E is the Young's modulus.

Model validation for VIV
It is necessary to validate the present model for simulating VIV. The experiment carried out by Lehn (2003) can be chosen for comparison. In the experiment, a model pipe was subjected to uniform incoming fluid flow by rotating in water. The experimental parameters are listed in Table 2. With these parameters and neglecting the seabed soil at shoulders, the present model to simulate VIV is carried out. The results from the experiment, CFD simulation (Wang and Xiao, 2016) and the present model are compared in Fig. 4. Fig. 4a plots the root mean square (rms) displacements along the pipe and it demonstrates that the first mode dominates this vibration. The results of the present model are larger than those from the experiment and the CFD simulation. At a point on the pipe, x/L=0.22, Figs. 4b and 4c show the vibration displacement and the response frequency, respectively. It can be seen that the results from the experiment, CFD simulation and the present model are in good agreement. Since in the present model, the in-line VIV is neglected and the coefficients of C D and C a are taken as empirical values, there are some discrepancies among these results. In general, the present model is reliable to predict the cross-flow VIV of flexible pipes.

Model validation for pipe-soil interactions
On the other hand, the pipe-soil interactions at the shoulders should also be validated. The experiment of Dunlap et al. (1990) is taken for comparison, in which a short ri-   Table 3. The soil reacting force with respect to the displacement of the pipe is plotted in Fig. 5. It can be seen that the results of the present model agree well with the results from the experiment and simulated by Aubeny and Biscontin (2009). Since the re-penetration and reverses of the pipe were not performed in the experiment, the results of the present model and those of Aubeny are not comparable. From Fig. 5, it can be known that the present model is reasonable to simulate the pipe−soil interactions.

Model validation for VIV of a free-spanning pipe
For further validation of the present model, the VIV in the middle span and the pipe−soil interactions at the shoulders are both taken into account. The studies of Larsen et al. (2004) and Ulveseter et al. (2016) can be chosen for comparison since there is lack of experiments to be compared. In their studies, the VIV was respectively simulated with VIVANA and the synchronization model (Thorsen et al., 2015), while the pipe−soil interactions were both modelled with springs and dashpots based on the RIFLEX software. The main parameters of the free-spanning pipe are presented in Table 4. With these parameters and the seabed soil listed in Table 1, the static configuration and dynamic responses of the pipe can be simulated with the present model.
For the static configuration, Ulveseter et al. (2016) used the RIFLEX to calculate the pipe deflection. The seabed soil was modelled with linear springs and the stiffness was taken as k s =40×10 3 N/m 2 . In the present model, however, three types of soft clay are calculated and the results are compared in Fig. 6. It can be seen that when the seabed soil is softer, the pipe deflection is larger. As for the lower range, the pipe's static configuration obtained from the present model is consistent with the result of Ulveseter et al. (2016). Hence, the present model is able to calculate the static configurations of free-spanning pipes.
When the VIV of the pipe is excited by the external incoming fluid flow in the middle span, the pipe−soil interactions will occur at the shoulders. The soil supporting force will decrease because of plastic deformation. Consequently, the pipe will penetrate into deeper soil, as shown in Fig. 7a. The vibrations of the pipe are around a new dynamic equilibrium which is below the static position. The bending stresses of the pipe will increase. In order to demonstrate the vibrations of the pipe and the soil reacting force, two typical points P and Q are selected from the pipe, as shown in Fig. 7a. The corresponding vibration displacements and soil forces are plotted in Figs. 7b−7c and Figs. 7d−7e, respectively. From Figs. 7b and 7c, it can be seen that the displacements of Points P and Q will become deeper and deeper and finally reach steady state with the vibration of the pipe. At the steady state, the displacement at Point P is aperiodic while that of Point Q is periodic. This is because Point P is far away from the middle span and it is less affected by the VIV and more constrained by the seabed soil. The soil forces at Points P and Q are also different. In Fig. 7d, it can be observed that the soil force at Point P is always positive, while Fig. 7e shows that the soil force at Point Q moves between positive and negative, which almost follows complete boundary loops.    In terms of the vibration amplitudes and bending stresses along the pipe, the results of the present model are compared with the results of Larsen et al. (2004) and Ulveseter et al. (2016), as depicted in Fig. 8. It can be seen that the vibration amplitudes are large in the middle span. While the bending stresses are large in the middle span and at the connections between the span and the shoulders. These places are dangerous because they are more easily to get fatigue damage. Fig. 8 demonstrates that the results predicted by the present model agree with the results from Larsen et al. (2004) and Ulveseter et al. (2016). Because the seabed soil at the shoulders is modelled differently, there are some reasonable discrepancies among these results. It can be known that the present model can simulate the vibrations of free-spanning pipes subjected to VIV in qualitatively.

Parametric study
The seabed soil at the shoulders may be different which can affect the vibrations of the free-spanning pipe undergoing VIV. Based on the present model, the influence can be analyzed. In this section, the influences from the suction factor f suc , separation factor f sep , mudline shear strength S u0 and strength gradient S ug are going to be analyzed. In these analyses, only the concerned parameter is changed while the other parameters are kept constant.

Influence of the suction factor f suc
In the subsea, the seabed soil may be different and it may have different suction forces. The suction force is larger when the soil cohesion is stronger. In Eq. (13), it can be seen that the maximum suction force is related to the max-  imum supporting force through a factor f suc . In order to analyze the influence of the maximum suction force on the VIV of a free-spanning pipe, the suction factor can be selected as f suc =0.1, 0.203, 0.3 and 0.4. Then, the vibrations of the pipe are simulated and the results are plotted in Figs. 9 and 10. With the increase of the suction factor f suc , the maximum suction force is getting larger and larger. Fig. 9a shows that the penetration of the pipe at Point P into the seabed is becoming smaller. Point Q on the pipe has the same characteristics, as shown in Fig. 9b. Figs. 9c and 9d show that the soil reacting force decreases with the increase of the suction factor f suc . In Figs. 10a and 10b, it can be seen that as the suction factor f suc increases, the vibration amplitudes and bending stresses both decrease. The reason for this may be that for a larger suction factor, the maximum suction force is larger and the pipe will be more strongly adhered by the soil.

Influence of the separation factor f sep
The suction factor f suc determines the maximum suction force of the seabed soil, while the separation factor f sep determines the width of the range of the suction force, as given in Eq. (17). With the increase of f sep , the suction range becomes wider. To analyze this influence on the vibrations of the pipe, the separation factor f sep is taken as 0.5, 0.6, Fig. 10. VIV of the pipe for different suction factor f suc : (a) vibration amplitudes and (b) bending stresses. 0.661, 0.7 and 0.8, and the calculation results are plotted in Figs. 11 and 12. Figs. 11a and 11b show that the penetration depths of the pipe at Points P and Q decrease with the increase of the separation factor f sep . In Figs. 11c and 11d, it can be seen that the soil forces also decrease as the separation factor f sep increases. It is more clearly in Fig. 11d because the suction force appears at Point Q. With respect to the vibration amplitudes and bending stresses, Figs. 12a and 12b show that they decrease with the increase of the separation factor f sep . This is because for a larger separation factor, the suction force is wider and the seabed soil is more cohesive.

Influence of the mudline shear strength S u0
The seabed soil may have different shear strength and it can affect the VIV of the free-spanning pipes. When the shear strength increases linearly with respect to the depth, at different depths it can be expressed as S u =S u0 +S ug z s . In this subsection, the influence of the mudline shear strength S u0 on the vibrations of the pipe is analyzed. The values of S u0 are taken as 0.6, 0.8, 1.2, 1.6 and 2.0 kPa, and the calculation results are plotted in Figs. 13 and 14. In Figs. 13a and 13b, it can be seen that as the mudline shear strength S u0 increases, the penetration of the pipe into the seabed soil largely decreases. Figs. 13c and 13d show that the vibration  region of the pipe in the soil also decreases with the increase of S u0 . This is because the seabed soil becomes harder with the increase of the mudline shear strength. In Fig. 14a, it can be observed that the vibration amplitudes of the pipe are slightly changed with the increase of the mudline shear strength S u0 . However, Fig. 14b shows that the bending stresses of the pipe are larger for a larger S u0 . The reason for this may be that as the mudline shear strength increases, the vibrations of the pipe will be more constrained.  Figs. 15c and 15d, it can be seen that the pipe's vibration region also decreases with the increase of S ug . This is because the seabed soil becomes harder as the strength gradient increases. In Fig. 16a, it can be seen that the vibration amplitudes of the pipe are slightly fluctuated with the increase of S ug . For the bending stresses along the pipe, Fig. 16b shows that they are increased with the increase of the strength gradient S ug . This is  because for a larger strength gradient, the pipe will be more constrained.

Conclusions
The free-spanning pipes are usually subjected to currents in the subsea and the vortex-induced vibrations would be excited. The pipe−soil interactions at the shoulders could significantly affect the VIV of free-spanning pipes. In this paper, the seabed soil reacting force is directly calculated with a nonlinear hysteretic soil model. The lift force caused by vortices shedding in the middle span is modelled with a classic van der Pol wake oscillator. Based on the Euler− Bernoulli beam theory, a dynamic model for a free-span-ning pipe undergoing VIV is developed. At first, the static configuration of the pipe is calculated and then the VIV is simulated. With the present model, the influence of the seabed soil on the vibrations of the pipe is analyzed. The conclusions of this paper are summarized as follows.
(1) When the maximum suction force of the seabed soil at the shoulders is larger, the penetration depth of the pipe into the seabed will be smaller and the pipe will vibrate with smaller amplitudes and bending stresses. The pipe is safer.
(2) When the suction range of the seabed soil is wider, the pipe undergoing VIV will penetrate into shallower seabed and the vibration amplitudes and bending stresses will decrease.  (3) With the increase of the shear strength at the mudline, the penetration depth of the pipe will become smaller, the vibration amplitudes will slightly change while the bending stresses will increase. As a result, the pipe is subjected to larger fatigue damage.
(4) When the seabed soil has a larger strength gradient, the penetration depth and vibration region of the pipe will become smaller, the dynamic bending stresses will increase and the pipe is more easily to become fatigue failure.
In practice, under the action of an incoming fluid flow, the free-spanning pipe will deflect in the horizontal direction and vibrate in this direction, i.e., the in-line VIV. Moreover, the seabed soil may be different along the pipe and the pipe may be multi-span. These effects should be considered in the future in order to get a better understanding of free-spanning pipes undergoing VIV.