Effects of Different Design Parameters on the Vortex Induced Vibration of FRP Composite Risers Using Grey Relational Analysis

Risers are indispensable components of offshore platform systems that connect the wellhead at the sea bottom and to the platform at the sea surface and are normally made of high grade steel. Nowadays, fiber reinforced polymer (FRP) composite has been recognized as an attractive alternative riser material. Similarly to steel risers, FRP composite risers are also vulnerable to vortex induced vibration (VIV), and the effects of the composite makeup of these risers on VIV are the subject of the present investigation. Three risers (the tailored design composite riser, the composite riser with orthogonal reinforcements and the steel riser), three current velocities (0.36 m/s, 1.22 m/s and 2.13 m/s) and three water depths (12.5 m, 25 m and 37.5 m) are considered. In total, 9 study cases using orthogonal array (OA) sampling are investigated to study the risers’ VIV characteristics. The computational fluid dynamics (CFD) simulations with coupled fluid–structure interaction (FSI) are used to obtain the risers’ natural frequencies, global displacements, global stresses and the stress distributions in each composite lamina. The effect of 5 parameters (Etension, Ebending, L/Douter, tension force and current velocity) on the VIV amplitude in a cross flow direction of the risers is analyzed using Grey relational analysis (GRA) and the Grey relational grade of these parameters are: r 05 ( Velocity ) > r 03 ( L / D outer ) > r 04 ( Tension   Force ) > r 01 ( E t e n s i o n ) > r 02 ( E b e n d i n g ) .


Introduction
Risers are a critical component of offshore platform systems, and over 80% of offshore oil and gas accidents are related to the failure of these risers. The risers connect the wellhead (sea bottom) and the platform (sea surface) and can be divided into two main types: A drilling riser (which guides the drilling equipment and conducts the drilling fluid) and a production riser (which transports and exploits the gas and oil).
Presently, these risers are typically made from high grade steel. However, because of their extreme weight and problems with corrosion, fiber-reinforced polymer (FRP) composite risers have become an attractive alternative. Compared to steel, FRP composites have lower density, better mechanical properties, better corrosion and fatigue performance, and better thermal insulation [1,2]. In addition, the FRP composites can be tailored-designed by optimizing the fiber and matrix combinations, fiber orientations, the thicknesses of the liner and composite lamina, and the stacking sequence, according to are included for each design variable, and the strength level number, t, and k the number of design variables [26]. More specifically, in this study, the OA program (bose) is employed as OA (q 2 , k, q, 2) [26]. Computational fluid dynamics (CFD) simulations with coupled fluid-structure interaction (FSI) are used to study the VIV characteristics of FRP composite risers', which is one of three basic different methods. The three basic methods are: Semi-empirical models, Navier-Stokes models (or CFD methods) and simplified wake models [27]. For semi-empirical models, the three main types were included: (1) The wake-oscillator models, in which common terms in equations were used to couple the body and the wake oscillations; (2) the single degree-of-freedom (SDOF) models, in which a single dynamic equation with aeroelastic forcing terms on the right-hand side was employed; and (3) the force-decomposition models which depended on the certain components of the forces on the structure from experimental results [28]. For CFD methods, the flow field around the cylinder is computed by solving numerically the unsteady Navier-Stokes equations [27]. For simplified wake models, the flow forces were modelled using the simplified models and a strip theory was used to reconstruct the 3-dimensionality [27].The significant sequence of different parameters on the risers' VIV characteristics are explored using the Grey relational analysis (GRA) which is a small sample, a small amount of calculation and no requirement for independency [29]. More specifically, GRA is utilized to determine the relational grade between the reference sequence (the VIV amplitude in cross flow direction) and a given set of parameters (E tension , E bending , L/D outer , tension force and current velocity). Then, the best one can be found by comparing the resultant relational grades, i.e., GRA can be regarded as a measure of similarity for finite sequences. The main procedure of GRA includes: (1) the Grey relational generation; (2) the Grey relational coefficient calculation; and (3) the Grey relational grade calculation. When a variable sequence has the highest Grey relational grade, it is the most important one [30].

