Development of submarine cable hydrodynamic loads calculation program using MATLAB

A MATLAB program was developed to estimate the hydrodynamic forces due to wave and current on an un-trenched submarine cable of diameter 9 in resting at the bottom of a 14-m depth sea, utilizing the Fourier decomposition method. The aim of this part of the research is to provide an accurate estimation of the hydrodynamic forces acting on the cable and to present a comparison between the Morison force model and Fourier decomposition force model. In addition, investigating the methods of fixation for submarine pipes and cables as the commonly used methods is of high cost for both manufacturing and installation and maintenance processes has many drawbacks. This research is concerned with the ability of the cable to be stable on the seabed using its own weight and, if required, estimating extra chain weight to be added helping stabilize the cable on the seabed. A validation of the developed program with UWAHYDRO program developed in the Western Australia University is introduced in this article.


Introduction
There is a need to utilize 70% of undiscovered land under seas and oceans to connect continents through power and communication cables, establish renewable energy farms in the open seas, or even excavate for offshore sources of hydrocarbon fossil fuels. Lately, pumps and compressors are placed at the bottom of the sea on the submarine pipelines for transferring fossil fuels and natural gas from offshore wells to onshore facilities. All the previous mentioned applications require the presence of cables at the bottom of seas or oceans with different depth, different meta-ocean properties, and different soil types. Therefore, the research for the stability of these cables is of huge importance due to its critical application and need. A MATLAB program was developed to estimate the hydrodynamic forces due to combined effect of waves and current per unit length acting on a 9-in cable in diameter, resting on a calcareous sand seabed with 14 m depth. The developed program was based on the theories and recommended practice presented in Ghebreghiorghis 1 based on DNV-RP-F109 ''On Bottom Stability Design of Submarine Pipelines'' published by DN Veritas 2 in 2004 and modified at 2014, 3 as there was no other known code or standard for such matter. The region meta-ocean data obtained is as follows: significant wave height H s of 4.54 m, associated peak period T p of 9.65 s, and a steady current velocity U c of 0.66 m/s measured 1 m above seabed. The above-mentioned wave properties were obtained from the meta-ocean data for a 100 years return period for both wave and current. Current velocity is assumed to be acting perpendicularly on the cable as shown in Figure 1.
The spreading function was not utilized to estimate the worst hydrodynamic forces that may occur as the wave is estimated to have a 90°angle of attack on the cable. First, the Joint North Sea Wave Project (JONSWAP) spectrum was used to simulate the wave elevation spectrum based on the wave properties in the frequency domain according to Branlard 4 where g is the peak enhancement parameter, a is the Generalized Phillip's Constant, and s is the spectral width parameter. They take the following forms where H s is the significant wave height, T p is the peak period. The peak enhancement parameter value varies approximately from 1 to 5 and with a mean value of 3.3 and variance of 0.62. The JONSWAP spectrum reduces to Pierson-Moskowitz spectrum when the value of g is 1.
Then, it is transferred at the cable level on the seabed at 14 m depth, as a velocity power spectrum using a transfer function suggested by DNV-RP-F101 2 and Azouz et al. 5 and Li 6 D being the water depth, and k is the wave number obtained from a Dispersion Relation assuming that tanh(x) ' x for small values of x following 7 Note that for more accurate prediction of the wave number, iterative method should be considered in equation (8) other than the above-mentioned assumption.
After obtaining the velocity power spectrum, Fourier series is used to transform the frequency domain spectrum to a time domain profile history of the wave particle velocity at the cable level where a i is the amplitude, u i is the phase angle which is random in nature and normally distributed from [2p to p], t is the storm/simulation duration.
Determining the values of the amplitude and phase for each frequency gives the irregular shape for the velocity profile.
Then the acceleration profile is calculated from the velocity profile through derivation of the velocity profile with time The Fourier decomposition method explained by Chakrabarti 8 is adopted to simulate the hydrodynamic drag, inertia, and lift forces acting on the cable. The general expression for the periodic quantity is where a i and b i are the Fourier coefficients driven from experimental measured forces. A more developed form, developed from the general expression, is as follows where C i and B are defined as follows To obtain non-dimensional force coefficient from the Fourier analysis, measured physical experimental forces were neutralized by ½ÁrÁDÁu(t) 2 , where u(t) is the wave velocity profile at the cable level. The horizontal drag force was obtained by subtracting the inertia component from the measured horizontal force where C a is the added mass coefficient with a value of 2.29 as determined from the potential flow theory 9 and a W (t) being the wave acceleration profile at the cable level.
The decomposition of the non-dimensional horizontal drag and vertical forces has the following expressions of Abdolmaleki 10 Also, for the purpose of comparison, the Morison force model was used to calculate the forces on the same cable using the same wave properties where U W is the wave velocity, U c is the current velocity, d is the pipe/cable diameter and C L , C D , and C M are the lift, drag, and inertia coefficients, respectively. Extensive physical tests were conducted to define the hydrodynamic forces coefficients of the Morison equation. However, as highlighted by Sorenson et al. 9 and Allen et al., 11 the Morison equation does not accurately describe the force time-series for combined flow.
The inertia force is calculated separately according to the equation proposed by Sorenson et al. 9 having the same formula for both models where a(t) is the wave acceleration. A discussion regarding both methods in terms of hydrodynamic loads estimation for irregular wave and associated current is presented in this article.
To adapt the Fourier method to an irregular wave profile, the nth-moment spectrum equation is used to obtain the significant wave velocity U s of the velocity power spectrum and its associated peak period T p where KC is the Keulegan-Carpenter number and M is the ratio of the current to wave velocities. Then using these two properties, the Fourier coefficients can be selected from the experiments held by Sorenson et al. 9 at the Danish Hydraulic Institute and was presented in the AGA PR-178-04405 12 report based on the Keulegan-Carpenter (KC) and the current to wave ratio M, with another two inputs being the cable and soil roughness. At first, a set of nine harmonics [9 coefficients + 9 phase angles] were assumed to provide an acceptable accuracy in the representation of the forces, then it was found that only the first five harmonics can be used, as they give a sufficient description of the forces because there was no significant change after the fifth harmonic. In the UWAHYDRO data input, a pipe of 250 m in length is to be simulated, while in the MATLAB program, a pipe of 1 m length was studied. The reason is that since the spreading function was not activated, therefore, the full forces load is assumed to be perpendicularly hitting the pipe. Therefore, the effect on one-node with length of 1 m is the same as for the full 250 m length except for the fixation points (start and end of pipe).

