Dynamics of risers analysed by means of the segment method, with consideration of bending and torsional stiffness

The paper presents the application of the finite segment method to the analysis of coupled bending torsional vibrations of risers. The method is formulated by means of joint coordinates using multibody methods for kinematics and dynamics. A new approach to calculating bending and torsion moments is presented. The mathematical model and computer program enable us to analyse both free and forced vibrations of risers caused by the motion of the base (vessel or platform) as well as hydrodynamic forces. The model is validated by comparing frequencies of free and forced vibrations calculated from the authors’ own models with the results presented by other researchers. Natural frequencies are also compared with analytical solutions. The influence of sea currents and of the initial twisting of the riser on its natural and forced vibrations is analysed.


Introduction
Problems concerned with static, kinematic and dynamic analysis of devices used for extracting and transporting hydrocarbons from the seabed are extremely complex. In the case of risers the difficulties are due to the fact that it is necessary to consider large deflections and the interaction of the sea, which lead to nonlinear equations of dynamics. In order to consider the influence of the sea, base motion caused by waves, currents, buoyancy forces, hydrodynamic drag as well as change in the added mass have to be taken into account.
The most popular methods used for modelling risers are the finite element method [1,2], finite difference method [3,4], lumped mass method [5,6], finite segment method [7][8][9] and the rigid finite element method [10][11][12][13][14][15][16][17]. The professional software packages both for general purpose (Abaqus, Ansys) and specifically for marine applications (Orcaflex) generally use the finite element method. A review of models and methods used for statics and dynamics of risers, including different excitation problems, is presented in [18]. More recent discussions of models and methods used for dynamics as well as control techniques can be found in [19].
Longitudinal flexibility is usually omitted in dynamic analysis of risers due to the fact that elongation of typical risers of 2 to 3 km length is usually less than 1 m. Torsional flexibility is also rarely considered [6,[20][21][22][23][24][25][26][27], but even then is usually examined only when the multi-layered structure of pipes is taken into account [22][23][24][25][26][27]. Geometric and bending stiffness are dominant although torsional vibrations may cause damage due to material fatigue. Some research is devoted to stability of risers under torsion and bending-torsion coupling [23,24] as well as structural analysis of the layer pipes and axial-torsional coupling [22,25,26]. Analysis of different modelling approaches to flexible pipes with layers under tension, torsion and bending loads is presented in [27].
The finite segment method (FSM) and the rigid finite element method (RFEM) are popular not only due to their simple physical interpretation. Customised formulations of these methods can be numerically effective and thus easy to apply for solution of optimisation problems [15][16][17][28][29][30].  The main idea of both methods is to divide a flexible link into rigid finite elements (rfe) reflecting mass features, and spring-damping elements (sde) reflecting flexible and damping features of the link. The motion of the subsequent elements is usually described in relation to the previous element [7,31,32]. The Euler angles are used in order to describe the motion, and continuity of displacements in connections of the elements is assumed. Such an approach means that the generalised coordinates defining the motion of element i depend not only on its own coordinates but also on generalised coordinates of all preceding elements and coordinates of a chosen point of the first segment. As a consequence the mass matrix of a segment not only depends on all these coordinates but is also full [32], whereas the stiffness matrix of the system has a simple form for linear physical and geometrical relations.
This paper further develops models described in [30,31]. Paper [30] describes the segment method using joint coordinates and its application to analysis of spatial bending vibrations, while torsional vibrations are omitted. The models presented there are numerically effective and they are used for the solution of the dynamic optimisation problem for stabilisation of forces and moments in the connection of the riser with a wellhead.
The model used in this paper is described in detail in [31] and it takes into consideration both bending and torsional flexibilities. The essential difference between the model presented here and the model presented in [31] is a different method of calculating bending and torsional moments acting at spring-damping elements. In paper [31], generalised forces reflecting bending and torsion are calculated by means of the energy of spring deformations, while in this paper special formulae for calculating moments are presented. Thus, the paper presents a nonlinear model of a riser including bending-torsional coupling, influence of the marine environment such as buoyancy, drag forces and sea current, and also motion of the platform or vessel that the riser is attached to. The model is validated by comparing the authors' results with those presented by other researchers. Finally, the influence of riser twisting on natural frequencies and on riser vibrations during the motion caused by the surface vessel is examined.

Riser model
The main idea of the finite segment method is to replace a continuous link with a system of rigid segments (rfe) connected by means of spring-damping elements (sde). The division is carried out in two steps: first the discretised link is divided into elements of equal length reflecting mass features of the link and then spring-damping elements are placed in the middle of primary elements (points A 1 … A n ) reflecting flexible features of the segment.
According to [32,33] coefficients characterising an sde for a homogenous beam can be calculated from the following formulae: where E , G are Young and Kirchoff moduli, respectively, J is the inertial moment of cross-sectional area with respect to the bending axis, J 0 is the inertial moment of cross-sectional area with respect to the torsional axis (lengthwise).
Thus a flexible link is replaced by a system of n + 1 rigid elements (segments S 0 ÷ S n ) connected at points ( A 1 ÷ A n ) by spring-damping elements ( e 1 ÷ e n ).
As a result, the motion of the divided link, presented in Fig. 1, is described by m = 3 + 3(n + 1) generalised coordinates, which are the components of the following vector: The model used in this paper is derived in [31] and thus its description here is limited to parts which are relevant for the present discussion. Coordinate system {i} � is assigned to segment S i . As shown in Fig. 2b, angles i , i , i , called heading, attitude and bank angles [34], are ZYX Euler angles, the name of which describes the order of rotations assumed for the model. Coordinates of point A of segment i in global coordinate system {} can be calculated according to the formula: where l j is the distance between points A j and A j+1 , The equations of motion are derived from the Lagrange equations of the second order. Having omitted detailed derivations, for which we refer the reader to [31], the first variational derivative can be written in the following block-matrix form: where rr , r , , r , are mass matrices and vectors referring to respective generalised coordinates Bending and torsional moments induced by spring damping elements are calculated as follows. Using stiffness coefficients defined in Eq. (1) we introduce subscript i to denote an sde, and thus c (b) i is the bending stiffness coefficient of sde i connecting segments i − 1 and i . The bending moments acting on segments i − 1 and i are then expressed in the form:  where cos Δ i = cos Δ i cos Δ i . Moments caused by torsion can be calculated according to the formula: Angles Δ i , Δ i , Δ i denote relative rotations of segment i with respect to segment i − 1 and for general discussion we assume that they are not small (Fig. 2). They can be calculated using formulae: where i,Δ i,j denotes element of the i-th row and j-th column of matrix i,Δ : Due to the deformation of sde i , the generalised forces in the form of the following vector will occur in the equations of motion: Gravity, buoyancy forces and hydrodynamic drag forces can be incorporated into the model by means of n the generalised forces as described in [31], and finally the equations of motion take the form: Reactions have to be introduced and equations of motion have to be supplemented with constraint equations.

Reactions and constraint equations
Numerical simulations are carried out for a riser (Fig. 3) whose motion at the top point (A 0 ) is known and whose coordinates are known functions of time: In this case the coordinates of point A 0 can be eliminated as generalised coordinates, and the equations of motion can be written as follows: where 0 is the vector of reactions at point A 0 .
Equation (12.2) takes the form: from which acceleration vector ̈ can be calculated. Then from (12.1) the vector of reactions at A 0 can be calculated in the form:

Fig. 3 Rigid connection of the riser with a wellhead
In the case when there are constraints imposed on point E = A n+1 (Fig. 2), for example when this end of the riser is connected with a wellhead, it is assumed that this connection is rigid (Fig. 3).
In such a case the influence of force E on the riser has to be considered and the respective constraint equations have to be formulated. Due to reaction E at point E , coordinates E are defined as follows: and the vector of reactions occurs in the equations of motion in the following form: The constraint equations can take the acceleration form: where E = E ̇ ,̇E ,0 ,̈E ,0 . Initial conditions are defined after assuming ̇=̈= 0 and solving equations of statics resulting from (15,16).The nonlinear set of algebraic equations is solved by the iterative Newton method.
Equations (13,17) are integrated using the Runge-Kutta method of the fourth order with a constant integration step. If it is necessary to calculate frequencies of free vibrations, we apply a procedure described in [14]. The Baumgarte method [35] is applied for stabilisation of the constraint equations.

Validation of the model
In order to prove the correctness of the models and computer programmes elaborated on the basis of these models, a validation is carried out. This involves comparing frequencies of free vibrations calculated from the authors' own models with the results presented in [36] obtained by the differential transformation method. It is important to note that the results of Chen et al. [36] are compatible with those presented by other researchers obtained by means of five different methods. As for large deflections and forced vibrations of risers, we compare our results with those presented in [37].
T Ë= E Natural frequency analysis is concerned with vibrations of the riser shown in Fig. 4, with parameters given in Table 1.
It is assumed that the riser is submerged in water, and vibration analysis takes into account the gravity force, drag force, buoyancy force and fluid inside the riser. Extensibility is omitted and bending in plane xy as well as torsional vibrations are analysed.
In order to calculate the free vibrations, we use an approach inspired by [38] and described in detail in [14] for a system with geometrical constraints. Table 2 contains the first five frequencies of the riser calculated for various numbers of segments n = 25, 50, 100 and 200 and results presented in [36]. Four different supports shown in Fig. 4b are analysed.   Table 2 shows that the segment method proposed in this paper gives results compatible with those obtained by means of different methods, especially those from [36]. The differences are not larger than 1.5% even for the number of segments as small as n = 25.

Analysis of the results presented in
We also compare results obtained by means of the segment method (SM) with analytical values calculated according to the following formula [39]: where L is the length of the riser, i is the mode, T is the top tension, m * is the mass per unit length. This formula is valid for H-H boundary conditions and buoyancy forces acting along the whole length of the riser. The comparison of the first five frequencies calculated for n = 100 is presented in Table 3. Very good correspondence of results is achieved. Table 4 presents torsional frequencies calculated by means of the segment method and using an analytical formula. The results refer to the riser with both ends fixed; rotation about axis x ′ at points O and E is eliminated. The EI m * analytical formula for bending natural frequencies is as follows [40,41]: where G is Kirchoff's modulus, A out =  Table 4 are obtained for the riser described in Table 1.
Analysis of the results implies the compatibility of the segment method with the analytical solution; the error  smaller than 1% for the first natural frequencies is achieved even for n ≥ 25 . One can also see stabilisation of the results as the number of segments n increases. It is important to emphasise that the values of torsional natural frequencies for the riser described in Table 1 are two orders of magnitude larger than frequencies of bending vibrations. Paper [37] presents a riser shown in Table 5

Numerical simulations
This section presents numerical simulations based on the method described above and concerned with both statics and dynamics of the riser with parameters presented in Table 6.
Both ends of the riser are fixed by means of springs with stiffness coefficients c (t) 0 , c (t) E equal to 10 8 Nm∕rad . Analysis carried out in this paper also includes the influence of different sea currents on the riser. Three cases are analysed: C0 -without current; C1 -current changing (Fig. 7a) Coordinate y E depends on the deflection of the riser due to the sea current acting in the direction of z axis. Figures 8 and 9 show how the static deflection of the riser is influenced by sea currents. It can be seen that the differences are considerable when different profiles of currents are taken into account. Figure 10 exemplifies how the profile of the sea current influences the torsion angle and bending moment along the riser in the static problem.
The influence of torsional loading on the riser is also examined by introducing initial twisting at the bottom end of the riser. The change of coordinate z with respect to different twisting angles when there is no current is shown in Fig. 11. The model presented enables us to carry out dynamic analysis for displacements occurring at discretisation points. The vibration analysis of a riser connected to a moving vessel is presented below. It is assumed that the hang-off point of the vessel moves in the horizontal xz plane along the trajectory shown in Fig. 12a, and its velocity changes in directions x and z in the way presented in Fig. 12b.
The functions assumed are similar to the harmonic functions: where a x = 5[m], a z = 4[m], = 2 ∕T and per iod T = 10[s].
These functions are changed for t = 0 and t = 2T in order to both fulfil homogenous initial conditions for the equations of motion of the riser and to stop the motion of the riser for t = 2T : This causes deformation of the trajectory of point A 0 , as seen in Fig. 12a.  The period of functions x(t) and z(t) is close to the most common frequencies of surface waves. The period of vibrations of the riser considered is about 0.4 s.
The calculations are carried out for a riser with the data from Table 6 over period t ∈ ⟨0, 30s⟩ with integration step h = 0.0005s and number of segments n = 25 . Figure 13 presents how the profile of the sea current influences torsion angle and torsion moment M T in the middle of the riser ( s = 250m ), while Fig. 14 shows this influence on displacements of the middle point of the riser.
The influence of additional loading on the riser either in the form of sea currents or twisting is important when stress analysis is concerned. Figure 15 shows the differences in values of the torsional moment in the middle of the riser under different current profiles and twisting angles acting at the bottom end of the riser.
It can be seen that initial twisting results in smaller increase in the moment than the sea current does. Analysis of the results from Figs. 10, 11, 13, 14 and 15 indicates that the intensity of sea currents and initial twist of the riser have an essential influence both on the displacements of the riser and on the loading. Thus, consideration of torsional stiffness in the mathematical model is important since this can change both bending and torsional stress in the riser.

Final remarks
The paper presents a spatial model of a riser for static and dynamic analysis with consideration of bending and torsional stiffness. The finite segment method is used for discretisation of a flexible link. The computer programme developed on the basis of the model presented can be used for nonlinear static analysis of displacements and loads as well as for analysis of both free and forced vibrations caused by the motion of the base (vessel or platform) to which the riser is attached. The authors' own results concerned with free and forced vibrations in water, when buoyancy, drag forces and added mass are taken into account, are compared with those presented in the literature, and good agreement is obtained. Static and dynamic analysis is concerned with the influence of sea currents and initial twisting of the riser on bending and torsional displacements. The results obtained show a considerable influence of the presence and intensity of the sea currents as well as initial torsional loading. Dynamic analysis deals with vibrations of a riser caused by the horizontal motion of the base along almost an elliptical path with diameters of 10 m and 8 m. The influence of the sea current and twisting is illustrated by vibrations of a point in the middle of the riser length. Time series of the torsional angle and torsion moment also show the considerable influence of the sea current and twisting. A strong sea current acting in a direction perpendicular to the riser causes a torsion of between ten and twenty degrees. Thus the stress due to torsion may be considerable.
Calculation time for simulation of dynamics of the riser over a time period of 30 s with an integration step of 0.0005 s are less than 3 min on an average PC computer. This means that the programs developed can be very useful in designing such devices.  Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.