Materials
The conventional steel riser (riser 3) is made from X80 steel as a benchmark. For the composite riser, AS4-epoxy was selected for the composite layers while titanium alloy was used for the inner liner. More specifically, AS4 is a type of high strength carbon fiber with elastic modulus of 235 GPa and an ultimate strength of 3590 MPa in the fiber direction. The thermoset epoxy which has good mechanical properties and well-suited for filament winding was selected as the matrix. The fiber volume fraction of AS4-epoxy composite is 0.6. The material properties of these materials are presented in Tables 1  and 2. Notably, for steel and titanium, the bilinear kinematic hardening model was considered while for AS4-epoxy, an elastic model was used.  Table 3 shows the geometries of the risers considered in this study. The internal diameters (IDs) of all three risers were fixed at 250 mm, while the thickness of each riser was determined by the previous studies [4,17,18,31] to accommodate all the local and global loads considered. More specifically, the thicknesses for risers 1, 2 and 3 were 30.5 mm, 39.5 mm and 25 mm, respectively. Both FRP composite risers, 1 and 2 (45.7 kg/m and 59.6 kg/m, respectively) achieved significant weight saving compared to the steel riser (169 kg/m) and the FRP composite riser 1 even led to a 23.3% weight reduction to the conventional FRP composite riser 2 [4]. These geometries of composite risers 1 and 2 have been verified by KIM [15] and Amaechi [3].

Riser Geometries
In terms of the length for the risers, 12.5 m, 25 m and 37.5 m are considered.

Flow Situation and Forces
In order to keep the risers' vertical position, the top tension forces, which are 1.5 and 2 times of the structural weight for steel and composite risers, are applied. The gravity and buoyancy loads are determined by the geometry and material density of each riser.
The 1 year winter storm (0.36 m/s), 100-year hurricane (1.22 m/s) and 100-year loop current (2.13 m/s) situations in Gulf of Mexico [15,32] were employed in this study based on the requirement of API Recommended Practice [33]. For the sea water, the density was 1024 kg/m 3 and the kinematic viscosity coefficient υ was 1.06 × 10 −6 . The Reynolds number (Re) and turbulent intensity (I) are calculated by Equations (1) and (2).

OA Sampling
The study cases in this paper were determined by OA sampling (bose) [26] based on the riser types, current velocities and water depths (Table 4). First, a pre-stress modal analysis was used to analyze the natural frequencies of all the risers. More specifically, the modules of Static Structural and Modal in ANSYS are required for the steel riser (riser 3) while the modules of ACP (Pre), Static Structural and Modal in ANSYS are involved for the FRP composite risers (risers 1 and 2).
The modules of geometry, ACP (Pre), transient structural, fluid flow, ACP (Post) and system coupling [34,35] were then employed to achieve the VIV simulation considering the fluid structure interface (FSI) for the composite risers 1 and 2. For the steel riser 3, ACP (Pre) and ACP (Post) were not needed.
The layered-structure method (LSM) [34], which could present the results in each composite layer, was employed for the FRP composite riser simulation. More particularly, in the module of ACP (Pre), the material properties, the thicknesses and layer numbers of composite laminate and liner, the fiber orientations and the stacking sequences of laminate were imported. In the ACP (Post), the stress and deformation response of every composite laminae can be obtained.
Next, the geometry and material data of FRP composite riser in ACP (Pre) were transferred to the module of transient structural. This is utilized to simulate the dynamic response of the risers. The gravitational acceleration, top-tension force, buoyancy, the end support condition (simple support at the top end, the fixed support at the bottom end) and the fluid solid interface were employed here and the time history of displacement and global stress distribution of the risers can be obtained as well in this module. The 12 elements/m in the axial direction and 20 elements in the hoop direction were employed for the risers. In the though-thickness direction, 18 elements, 22 elements and 1 element were employed for risers 1, 2 and 3, respectively ( Figure 1a). The number of elements in though-thickness direction was set automatically as the number of layers of the structure from the ACP (Pre) module.
The simulation of flow conditions was achieved in the module of fluid flow (Fluent) with the input of flow density, kinematic viscosity coefficient, flow velocity, turbulent intensity, etc. The flow dimension was set as 20D (X-direction) × 10D (Y-direction) × water depth (Z-direction), where D is the outside diameter of the risers, X, Y and Z directions are in-line flow direction, the cross flow direction and the depth direction, respectively. The large eddy simulation (LES) with bounded second order implicit transient formulation, which reduces the computational cost with acceptable accuracy, was chosen as the mathematical model for turbulence in computational fluid dynamics (CFD). The LES was used here to give the numerical solution to the Navier-Stokes equations for turbulent flow. This allows the largest and most important scales of the turbulence to be resolved and is suitable for the riser with medium length [27], which allows larger time-steps and coarser mesh sizes than direct numerical simulation (DNS) but needs more computational costs than the Reynolds-averaged Navier-Stokes equations (RANS). Also, the fluid grid of 5D × 5D area around the riser was refined and set as the dynamic mesh zone with smoothing and re-meshing methods. In this study, the number of elements for fluid zone varied from 0.56 million to 2 million depending on the mesh size and flow zone size ( Figure 1b). The module of system coupling was employed to achieve the data transfer between the transient structural and fluid flow, i.e., the two-way fluid structure interaction (FSI).
This analysis procedure and models selected have been proved by previous studies [34,35]