Validation of the proposed MATLAB program
Results of velocity, acceleration profile estimations of MATLAB Program, versus UWAHYDRO are presented in Figure 2. The Matlab M-file and function file are presented as supplemental material 1 and 2 respectively.
The velocity and acceleration profiles show near perfect match except for a small shift in phase; this is due to the random phase angle used to simulate an irregular sea state, which makes it hard to achieve a 100% match. As for the hydrodynamic forces' profiles, the drag, lift, and inertia forces presented in Figure 3 show a great agreement in terms of profile and values with the estimated forces from UWAHYDRO. The smaller time increment 0.1 s chosen showed more regressions,  even higher peaks than the 1-s time increment used in UWAHYDRO program, which helps in increasing the accuracy of the estimated forces and anticipating even higher loads than smaller time increments utilized in most programs.

Case study
A case study of an un-trenched submarine cable of diameter 9 in resting at the bottom of a 14-m depth sea has been performed to estimate the hydrodynamic forces due to wave and current. Data input to the proposed MATLAB program are presented in Table 2.

Results and discussion
The time increment was decreased so to detect peaks that can be missed using normal time steps as 1 s, therefore increasing the accuracy of the calculated hydrodynamic forces. To overcome the random parameter    used, seven runs were conducted and the average mean of all seven runs was calculated. Figure 4 presents JONSWAP wave elevation spectrum at sea surface and velocity power spectrum at cable level, while the transfer function results are shown in Figure 5. The calculated velocity and acceleration profiles at cable level is shown in Figure 6. Concerning the calculated drag, inertia, and total horizontal forces acting on the cable is presented in Figure 7, while the calculated Lift force acting on the cable is shown in Figure 8.
Results showed the superiority of the Fourier model over Morison's model in estimating the hydrodynamic forces. In the case of the drag force, the Morison model overestimate the forces, which means it anticipates that the cable will move a much longer distance than that in reality, resulting in over-weighting the cable causing a higher cost in manufacturing and installation. The inertia force for Fourier and Morison models is the same as it takes the same formula proposed by Sorenson et al. 9 and AGA. 11 The lift force estimated by Morison model shows a good agreement with the Fourier model in terms of profile shape, without taking into consideration the down force acting on the cable due to its small diameter. Which indicates that the hydrodynamic uplift force somewhere along the force time history will help fixing the cable into the seabed and increase its embedment. The Fourier model results in the lift force profile shows peaks that Morison model comes short to anticipate. Therefore, in terms of drag and lift forces, the Fourier model is a much better option in estimating a more realistic and accurate hydrodynamic forces than Morison.

Conclusion
A MATLAB program was developed to estimate the hydrodynamic forces due to wave and current on an un-trenched submarine cable of diameter 9 in resting at the bottom of a 14-m depth sea, utilizing both Fourier decomposition method and Morison force model. JONSWAP spectrum was used to simulate the wave elevation spectrum based on the wave properties in the frequency domain. The developed program was validated using UWAHYDRO program developed in the Western Australia University and results showed very good agreement. Furthermore, a case study was presented to provide a comparison between force calculations obtained from Fourier decomposition method and Morison force model. Thus, in terms of drag and lift forces, the Fourier model is a much better option in estimating a more realistic and accurate hydrodynamic forces than Morison. The next phase of the research will utilize the resultant hydrodynamic forces from the Fourier model to be an input to the cable-soil interaction model, also estimating the correction factors for the drag, inertia, and lift forces due to cable movement, penetration, and lifting of the soil.

Declaration of conflicting interests
The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Funding
The author(s) received no financial support for the research, authorship, and/or publication of this article.

Supplemental material
Supplemental material for this article is available online.