Results and Discussion
In this section, the VIV responses of all 9 study cases are presented and discussed, and the Grey relational grade of the 5 parameters (Etension, Ebending, L/Douter, tension force and current velocity) on the amplitude of VIV in cross flow direction for the risers is determined.

Natural Frequencies
When the natural frequency of a riser is close to the vortex shedding frequency, lock-in might occur, which would lead to fatigue failure of the riser. Using Equation (3), the vortex shedding frequency can be calculated.
where, St is Strouhal number, D is the outside diameter of riser and U is the flow velocity. According to the relationship between Re and St [36], the St are 0.2 for study cases 1, 4 and 7, 0.21 for study cases 2, 5 and 8 and 0.22 for study cases 3, 6, and 9.

Results and Discussion
In this section, the VIV responses of all 9 study cases are presented and discussed, and the Grey relational grade of the 5 parameters (E tension , E bending , L/D outer , tension force and current velocity) on the amplitude of VIV in cross flow direction for the risers is determined.

Natural Frequencies
When the natural frequency of a riser is close to the vortex shedding frequency, lock-in might occur, which would lead to fatigue failure of the riser. Using Equation (3), the vortex shedding frequency can be calculated.
where, S t is Strouhal number, D is the outside diameter of riser and U is the flow velocity. According to the relationship between Re and S t [36], the S t are 0.2 for study cases 1, 4 and 7, 0.21 for study cases 2, 5 and 8 and 0.22 for study cases 3, 6, and 9.
Also, the possibility of a riser's lock-in can be estimated using reduced velocity U r (Equation (4)). Lock-in may occur, when U r ranges 4-8 [36].
The natural frequencies, vortex shedding frequencies, and reduced velocities are listed in Table 5. From Table 5, it can be seen that the natural frequencies of risers 1 and 3 in cases 3 and 9 are close to their vortex shedding frequencies, and similarly the reduced velocity U r in these cases are close to the range of 4-8. Hence, lock-in might occur for these two cases. However, Equation (3) for vortex shedding frequency and Equation (4) for reduced velocity are based on the rigid cylinder without deformation and fluid structure interaction.

Time History of Displacements and von Mises Stresses for All Study Cases
The time history of displacements in the in-line flow direction (X) and cross flow direction (Y) for all study cases are plotted in Figures 2 and 3, respectively. For all the study cases, the results of a typical location, i.e., the middle point of each riser, are presented. Also, the possibility of a riser's lock-in can be estimated using reduced velocity Ur (Equation (4)). Lock-in may occur, when Ur ranges 4-8 [36].
The natural frequencies, vortex shedding frequencies, and reduced velocities are listed in Table  5. From Table 5, it can be seen that the natural frequencies of risers 1 and 3 in cases 3 and 9 are close to their vortex shedding frequencies, and similarly the reduced velocity Ur in these cases are close to the range of 4-8. Hence, lock-in might occur for these two cases. However, Equation (3) for vortex shedding frequency and Equation (4) for reduced velocity are based on the rigid cylinder without deformation and fluid structure interaction.

Time History of Displacements and von Mises Stresses for All Study Cases
The time history of displacements in the in-line flow direction (X) and cross flow direction (Y) for all study cases are plotted in Figures 2 and 3, respectively. For all the study cases, the results of a typical location, i.e., the middle point of each riser, are presented.   As can be seen in Figure 3, the cross-flow direction (Y) displacements in all cases demonstrate continuous vibrations. Similarly to the displacements in the in-line flow direction, the amplitude of the middle point of riser 1 in case 3 in the cross flow direction was the maximum (9.01 × 10 −2 m). Then, the amplitude of riser 2 in case 5 was 1.28 × 10 −2 m. For all the other cases, the maximum amplitudes of the middle points of risers were well below 7.1 × 10 −3 m.
The time history of von Mises stress for all study cases is plotted in Figure 4. For all the study cases, the results of a typical location, i.e., the fixed bottom end of each riser, are shown. Note that, for cases 1 and 7, which are dominated by tension forces, the von Mises stresses are not the maximum at the bottom end.  As can be seen in Figure 3, the cross-flow direction (Y) displacements in all cases demonstrate continuous vibrations. Similarly to the displacements in the in-line flow direction, the amplitude of the middle point of riser 1 in case 3 in the cross flow direction was the maximum (9.01 × 10 −2 m). Then, the amplitude of riser 2 in case 5 was 1.28 × 10 −2 m. For all the other cases, the maximum amplitudes of the middle points of risers were well below 7.1 × 10 −3 m.
The time history of von Mises stress for all study cases is plotted in Figure 4. For all the study cases, the results of a typical location, i.e., the fixed bottom end of each riser, are shown. Note that, for cases 1 and 7, which are dominated by tension forces, the von Mises stresses are not the maximum at the bottom end.  As can be seen in Figure 3, the cross-flow direction (Y) displacements in all cases demonstrate continuous vibrations. Similarly to the displacements in the in-line flow direction, the amplitude of the middle point of riser 1 in case 3 in the cross flow direction was the maximum (9.01 × 10 −2 m). Then, the amplitude of riser 2 in case 5 was 1.28 × 10 −2 m. For all the other cases, the maximum amplitudes of the middle points of risers were well below 7.1 × 10 −3 m.
The time history of von Mises stress for all study cases is plotted in Figure 4. For all the study cases, the results of a typical location, i.e., the fixed bottom end of each riser, are shown. Note that, for cases 1 and 7, which are dominated by tension forces, the von Mises stresses are not the maximum at the bottom end.   Figure 4 indicates that the time history of von Mises stresses has the similar change trend to displacement in the cross flow direction for each case. In case 3, the maximum von Mises stress of the bottom end of the riser 1 was 1.3 × 10 8 Pa, followed by riser 3 in case 9 (3.07 × 10 7 Pa). In cases 5 and 2, the maximum von Mises stresses were approximately 2.37 × 10 7 Pa and 2.24 × 10 7 Pa, respectively. For all the other cases, the maximum von Mises stresses were below 8.1 × 10 6 Pa.

Maximum Displacements and Von Mises Stresses for All Study Cases
Similarly, the time history of displacements in the in-line flow direction (X), the cross flow direction (Y) and von Mises stresses at all locations of risers can be obtained. Based on these data, the maximum displacements and stresses at different locations of risers are presented in Figures 5-9.  Figure 4 indicates that the time history of von Mises stresses has the similar change trend to displacement in the cross flow direction for each case. In case 3, the maximum von Mises stress of the bottom end of the riser 1 was 1.3 × 10 8 Pa, followed by riser 3 in case 9 (3.07 × 10 7 Pa). In cases 5 and 2, the maximum von Mises stresses were approximately 2.37 × 10 7 Pa and 2.24 × 10 7 Pa, respectively. For all the other cases, the maximum von Mises stresses were below 8.1 × 10 6 Pa.

Maximum Displacements and Von Mises Stresses for All Study Cases
Similarly, the time history of displacements in the in-line flow direction (X),the cross flow direction (Y) and von Mises stresses at all locations of risers can be obtained. Based on these data, the maximum displacements and stresses at different locations of risers are presented in Figures 5-9.   As can be seen in Figure 6, it is also clear that the maximum amplitude in the cross flow direction (Y) occurs above the middle point in all the study cases, which indicates the simple support (top end) leads to larger displacement compared with the fixed support (bottom end). According to the procedure of creating Figure 6, the maximum amplitudes in the cross flow direction and the location they occur for study cases 1-9 were 7.01 × 10 Comparing the results in Figures 5 and 6, the in-line displacement is much larger than the cross flow displacement. Based on the time when the maximum in-line displacement occurs for all cases, Figure 7 presents the total displacement distributions of all cases at those specified times. More specifically, the maximum total displacements for cases 1-9 were 2.  As can be seen in Figure 6, it is also clear that the maximum amplitude in the cross flow direction (Y) occurs above the middle point in all the study cases, which indicates the simple support (top end) leads to larger displacement compared with the fixed support (bottom end). According to the procedure of creating Figure 6, the maximum amplitudes in the cross flow direction and the location they occur for study cases 1-9 were 7.01 ×    Comparing the results in Figures 5 and 6, the in-line displacement is much larger than the cross flow displacement. Based on the time when the maximum in-line displacement occurs for all cases, It is noted that only the failure of steel risers (cases 7, 8 and 9) can be verified by the von Mises stresses. Based on the American Bureau of Shipping standard [37], 67% of the yield stress of X80 steel (371 MPa) is the allowable stress. Therefore, the steel risers studied in this paper would not fail for all cases 7 (6.6 MPa), 8 (3.5 MPa) and 9 (30.7 MPa). In contrast, von Mises stresses of FRP composite risers 1 and 2 cannot be utilized to verify the failure because the global stresses cannot represent the stress distributions in each composite laminae. The global stress distributions of composite risers 1 and 2 in Figure 8 are only used to determine the locations and time of the maximum stresses and these data can be used to verify the stress failure in each composite layer in Section 3.4.
For all study cases, the maximum von Mises stresses occur at bottom end of all risers, except in case 1 (−2.5 m) and case 7 (top end). In the middle parts above the middle point are the locations where second maximum von Mises stresses happen. In the point which is approximately 25% of the riser length to the bottom end, the minimum von Mises stresses occur. This is because the fixed support at the bottom leads to maximum stress of the bending situation, and the middle parts have the largest deformation which also results in large stresses. In case 1, the tension force and vertical deformation together lead to maximum stresses at approximately −2.5 m, while in case 7, the tension force is dominant and the maximum stresses occur at the top end. For the part which is approximately 25% of the riser length to the bottom end, a comparatively small deformation happens and also relatively far from the fixed support, therefore, the minimum von Mises stresses occur. It is noted that only the failure of steel risers (cases 7, 8 and 9) can be verified by the von Mises stresses. Based on the American Bureau of Shipping standard [37], 67% of the yield stress of X80 steel (371 MPa) is the allowable stress. Therefore, the steel risers studied in this paper would not fail for all cases 7 (6.6 MPa), 8 (3.5 MPa) and 9 (30.7 MPa). In contrast, von Mises stresses of FRP composite risers 1 and 2 cannot be utilized to verify the failure because the global stresses cannot represent the stress distributions in each composite laminae. The global stress distributions of composite risers 1 and 2 in Figure 8 are only used to determine the locations and time of the maximum stresses and these data can be used to verify the stress failure in each composite layer in Section 3.4.  For all study cases, the maximum von Mises stresses occur at bottom end of all risers, except in case 1 (−2.5 m) and case 7 (top end). In the middle parts above the middle point are the locations where second maximum von Mises stresses happen. In the point which is approximately 25% of the riser length to the bottom end, the minimum von Mises stresses occur. This is because the fixed support at the bottom leads to maximum stress of the bending situation, and the middle parts have the largest deformation which also results in large stresses. In case 1, the tension force and vertical deformation together lead to maximum stresses at approximately −2.5 m, while in case 7, the tension force is dominant and the maximum stresses occur at the top end. For the part which is approximately 25% of the riser length to the bottom end, a comparatively small deformation happens and also relatively far from the fixed support, therefore, the minimum von Mises stresses occur.
Based on the time when the maximum von Mises stresses occur for all cases, Figure 9 presents the von Mises stress distributions of all cases at those specified times. More specifically, the maximum von Mises stresses for cases 1-9 were 2.68 MPa at −2.

Lock-In Phenomenon for Study Case 3
As mentioned in Section 3.1, the lock-in likely occurs for study case 3. Further, it is confirmed by riser 1 in case 3 which has the largest vibration amplitudes and stresses. In this section, the fast Fourier transform (FFT) is utilized to find out the vortex shedding frequency and the vibration frequency of riser 1 in case 3 to analyze the lock-in" phenomenon. Figure 10 shows the movement track of riser 1 at 18.8 m under the water surface in case 3. From Figure 2, the stable vibration occurs after 10 s (the displacement in flow direction is comparatively stable at 0.44 m (Figure 5), while the maximum amplitudes in cross flow direction is about 0.094 m (Figure 6)), and the dominant vibration frequency 0.996 Hz of riser 1 can be obtained using FFT (see Figure 11).

Lock-in Phenomenon for Study Case 3
As mentioned in Section 3.1, the lock-in likely occurs for study case 3. Further, it is confirmed by riser 1 in case 3 which has the largest vibration amplitudes and stresses. In this section, the fast Fourier transform (FFT) is utilized to find out the vortex shedding frequency and the vibration frequency of riser 1 in case 3 to analyze the lock-in" phenomenon. Figure 10 shows the movement track of riser 1 at 18.8 m under the water surface in case 3. From Figure 2, the stable vibration occurs after 10 s (the displacement in flow direction is comparatively stable at 0.44 m (Figure 5), while the maximum amplitudes in cross flow direction is about 0.094 m (Figure 6)), and the dominant vibration frequency 0.996 Hz of riser 1 can be obtained using FFT (see Figure 11).     Similarly, the vortex shedding frequency (0.996 Hz) can be calculated using FFT as well. The result is shown in Figure 12. The data used here for FFT is collected using the magnitude of the velocity at 18.8 m under the water surface in the cylinder wake. Similarly, the vortex shedding frequency (0.996 Hz) can be calculated using FFT as well. The result is shown in Figure 12. The data used here for FFT is collected using the magnitude of the velocity at 18.8 m under the water surface in the cylinder wake. From Table 5 and Figures 11 and 12, the natural frequency of riser 1 in case 3 is 0.94 Hz. The vibration frequency of riser 1 in case 3 is 0.996 Hz, which equals the vortex shedding frequency and thereby indicates the occurrence of the lock-in phenomenon. Here, it is noted that the vortex shedding frequency obtained using FFT is different from that in Table 5, because in Table 5 the vortex shedding frequency calculated by Equation (3) is based on the rigid cylinder without deformation and fluid structure interaction.

Stress Distributions in Every Layers for Study Case 3
For each case, there are two significant instances in time, i.e., the times at which the maximum From Table 5 and Figures 11 and 12, the natural frequency of riser 1 in case 3 is 0.94 Hz. The vibration frequency of riser 1 in case 3 is 0.996 Hz, which equals the vortex shedding frequency and thereby indicates the occurrence of the lock-in phenomenon. Here, it is noted that the vortex shedding frequency obtained using FFT is different from that in Table 5, because in Table 5 the vortex shedding frequency calculated by Equation (3) is based on the rigid cylinder without deformation and fluid structure interaction.

Stress Distributions in Every Layers for Study Case 3
For each case, there are two significant instances in time, i.e., the times at which the maximum stress and the maximum displacement occur. Comparing the stresses in FRP composite layers in these two time points, the maximum stresses in every lamina of FRP composite risers (risers 1 and 2) are obtained and the factor of safety (FS), which is defined as the allowable strength/stress, are illustrated in Figure 13. In this paper, the maximum stress failure criterion [38] considers failures in the fiber direction, transverse direction and in-plane shear separately for all composite layers. The minimum FS required are 1.53 for composite layers and 1.68 for titanium liners [39].

Grey Relational Analysis of Multiple Parameters
The parameters considered in this paper that affect the VIV amplitude in cross flow direction of the risers are Etension, Ebending, L/Douter, tension force and current velocity, listed in Table 6. It is noted that, for a composite laminate, there can be a significant difference between its effective moduli in tension and bending (Etension and Ebending). In this paper, Etension is calculated based on the 3D effective properties of the composite tube using the 3D laminate property theory [40,41]. Ebending is evaluated using static analyses of the FEA models of the selected lay-ups with Solid186 (layered brick) under bending situations. More specifically, the FEA model for calculating Ebending is a cantilever pipe (30 m) under a transverse force of 1000 N and therefore, the ≈ ∆ , where, P is the transverse force, L is the pipe length, I is the moment of inertia, and ∆ is tip displacement of the pipe.  , respectively. Therefore, in case 3, the in-plane transverse stresses are the most critical stresses even though they are still far away from the failure stresses.

Grey Relational Analysis of Multiple Parameters
The parameters considered in this paper that affect the VIV amplitude in cross flow direction of the risers are E tension , E bending , L/D outer , tension force and current velocity, listed in Table 6. It is noted that, for a composite laminate, there can be a significant difference between its effective moduli in tension and bending (E tension and E bending ). In this paper, E tension is calculated based on the 3D effective properties of the composite tube using the 3D laminate property theory [40,41]. E bending is evaluated using static analyses of the FEA models of the selected lay-ups with Solid186 (layered brick) under bending situations. More specifically, the FEA model for calculating E bending is a cantilever pipe (30 m) under a transverse force of 1000 N and therefore, the E bending ≈ PL 3 3I∆ , where, P is the transverse force, L is the pipe length, I is the moment of inertia, and ∆ is tip displacement of the pipe. In this paper, the Grey relational analysis (Equations (5)- (12)) is utilized to determine the Grey relational grade of the 5 parameters on the amplitude of VIV in cross flow direction for the riser.
The detailed values for X 0 to X 5 are presented below: For X 0 to X 5 , the units are different and their values have a very large range. Hence, these sequences are normalized to comparability sequences using a process of Grey relational generation (Equation (6)). Here, 9 n=1 X i (n) is the mean value of X i (n), n ∈ (1,9).
, . . . . . . , From Equation (6), Y 0 to Y 5 are:  (7) is utilized to calculate the absolute difference between the original data and the comparative series.
where i ∈ [1,5], n ∈ [1,9]. The identification coefficient β ∈ (0, 1) is used to adjust the distinction between normalized reference series and comparative series. From the former study, the rank of the Grey relational grade cannot change with the change of β value and β only has an effect on the magnitude of the relational coefficient [42]. In this paper, β = 0.5, which leads to medium distinguishing effect and stability, is utilized [43]. (2) From Equation (11), Z 01 to Z 05 are:  (12) is used to calculate the Grey relational grade (GRG), which is defined as the numerical measure of the relevancy between the reference sequence (X 0 ) and the comparability sequence (X 1 to X 5 ) and indicates the degree of similarity between the comparability sequence and the reference sequence [44]. The higher the GRG of the comparability sequence of a parameter, the more similar this comparability sequence is to the reference sequence, i.e., the parameter from this comparability sequence affects reference sequence (VIV amplitude in crossflow direction) more significantly.

Conclusions
The study of vortex induced vibration of composite risers and steel riser is presented in this paper. Three types of risers, three current velocities and three water depths are considered which lead to five design parameters, i.e., E tension , E bending , L/D outer , tension force and current velocity. Based on OA sampling, nine study cases are conducted using CFD with coupled FSI to investigate risers' VIV characteristics including natural frequencies, global displacements, global stresses, and the stress distributions in each layer of the FRP composite risers. The effect of five parameters on the VIV amplitude in cross flow direction of the risers is analyzed using GRA and the Grey relational grade of these parameters is obtained.
The results show that: (1) Lock-in might occur for study cases 1 and 3 according to the relationship between natural frequency and vortex shedding frequency and the range of reduced velocity. (2) After an initial growth, the in-line flow displacement and von Mises stress maintain a relatively stable value and vibrate with a modest amplitude, while a significant vibration in the crossflow direction is obvious. (3) In general, deeper water depth and faster current velocities lead to a more severe VIV response. (4) Simple support leads to larger displacement compared with the fixed support.