DYNAMIC ANALYSIS OF A MARINE RISERByYongzhi GuoB.A.Sc. Taiyuan University of Technology, Taiyuan, China, 1982M.A.Sc. Taiyuan University of Technology, Taiyuan, China, 1985A THESIS SUBMITTED IN PARTIAL FULFILLMENT OFTHE REQUIREMENTS FOR THE DEGREE OFMASTER OF APPLIED SCIENCEinTHE FACULTY OF GRADUATE STUDIESDEPARTMENT OF MECHANICAL ENGINEERINGWe accept this thesis as conformingto the required standardTHE UNIVERSITY OF BRITISH COLUMBIADecember 1992Yongzhi Guo, 1992In presenting this thesis in partial fulfilment of the requirements for an advanceddegree at the University of British Columbia, I agree that the Library shall make itfreely available for reference and study. I further agree that permission for extensivecopying of this thesis for scholarly purposes may be granted by the head of mydepartment or by his or her representatives. It is understood that copying orpublication of this thesis for financial gain shall not be allowed without my writtenpermission.(Signature)Department of mechanical engineeringThe University of British ColumbiaVancouver, CanadaDate^C .^2 CI^q c1 zDE-6 (2/88)ABSTRACTThe problem of the dynamic analysis of a marine riser is considered in this thesis.The equations of motion for a marine riser undergoing large three dimensional de-flections and rotations are derived using modal discretization approach. The nonlinearexpressions for the bending and torsional deformations are obtained for a beam undergo-ing finite rotations. The elastic strain energy is expressed in terms of these deformationsand the Lagrange equations are used to derive the equations of motion valid for largedeformations. The three dimensional shear effect based on nonlinear elastic theory isincluded in the formulation.The dynamic response of a typical marine riser system due to excitations of the oceanwaves, current and surface vessel motion is studied. The hydrodynamic loading on theriser is represented by a general form of the Morison equation. The nonlinearity ofstructure and drag force are completely retained and a numerical procedure in the timedomain is developed. The effects of structural and hydrodynamic parameters on motionamplitude in two dimensional case are evaluated.Alternatively, the dynamic response of the simplified riser system is studied by ananalytical approach based on Butenin's method. The Butenin's method involves lesscomputational effort and the solution based on this method matches the numerical solu-tions well when the nonlinearities are small.iiCONTENTSABSTRACT^ iiLIST OF FIGURES^ viLIST OF SYMBOLS ixACKNOWLEDGEMENT^ xii1 INTRODUCTION 11.1 Preliminary Remarks ^ 11.2 A Brief Review of the Relevant Literature ^ 31.3 Scope of the Present Study ^ 82 FORMULATION OF THE PROBLEM 102.1 Model and Reference Coordinates ^ 112.2 Kinematics of the System 202.3 Shear Effects in Three Dimension ^ 242.4 Discretization Methodology 332.4.1^Kinetic energy ^ 332.4.2^Potential energy 362.4.3^Lagrange equations of motion ^ 402.4.4^Generalized force ^ 462.5 Particular Case: Planar Motion of an Eulerian Beam with Tip Mass 502.6 Modal Functions of the System ^ 53iii2.6.1 Modal functions of an Eulerian beam ^532.6.2 Modal function for a Timoshenko beam ^ 573 PARAMETRIC NUMERICAL ANALYSIS 633.1 Computational Considerations ^633.2 Effect of System Parameters ^673.2.1 Flexural rigidity ^673.2.2 Length of the riser ^693.2.3 Tip mass ^753.2.4 Ocean current ^753.2.5 Ocean waves ^803.3 Initial Conditions ^934 CLOSED-FORM SOLUTION USING THE BUTENIN METHOD 1044.1 Simplified Equations of Motion ^ 1044.2 Butenin's Approach-A Variation of Parameters Method ^ 1084.2.1 The eigenvalues and the solutions corresponding to linear system^1084.2.2 Solutions corresponding to nonlinear system ^ 1124.3 Typical Response Results: Comparison Between Analytical and NumericalData ^ 1195 CONCLUDING REMARKS^ 1235.1 Summary of Results 1235.2 Recommendation for Future Study ^ 125BIBLIOGRAPHY^ 127APPENDICES 130ivA COEFFICIENTS OF ANALYTICAL SOLUTION^ 130LIST OF FIGURES1-1 A schematic diagram of a drilling riser ^21-2 A typical flexible riser production system. ^42-1 Mechanic model ^ 122-2 A schematic representation of the undeformed and deformed beam . .^132-3 Eulerian angles /./.,, 9 and q5 ^152-4 The coordinate relationship ^162-5 The coordinate relationship ^172-6 The coordinate relationship ^182-7 Vector relations between inertial coordinate I and rotary coordinate e^212-8 The element change its orientation in space due to bending ^252-9 Shear angles in the plane consisting of e l and e3 ^262-10 Shear angles in the plane consisting of e 2 and e3 ^282-11 Shear angle relations ^292-12 Curvature components of deformed beam ^ 372-13 Definition sketch for linear wave theory ^482-14 Simplified model ^513-1 The trend of the first and the second natural frequency ^693-2 Change flexural rigidity—the first resonance case ^703-3 Change flexural rigidity—FFT analysis of the first resonance case^713-4 Change flexural rigidity—the second resonance case ^723-5 Change flexural rigidity—FFT analysis of the second resonance case^73vi3-6 The trend of the first and the second natural frequency ^743-7 Change length of the riser—the first resonance case ^763-8 Change length of the riser—FFT analysis of the first resonance case . .^773-9 Change length of the riser—the second resonance case ^783-10 Change length of the riser—FFT analysis of the second resonance case^793-11 The trend of the first and the second natural frequency ^803-12 Change tip mass—the first resonance case ^813-13 Change tip mass—FFT analysis of the first resonance case ^823-14 Change tip mass—the second resonance case ^ 833-15 Change tip mass—FFT analysis of the second resonance case ^ 843-16 Change the amplitude of current(negative)—the first resonance case .^853-17 Change the amplitude of current(negative)—FFT analysis ^ 863-18 Change the amplitude of current(negative)—the second resonance case^873-19 Change the amplitude of current(negative)—FFT analysis ^ 883-20 Change the amplitude of current(positive)—the first resonance case . .^893-21 Change the amplitude of current(positive)—FFT analysis ^ 903-22 Change the amplitude of current(positive)—the second resonance case .^913-23 Change the amplitude of current(positive)—FFT analysis ^ 923-24 Change the wave height—the first resonance case ^ 943-25 Change the wave height—FFT analysis of the first resonance case^953-26 Change the wave height—the second resonance case ^ 963-27 Change wave height—FFT analysis of the second resonance case ^ 973-28 Change the wave frequency ^ 983-29 Change the initial conditions—the wave frequency ^ 1003-30 Change the initial conditions—the first resonance case 1013-31 Change the initial conditions—the second resonance case ^ 102vii3-32 Change the initial conditions—the wave frequency ^ 1034-1 Comparison of numerical solutions and analytical solutions ^ 1204-2 Comparison of numerical solutions and analytical solutions ^ 122viiiLIST OF SYMBOLSa 13^Coefficients of nonlinear terms.c Wave speed.e or (e i ,e 2 ,e3 )^Rotary coordinate basic vectors.f^Frequency of oscillation(Hz).fi, f2^The first and second natural frequencies of oscillation(Hz).fx) fy In-line forces exerted on riser.k^Wave number.ki , k2^Eigenvalues of the system.1 Wave length.n i , n 2^Dimensionless natural frequencies.nx , ny , nz^Direction normals.qUI qt/1 qW1^ Generalized coordinates.t^Time.u, v, w^Deformations of beam along the directions of x, y, z, respectively.u', v', w'^Derivatives of u, v and w, respectively.up , vp , wp^Deformations of platform along the directions of x, y, z, respectively.vu,^The velocity of ocean waves.vc^The velocity of ocean current.x, y, z Inertial coordinates.xi, m , z i^Translational coordinates.ixA^Area of the beam cross section.A 1 , A:, B1 ,^Coefficients of integrations.CD^Drag coefficient.Damping coefficient of the system.CM^Inertia coefficient.Cs, Cy) C.^Force coefficients based on Froude-Kryloy theory.Dt^Riser pipe inner diameter.Do^Riser pipe outer diameter.D, Sphere diameter.E Yong's modulus of the beam.Fx , Fy , F,^Wave forces exerted on platform.Shear modulus.H Wave height.I or (4, i 2 , i3 )^Inertial coordinate basic vector.II or (i 1 ,i 2 ,i3 )^Translational coordinate basic vector./Cc , Ian , /cc^Inertia moments of the riser../&^Inertia moments of the platform.K 1 , K2^Numerical factors of the beam cross section.K„ Modal stiffness of the system.• Kn , Kc^Bending curvatures and torsional curvature of the riser.L Length of the beam.M12^Modal mass of the system.Mp^Mass of the platform./4^Mc^Bending moments and torsional moment of the riser.P Transform matrix.Qui, Qvi, Q.,^Generalized forces.T^Wave period.T j Fluid Kinetic energy in the riser.Tp^Kinetic energy of the platform.T, Kinetic energy of the riser.U1 , 172 , W,^Modal functions corresponding to u,v and w, respectively.al, a l^Complex modes.01, 02^Shear angles.Axial strain.A i , A2^Dimensionless damping coefficients.w Frequency of Oscillation(rad/sec).Dimensionless coefficient.w4, Lon , wC^Angular velocity of coordinate e w.r.t. coordinate I.0, 0 , 0^Eulerian angles.ei,^Modal functions corresponding to 0, 0 and 0, respectively.p Mass density of the beam.Dimensionless time.Rotary coordinates.4.1 4. 2^Complex modal coordinates.xiACKNOWLEDGEMENTI would like to express my gratitude to my supervisor Dr. A.S. Atadan for his inter-ests and valuable suggestions. I also want to express my appreciation to co-supervisorsDr. S.M. Calisal and Dr. V.J. Modi for their helpful suggestions, comments and dis-cussions. The investigation was supported by the Natural Sciences and EngineeringResearch Council of Canada, Grant No. A-2181 and OGPIN013, held by Dr. Modi andDr. Atadan, respectively.xiiChapter 1INTRODUCTION1.1 Preliminary RemarksThe interest in offshore floating platforms, for oil exploration and production facilitieshas increased during recent years. This is due to the fact that petroleum exploration andproduction is progressively moving into deeper waters.The marine riser is an important structural subsystem of any offshore platform usedfor oil and gas recovery. A riser is essentially a pipe which connects the wellhead atthe seabed to a fixed, a floating platform, or a vessel. It may be categorized as eithera production riser or a drilling riser. In the case of a simple riser (Fig. 1-1), a steelpipe extends from the ball joint at the top of the blowout-preventer to the surface ves-sel. Additional components of the riser system may include guide lines, choke-and -killlines, and buoyancy and tensioning devices. When the riser forms part of the producttransportation system, it may comprise of several pipelines (`import', 'export', test andcontrol lines). The system may have bundled or individual risers with separate tensioningdevices.As exploration and production moves into deeper waters, application of flexible risershas started to receive considerable attention[1]. In 1978, the first flexible riser was success-fully introduced at the Enchova Field in Brazil, and subsequently, in 1979, in Spain. Sincethen, flexible risers have found application in various parts of the world. More recently,1DrillingPlatform^/TensioningSystemLower Ball JointBlowout PreventerWell HeadSea FloorFigure 1 - 1: A schematic diagram of a drilling riser.2several companies have proposed floating production system concepts using flexible ris-ers for future oil and gas production in the North sea. Flexible riser systems have twomajor advantages over the conventional "straight" steel risers[2][4] in that they allow formore mobile positioning of the surface platform or vessel. Moreover, the flexible risersdo not require the expensive fixation and heave compensation systems that conventionalpretensioned risers must have. Figure 1-2 shows a flexible riser configuration.Regardless of the type of riser configuration under consideration, its design wouldinvolve a slender, structural member, in general a beam, in the ocean environment. Theanalysis, if developed with due consideration, can be equally applicable to a flexiblejacket platform, individual tendons of a tension leg platform, mooring lines of a catenaryplatform, submarine pipelines, etc.The risers used for transportation of fluids from the ocean floor to vessel mountedequipment is of particular importance for floating production systems. A productionriser system must be able to withstand severe storms(as per specification) in order tominimize the need for emergency riser disconnection. Thus the risers must have sufficientstructural strength to survive extreme environmental conditions. This calls for accuratedynamical analysis to assess response in the present complex and diverse loading. To laya sound foundation for such structural analysis is the subject of this study.1.2 A Brief Review of the Relevant LiteratureThe riser analysis has attracted considerable attention partly due to the fact that itis, and will continue to be, a vital link between the floating platform and the subsea borehole. Furthermore, the analysis poses several challenging problems.Marine analysis of the riser system may be static or dynamic. The static analysis isgenerally concerned with the maximum riser response in a vertical plane and does not3ProductionPlatform1M •SubsurfaceBuoySea Floor•^ MN\ 1 MA ^1■NNNWRNMAN^ NNNNNNNNNNNFigure 1-2: A typical flexible riser production system.4take into consideration the time-varying effects of waves, vessel motion and the inertia ofthe system. The vessel motion (one of the most dominant factors as far as the bendingstress is concerned) and the waves and current are advanced at suitable time intervalsand the maximum values of the critical parameters (position, stress, tension, etc.) arecalculated over an appropriate period (e.g., wave period). Such static analysis can beaccomplished in several different ways: finite difference approach, direct integration usinga fourth-order Runge-Kutta method, assumed deflection shapes of an elastic catenary orpower series, and finite element method.Obviously, a dynamic analysis of a marine riser system would be more complex. Thisis due to the fact that it would require information about the variation of the kinematicsof the waves and current as a function of x, y, z, and t, In addition to the rig motion.Furthermore the riser may move in various directions in the x-y plane at various elevationsdue to the omnidirectionality of the waves and current.The first step in marine riser dynamic analysis involves mathematical modelling of thesystem and the evaluation of the hydrodynamic loading. The mathematical simulationshould account for the geometric nonlinearities caused by large rotations, nonlinear termsin the strain-displacement relations and suitably defined boundary conditions. Physicalnonlinearities, due to the relations between the components of stress and strain, is gener-ally not considered because in the riser application the material remains elastic. Irani[29],and Vogel and Natvig[11] have presented models for studying large three-dimensional mo-tion of compliant risers based on the Eulerian beam theory.The system modelling also requires the consideration of wave and current forces dur-ing severe as well as nominal sea states, accounting for water depth, the rig motion, risergeometry, etc. The waves and current may give rise to in-line as well as transverse forces,and the dynamical response which is not necessarily planar due to the omnidirectionalityof the surface and internal waves, and the current over the length of the riser. Generally5speaking, the hydrodynamic loads on the marine riser are calculated by using Morison'sequation, modified to account for the velocities and acceleration of the flexible pipe. Acareful choice of the drag coefficient has to be made (the inertia coefficient is usuallynot of importance[3]). The selected drag coefficient should reflect a degree of uncer-tainty in Reynolds and Kvelegan-Carpenter numbers along the riser, possible transverseoscillations due to vortex shedding, and relative fluid velocity in the wave zone.Once the mathematical representation of the system is established, the next logicalstep is to select an appropriate solution procedure.Available methodologies may be classified into two major categories: frequency andtime domain analyses.Typically, frequency domain methods are based on a linearized drag load and struc-tural representation. They are quite efficient with respect to computational time andeffort, however they are, in general, approximate. Various refinements have been intro-duced to incorporate effects of finite wave amplitude and instantaneous load, but theseextensions are based on linearized loading and apply only to the regular wave case. Theseprocedures are widely used for the analysis of drilling risers and offshore platforms withonly limited nonlinearities, like those posed by nonlinear hydrodynamic forces. Theiraccuracy may be unacceptable in problems involving large nonlinearities or large changesin structural geometry. In such situations, one has to resort to numerical integration ofthe nonlinear governing equations to assess the structural response.Time domain methods on the other hand, allow for a nonlinear drag force description,and they are based on a stepwise calculation of the riser response. They may also includea fully nonlinear representation of the riser stiffness properties. At each simulated step,the stiffness terms can be recalculated according to the instantaneous response on thebasis of the nonlinear theory of elasticity. This implies that the displacement dependentstiffness terms are updated during the analysis. Obviously, the time domain methods6usually require considerable amount of computer time. This particularly true for thehighly nonliner representation of the system.Approach to the time domain analysis may include the finite difference or Galerkintechniques to solve the governing differential equations. They are, however, often limitedby the assumptions made to simplify the differential equations. Each of the methods[19][21] neglects torsion, and either ignores elongation or utilizes Lagrange multipliersto determine the riser tension separately as an auxiliary problem. These assumptions areintroduced to uncouple the equations in the two transverse principle directions. However,such simplification is not valid in problems involving redundant structures, applied mo-ments, longitudinal vibration and elements with noncircular cross sections. When used,the Lagrange multipliers lead to additional degrees of freedom and increase the numberand bandwidth of the algebraic equations to be solved at each time step. Furthermore,application of the finite difference techniques is frequently limited to structures withsimple end conditions.Another approach to the time domain analysis, which is more extensively used, isthe finite element method. Although finite element based methods are capable of mod-eling three-dimensional structures, including the effects of torsion and elongation, it hasbeen frequently applied to the two-dimensional problems. The suitability of this method[23][24][25] for highly nonlinear problems is also limited by formulations which utilize con-stant stiffness, damping and mass matrices. Of the finite element based methods citedabove, only one [22] appears suitable for general highly nonlinear structures. The methodis only briefly discussed in [22]. Irani[29] retained full nonlinearity of the hydrodynamicloading but left the rest of the equations of motion in a linearized form.71.3 Scope of the Present StudyThe objective of this thesis is to develop tools that are useful for the design and anal-ysis of marine risers connected floating platform system. A relatively general character ofthe Lagrangian formulation makes it applicable to a wide range of system configurations.The three dimensional shear effects based on nonlinear elastic theory are formulatedand the Lagrange equation is used to derive fully nonlinear equations of motion.Chapter 2 presents details of the formulations and modal discretization approach.The Lagrangian formulation of the equations of motion and the finite element methodare surprisingly similar: they need the establishment of kinetic and potential energies ofthe system, need shape functions to discretize the continuous system, and use Lagrangeprocedure to establish the ordinary governing differential equations of motion. However,they have two major differences: Finite element methods use interpolation polynomialas a shape function, while modal discretization method use admissible functions fordiscretization. The shape function of the finite element method is defined piecewisealong the beam, but the admissible function is expressed over the entire domain.The dynamic response of a two-dimensional riser-platform system due to excitationsof ocean waves and current is studied numerically in Chapter 3. The hydrodynamicloading is represented by a general form of the Morison equation. Nonlinearities of thestructure and drag force are completely account for and a numerical procedure in the timedomain developed. The objective is to assess the effects of structural and hydrodynamicparameters on the system response. The results obtained are prerequisite for efficientdesign of this class of system.Finally, the dynamic response of a simplified riser system is obtained analytically,through a closed form solution, based on the Butenin method. The thesis ends with asummary of results and some thoughts concerning direction of the future research which8may prove fruitful.9Chapter 2FORMULATION OF THE PROBLEMThis chapter begins with a discussion on the coordinate system, followed by thedetermination of the kinematics of the system. The kinetic and potential energies of thesystem are then derived. Finally, using Lagrange's formulation, the equations of motionfor a riser which undergoes three large dimensional deformation are determined.The mathematical modelling of the structure requires derivation of the equations fora tensioned beam undergoing large three-dimensional deflections. In order to developthe equations of motion for a tensioned beam undergoing large deformations, the gen-eral nonlinearity theory of elasticity needs to be considered. The nonlinearities are bothgeometrical and physical in origin. Geometric nonlinearity is associated with the neces-sity to consider the deformed configuration and the need to include nonlinear terms inthe strain-displacement relations. Physical nonlinearity is due to the relations betweenthe components of stress and strain being nonlinear. In riser applications the materialremains elastic so it is adequate to consider geometric nonlinearities only[29].Literature review reveals that the formulation developed so far is restricted to theassumptions inherent in the Euler beam theory which does not including the sheareffects[5][14][18]. Derivation of Timoshenko's beam is far more complex than that ofEuler's beam. Here a general derivation is presented, for the first time, which considersnot only the internal flow and arbitrary formed platform, but also the shear effects basedon three dimensional nonlinear elastic theory. Finally, the hydrodynamic loading due toocean waves and current is evaluated using the Morison equation.102.1 Model and Reference CoordinatesThe marine riser in its simplest form is a pipe which connects a wellhead at theseabed to a fixed or floating platform or a vessel. It may be a production riser or adrilling riser. In the case of production riser, the internal flow, i.e. gas or oil, can betransported simultanneously to another vessel which is used to transport the gas or oilfrom sea production to inland. Thus, the mass of the platform under such a considerationwill remain constant despite of the internal flow from sea bed.The corresponding mechanic model is a cantilever beam with internal flow, which isconnected to a platform and subject to wave and current forces, see Fig. 2-1.There are three coordinates used in this paper. They are: inertial coordinate ( o,x, y, z) with basic vector I or (i 1 , i 2 , i3 ), translational coordinate ( c, xi, yi,^) withbasic vector I I or (^) and rotary coordinate (c, 6,77, C ) with basic vector e or(e 1 , e 2 , e3 ).All deformations of the beam are referred to the inertial coordinate system x, y, z. Thez axis is taken to be aligned with the centroidal axis of the undeformed beam. When thebeam deforms, the centroidal axis at an arbitrary section translates linearly by amountsu, v, w in the x, y, z directions, respectively, and the section twists by an angle aboutthe centroidal axis.The coordinate axes x o , yo , 20 are fixed at an arbitrary point on the centroidal axis ofthe undeformed beam. Before deformation, x o , yo , zo are parallel to x, y, z, respectively.Deformations u, v, w and displace the x o , yo , zo to x 1 , y i , 2 1 and rotate x 1 , y i , 2 1 towhere the axis c is tangential to the deformed axis (see Fig. 2-2). The orientationof c with respect to x 1 , y i , z i may be expressed in terms of Eulerian angles. Fig. 2-3shows the Euler angles 0, e and where the rotations are taken in the following order:a positive rotation b about the y i axis resulting in x 2 , y 2 , z2;1 1WavesC.G.SurgeRollPitch pro:/SwayInternal FlowCurrentElementCoordinateSystemGlobalCoordinateSystemSea FloorHeaveYawFigure 2-1: Mechanic model12Figure 2-2: A schematic representation of the undeformed and deformed beam13a positive rotation B about the x 2 axis resulting in x3 , y3 , z3 ;a positive rotation about the z 3 axis resulting in 4, 7/ ;The relationship between the translational coordinates and the rotational coordinatescan be expressed explicitly by mathematic formulations.From the theory of linear algebra, we know that if i 1 , i 2 , i3 are coordinate basic vectorsin three dimensional space, any other vector can be expressed as a linear combination ofthe coordinate unit vectors:r^i + X2 i2 X3i3If e l , e 2 , e3 are another coordinate unit vector different from i 1 , i 2 , i3 , there must exista linear relationship between these two coordinate unit vectors:e = PIor by matrix notation:e1 P11 P12 P13e 2 P21 P22 P23 i2e3 P31 P32 P33 13where the matrix:P11 P12 P13P P21 P22 P23P31 P32 P33is called transform matrix representing the relationship between two coordinate systems.As showed in Fig. 2-4, a positive rotation IP of coordinate x l , y i , z 1 with unit vector1 2^Theabout the y, axis results coordinates x 2 , y 2 , z 2 with unit vectors i, 1 2 1 23relationship between the two coordinate systems is:1422 , 23Y1 / Y2Figure 2-3: Eulerian angles ii), 8 and 015z 0xFigure 2-4: The coordinate relationshipX2 co.s0^0^— sin0 X1Y2 0^1^0 Y1 (2.1)Z2 siroP^0^cos0 Z1or in particular,i2i cos ^0^— sin0 1 11; 3 = 0^1^0 l• 21 (2.2)•213 simb^0^cos0- -1•31Similarly, a positive rotation 8 of coordinates x 2 , y2 , z 2 with unit vector ii , iz , i3 aboutthe x 2 axis results coordinates x 3 , y3 , z3 with unit vector 4. , i3, i3. From Fig. 2-5, wehave,16ZX2 , X30/Figure 2-5: The coordinate relationshipS3 1^0^0 X2y3= 0^co38^sine Y2 (2.3)Z3 0^—sin°^cos8 Z2Or,•31 1 1^0^0 21 1•31 2= 0^cos8^sine •21 2 (2.4)•31 3 0^—sine^cosi) 21 3Finally, a positive rotation 0 of coordinates 5 3 , y3 , z3 with unit vector 4, i 23 , il about thez3 axis results coordinates 4", ri, C with unit vector e l , e 2 , e3 . Fig. 2-6 indicates:cosh^sincb^0 S3= —sin0^cosh^0 Y3 (2.5)0^0^1 Z3yx71C17z3 ,zx30x yFigure 2-6: The coordinate relationshipcosq sing—sin0 cosq'0^00^i30^1 231^i33_^- Or ;ele2e3Combination of equations 2.2, 2.4 and 2.6 results:e l cosq5^sing^0 1^0^0 cosi')^0^— sin011^-e2 —sink^cosq5^0 0^cost)^sine 0^1^0 12e3 0^0^1 0^—sine^cost) sinik^0^cos0 i cosOcos11) sinOsinesin0 cosesin0 —sinOcos0 cosOsinesin0—cosOsin0 sinIksinecos0 cosecos0 sinOsin0 cosOsinecosq5sinOcose^—sine^cosOcos8or,e le2e318Equation 2.7 can be expressed concisely:where,cosOcos0e = PI 1sinOsinesinck^cosesin0^—simbcos0 cosOsinesing5(2.8)P —coszksin0 simbsinecos0^cosecos0^sinOsin0 cosOinecos0 (2.9)sinikose^—sine^coslk cos()Transformation matrix P represents the relationship between translational and rotationalcoordinates and it has two properties:1. P -1 = PT (2.10)that is P an orthogonal matrix.0—we 0Wn —we 02.^PPT = Si= (2.11)where, WTIWC0.cose.sin0 Ocosch0.cose.cos0 — Osin0—(2.12){w} or Si is an angular velocity matrix. It represents the angular velocity of coordinatee with respect to inertial coordinate I (or translational coordinate I ' ) at arbitrary timet.192.2 Kinematics of the SystemFigure 2-7 shows an element taken from beam or platform. The position of arbitrarypoint in the element at arbitrary moment can be expressed as:r = rc + rewhere rc denotes the position of centroid of the element and rt denotes the position ofarbitrary point of the element in coordinate e:^=^ne2 Ce 3Because the coordinate (,77, 0 is fixed to the body of the element, we get,4^d( 0dt^dt^dtthe velocity of arbitrary point in an element is:^= + = +^+^+ Ce3 = + { 6 }T ewhere,rc = it(z, t)i i + 1)(z, t)i 2^zi)(z, t)i3foT= H, 77,e le263From equation 2.8, we have6 =^+ P11Sinced= dt = 01120Figure 2-7: Vector relations between inertial coordinate I and rotary coordinate ewe obtainé = h i^(2.13)Substitution of the inverse of equation 2.8II , p-i einto equation 2.13 gives rise to,e = PP -l eAccording to the second property of P,e = PP -l e = PPTe= ftewhere ft is given by equation 2.11. Hence,ic = ei + net + Ce3 = {We = {Wile21or, rE = (w1 — wo)e i + (wg — we C)e2 + (we n — wn0e 3 =el e 2 e3wC (.4.41e^7l^C= wxr er = r, + re = r, + w x rewhere,w = we + wil e2 + wce3re =^+77e2 + (e 3The kinetic energy of an infinitesimal mass separated from the element can be writtenas:dT = 1/2i.idmThen the kinetic energy of the element becomesT = 1/2fAf= 1/21m (i, + i e )(ic + te )dm= 1/2fm i(Cc c +^+^ (2.14)where,I22T = 1121_me(u2 + 7) 2 +12,2) 4_ 1/2{w} T .[Iie.{1.0}2(2.15)is the relative velocity of centroid with respect to centroid. Namely,Ve = 0following which equation 2.14 becomes,T = 1/2fm (ic .ic te .it )dmwhere,= ii 2 (z,t)^1) 2 (z t)^11) 2 (z, t)and,it .4^(wn.0^wC. 11) 2^(wC.e^we.C) 2^(we. 71^wn.0 2Finally,(71 2 + C 2 )^— e.n( 2 + C 2 )^— p.0_^—C.71^(e 2 + 77 2 )where {w} represents the angular velocity of coordinate e with respect to the inertialcoordinate I and its expression is given by 2.12,— 11;.sin8= I dm111(2.16)We{w} = wn0.case.sin0 Bcos00.cos0.cos0 — Bsin023[l]e = IM(7) 2 + C 2 )n6•(—6.71(6 2^( 2 )—6.(7/•((62 + n 2)dm (2.17)Equation 215 is the kinetic energy expression derived for any element. The integrations2.16 and 2.17 are carried out on the entire element. It can be used both for platform andriser.2.3 Shear Effects in Three DimensionIt is well known that the shear angle in two dimensional problem is determined by,(2.18)where is the slope due to bending and 1/) ]. is the slope of center line of an element whichundergoes both bending and shear effects.^is determined by the following equations :1,/) ]. = (2.19)In space, due to bendings in two directions and torsion, the relationship betweenshear angle and displacement becomes very complicated.Let us take out an infinitesimal element with rectangular shape from the elastic body.The central line segment of the element MN is a vector in space. Its magnitude MN issupposed to be dz in length before deformation. We fix the rotational coordinateon to the element and make the axis C tangential to the central line MN and parallelto the z axis before deformation. After deformation, the element displaces to a newposition with displacement u, v, w in the direction of x, y, z, respectively. Meanwhile,the length of the element will change from Al N to M* N* and its orientation in space dueto bending, but its rectangular shape will not change if the shear effect is neglected[16]as showed in Fig. 2-8.2411Figure 2-8: The element change its orientation in space due to bendingbut keep its rectangular shapeIn this case, the direction of e3 is the direction of M*N* and the projection of M*N*in the direction of el and e 2 are zero,cos(e3,M*N*) 1cos(e i ,M*N*) = 0cos(e2,M*N*) 0V. V. Novozhilov stated in his "Foundations of the Nonlinear Theory of Elasticity" [16]that the element would not keep its rectangular form after deformation and there wouldbe nonzero projections of M*N* in the directions of e l and e 2 if the shear effect isconsidered.In the plane consisting of el and e3 (see Fig. 2-9), ab will be perpendicular to ac ifthe shear effect is neglected. If the shear effect is taken into account, a will translate toa*, and c will translate to c*. a*c* will not be perpendicular to a*b and there will be an25zIIeibYe2e3ieiFigure 2-9: Shear angles in the plane consisting of e l and e 326angle (32 between e3 and e31 where e31 denotes the direction of the projection of centerline 1\4"1■1* in this plane.Similarly, in the plane comprised of e2 and e3, there will be an shear angle 13 1 due toshear effect as showed in Fig. 2-10.Suppose the projections of line IVI*N* in the directions of e l , e2 and e3 are p i .dz ,p2 .dz and p3 .dz respectively. Then,pi.dz = M* N* .cos(M*N* , e l )p2 .dz = APN*.cos(M*N*, e2)p3 .dz = M* N* .cos(M*N*, e3)It is clear,10- *11* == ds = 69 1 2 -I- 292 2 dhp3 2 . dz (2.20)Figure 2-11 shows the relationship between vector 1\1*N*, el , e2 ,e3 ,e31 and e32, wherethe segment M*N*3 1 is the projection of \I*N* in the plane consisting of el and e3 :M* N*31 = \APi.dz) 2 + (P3.dz) 2 = '10)0 2 + (133 ) 2 .dz02 is determined by one of formulas :(2.21) pi .dzsin/32p? + p3.dz^ +P3p3 .dzcos/32 =+ pi.dzP3N/Pi27Y...ze2e3e32e2Figure 2-10: Shear angles in the plane consisting of e 2 and e328e 3^e32Figure 2-11: Shear angle relationstan/32^/p3Similarly, the projection of M*N* in the e2, e3 plane is the segment M*N*32.M* N*32 OP2 • dz )2 (p3 dz )2 Op2 )2 (p3 )2 . dz01 is determined by one of formulas:sin/j1 =p2 .dz^p3 p3.dzP2^+ Pip3 .dzp3.dzP3 Jpi +29tanI32 = — p2/p3In order to get the expressions of p 1 , p 2 and p3 in terms of displacement u, v, w, wewill make use of the relationship between e and I.We have shown that an element which is dz in length before deformation becomesds=M*N* in length after deformation. If the projections of M*N* are dx*, dy* and dz*in the coordinate of x 1 , y i and z 1 respectively, Then,dx*^= u' .dzdy*^v'.dzdz* = (1 w').dzand,M*N* dy*i 21^dz*i i3+ dzi 21 H- (1 (2.22)where u', v' and w' are derivatives of the displacements u, v and w, respectively. Werewrite equation 2.22 as,M*N* = {dX*} TI i^(2.23)where,{dX*} T = [ dx*, dy*, dz*1 1=i31By substituting the relationship e = PI 1 into equation 2.23 , we have,M*N* = {dX*}T Pe^ (2.24)30or,^M*N*^(Pli.dx* + P12 .dy*-F P13.dz * )ei= (P2i.dx* +• P22.dy* +• P23 .dz*)e2= (P3i.dx * + P32 .dy* P33 .dz*)e3^(2.25)or,M*N* =^P12.1/ + P13.( 1 w 1 )1.dz.ei= [P21.u i + P22.v' + P23 .(1^w')] .dz.e2= [P3i .u 1 + P32.v' + P33 . 1 + w 1 )].dz.e3Formula 2.25 expresses the projection of line M*N* in the direction of e l , e2 and e3respectively. Comparing 2.23 with 2.25, we find,^p l =^/312.vi^P13.(1 + w')]P2 = [P21•v i + P22•v' + P23.(1 + w')]p3 = [P3 i .u 1 + P32.v ' + P33 .( 1 + w')]pl P11 P12 P13P2 P21 P22 P23 V (2.26)P3 P31 P32 P33 (1 + w')in which P13 is given by 2.9, i.e.P11 P12 P13 cosOcos0 + sinOsinOsincb cosOsin0 —sinOcos0 + cosiksin8sin0P21 P22 P23 —cosOsin0 + sinOsinOcos0^cosBcos0 simPsin0^cosOsinOcosqSP31 P32 P33 sinOcose —sin() CO371COSThe shear angles are determined by:41 1 =tan -1 ( P2 )P3(2.27)31132 = tan -1 (^ ) (2.28)P3Equations 2.26, 2.9, 2.27 and 2.28 express the relationships among shear angles, Eulerangles and the displacements in three dimensional problem.As an example, we use the formulas in two dimensional problem, where only thebending and shear effect in the plane of xoz are considered. Thus,= 0= 0v' = 0= 0p 1 = coszk.u l — sin0p3 = sin0^f cos,L'note that,t an 82 P1P3SiTLIP - coszlsin0 .u' cos0sinlkcosi)sty* u , + 1cost',= Jamb ian — tan21,1tan,32 = + tan' I, .tanzP 1tan(0 — 0 1 )that isQ2= 1/1- 01The three dimensional shear angles are reduced to the shear angles in a plane.322.4 Discretization MethodologyThe first step in riser analysis is the derivation of governing equations of motion fora pipe connecting to a platform. The pipe is an elastic body with infinite degree offreedom in general. In order to get a numerical model calculable by computer, a properdiscretization method should be selected. The modal discretization method is one of thesimple and effective methods and therefore is used here.2.4.1 Kinetic energyThe expression of kinetic energy (Eq. 2.15) derived for any elements is,1T = 2— e (ii 2^th2)^1/2{w} .[i] e .{W}where,Me = I dmAl= IM(71 2 + C 2 )^-.77-01^( 2 +^-77.0-(.77^( 2 + 71 2 )dm For a platform,Ale= [ 1p]{w} = {wp}^{4.,)(L,u = uP = L , t)= p = 1)(L , t)33tb = tbp^11)(L, t)The expression of kinetic energy for the platform becomes,Tp = 2Mp(lip2 75p 2 '142) 1/2{WP}TIIPLIWPIwhere,{wp} =For the riser,wCp WE(L,t) Op.COSOp.SiTLOp + O.COSOP= =col7P Wr)(L,t) Op .SinOpOp . COSOp .COSOp —wCp WC(L,t) 0. 19^0p.SiThOPm[/], =m and [I] are the mass and inertia moments of unit length of the riser respectively. Thekinetic energy of the riser is,Tr = 1/2J L [m(u 2 + 1) 2^2L72) --I- {W} T . [I].{w}]dz^T0where, {w} =WEwri7,b.cos6.sin0 9cos0,b.cosO.cosq — 9sin0— OsineTf is the fluid kinetic energy in the riser[29],T f = 1/21 0 [p.A(iiu,'.U`) 2 p.A(i)^.U*) 2^p.A1.1` 2 ]dzThe principal axes can be selected so that,ItE000Inn03400ISM(2.29)[4] =Ice000Inn000Icctherefore,T.tipl-{wp} 2 I UP (1,1)P .cos8p .sinOp 6C.cos4) 221• — / 13 (74.COSOp .COSOp — dp .Si7/4) 221-/PC (0 -^.Sin0p)2C12- {w}^'{W}12-/ce (0.cosO.sin0^O.cos.0) 2 1^•• Inn (7b.cosB.cos0 - iLsin0) 21• 2 Icc( - 7/).sin0) 2 (2.30)Because of the axial symmetry of the pipe,/ec =^= /equation 2.30 can be reduced to a simple form:1-{w}T^1.[ 1 1 . { w }^I ( 11)2 CO 32 + 8 2 )+ 2 1cc(:k - .sin0) 2The total kinetic energy of the system is,T = Tpm ti 2 + I; 2^II) 2 \2 pkp^p p)-/P(721 U ,bP .cos8P .sin0P dp .cos0p ) 2352gii (Op .cos8p .cosOp — Op .sin0p ) 22-1 /1‘). c (Op — Op .sin0p ) 2_1 I L ^+ 7)2^w. 2 )/(74 2 cos 2 0^e 2 )^/cc (ciS — 1P.sin0) 2^1.1 1 .U*) 2^p.A(i)^71.U*) 2^p.A.U* 2 ]dz2.4.2 Potential energyWe have derived the general form of the angular velocity of coordinate e with respectto inertial coordinate 2.12 as, {w} = can0.cos8..sinqS^()coscb0.cos8.cos0 — Osin0—(2.31)wCAssume that the origin of the ,^ moves along the deformed centroidal axis ofthe beam with unit velocity such that the^77 and C axes are always directed alongthe principal torsion-flexure axes of the beam[15]. Kirchoff's analogy states that theinstantaneous angular velocity components of the ,71, C system with respect to the x, y,z system correspond to the components of curvature of the deformed axes[15]. Figure2-12 illustrates this principle, where the space derivative with respect to the curvilinearcoordinate, s, along the deformed centroidal axis corresponds to the time derivative forthe angular velocities.Thus, the bending curvatures^KT, and the torsional curvature IQ obtained by36Y 1 ) Y2Z3 5 CFigure 2-12: Curvature components of deformed beam37projecting aeas 7 as and 2 along the 4 . , i and C have a similar form as formula 2.31,as^as+ - cos0.COSO.COSO — sin0asaA' — a- sin()as^asKinKcThe bending strain energy and torsion strain energy can be obtained by one of thefollowing forms:L M2 M2 m2U1 + U2 = 1/2 o (^ H- ^ ).dz E J4 E Jr, GJcOr,Ul + U2 = 1/2 (E.J.1q EJn .K,72 GBecause of the axial symmetry of pipe we get,= 1,1 JThus the bending strain energy has a form:21 foL (Ele .lq-12L foL EJ(K1 1C2 ).dz21 foL EJ {(7:s cos0) 2 + ( a:se ) 2 1.dzThe shear energy is given by,U3 = 1/2f (K 1 AG13 K2AGOD.dzwhere K1 and K2 are numerical factor depending on the shape of the cross section, A isthe cross-sectional area, and G the modulus of rigidity. i3 and /32 are determined by,tan -1 (- 131 )P3tan -1 (—E1 )P338P1P2P3P11P21P31P12P22P32P13P23P33 (1V+ w')andPi t P12 P13 cosOcos0 siniPsinesincb cos8sin0 — sinOcos0 cosOsinesin0P21 P22 P23 —cosiksin0 sinOsin8cosck^cosecosq5 siniksin0 cosiksinecosOP31 P32 P33 sinOcose —sine cosOcoseThe tension strain energy of the system is,T2U4 = 1 IL ^ .dz2 o EAor,LU4 = 1 fJo EA.E2.dzwhere the axial strain c may be expressed by a nonlinear form in terms of derivatives ofdisplacements u, v and w [16][29],1+ 2- (1.1 /2 + V/2 + 21) /2The potential energy caused by internal pressure is,1 12 oU5 = -^Tou l2 dzThe total energy of the system is,U = + U2 + U3 + U4 + U5392.4.3 Lagrange equations of motionThe energy expressions obtained are:Kinetic Energy:T = Tp TT1 M^2 ti 2^2)2 Pk PIP^u) 2'^Pl+ 2143j7,bp.COS9p.,S7T/Op + 9p .0O390) 21+ - 1-fm (14.COSOp .COSOp — OP.sin0p ) 21.ZC)C ( -^. S in Op ) 2+ -1 I L [M(11 2 + 1) 2 + 2// 2 )o./(ti) 2 cos 2 0^2 ) + ./cc(4 — IP.sin0) 2u'.U*) 2^p.A(i) v' .U') 2 p.A.U"JdzPotential Energy :U^Ul + U2 + U3 + U4 + U51 1^alk^2 89— 0 [EJ{( as cose) + (—as )2 }DO a0+ Gfc(^assznO) 2▪ K i AGO K2AGP▪ E A.{w'^(u,12^w'2)}2+ ToU 12 1dzin whicha-s = {u12 + 1/2 + ( 1 + w ' ) 2 } - 2 . 5Z40We transform the natural coordinates into modal coordinates as follows:u(z, t)^E Ui (z).qui(t)it(Z, t)^E ut (z).4„,(i)u/(z,i)^E u:(z)• qui(t)v(z, t)^E 17,(z).q,„(t)t=i1.)(Z, t)^E vi (z).4,,(t)2=1vi(z,^E2=1W(Z, t)^E lirt (z)•q„„(t)E wi (z).4„„(t)W I (Z, t)^E lic(z)•qw,(t)i=1ip(z,^E c(z). qui(t)P(Z, t) = E 4),(z).4„,(t)4111) 1 (z, t) = E Taz)• q„,(t)8(z, t) = E z ( z). qvi (t)B(z ) t) = E o z ( z).4,, i (t)i=1n8' (z , t) = E 0:(z). qvi(t)2.1Ch(z, t) = E 4),(z)• goi (t)0(z, t)^E cp i ( z).4„,(t)n01 (Z 2 t) = E 4:(z). q4„(t)i=1where Ui (z), 171 (z), Wi (z), 111 1 (z), 0,(z) and 1 1 (z) are modal functions of the systemand quz, Qvi, qwi and q,, are modal coordinates of the system. According to Lagrangeequations in matrix form:d aT^aT^audt^ „} ) (9{g,,}^a^ = {Q.}we obtain the equations of motion in qu-coordinate:Miiip{U(L)} If.,„01,{T(L)}+ I Km pA)ii{U} /„,2 1 .4{111}1dz+^[E J 1P 1{W'} (T0 — p AU " )u 1 {U^K AG (11) — u 1 )({T} — {U'})] dz42foLp.A.U*(i.Z{U} — U.{U'})dz.• 11)4 [(cosepcosOp • ck p — sinOpsint9p • O p)(4) pcosOpsinOp + pcos0p)+ cosBp s in0p(OpcosOp s inOp + OpcosOpsinOp • ckp- 7419psin8psinOp + iipcos Op — tjp cbpsin(pp )101(L ) 1• InPn [—( i)psin9pcos cbp — Opcos8psin0p) (1,GpcosOpcosOp — i)ps in0p)cos8pcosOp(OpcosOpcos q5p — 1,bpOp sinOpcosOp- OpOp cosOpsinOp — OpsinOp — Op Opcos0p) — Ikp ] Pli(L )1- IcPc [pcosOp(Op — Opsin9p) + sinOp( Op — OpsinOp — Opt9pcosep)H4(L)}+ fo [I(1/)cos 2 0 — Osin28 —-^ — ti4cosO)sin.0 -+ &:os6(:k — 11;sin0))d{xli}dzL ' AG . 131 r D[ 3 (simb.sinO.cos0 — coszk.sin0){11 1 }- Jo p2 p^2P3 (u'cosO.sine.cos0 + u /sinO.sin0• (coszk.sin0 — simk.sine.cos0)(1+ w')){W}- P2 sin.O.cos9{U'} — P2 (u'coslk.cosO — sinO.co.s0(1+ u1)){T}]dz-^1L K'AG.132 [P3 (cosO.cos0 + simk.sinO.sin0){1110 Pi 4- P32+ P3 {—u'sirrik.cos0 + u'cosO.sin0.sint9+ (1 + w')(—co.50.cos0 — sinik.sine.sin0)}{T}Pi simb.cosefU l l — Pi (u lcaslii.cos0 — sinO.cos8(1 + w i )){W}]clzJo K'AG(.4) — u')({W} — {U'})dzL^Vcos20+fo^u'2 + v'2 + (1 + .0) 2EJ[ ^ti]ftlilelzfL (uV2cos26 + Oa) EJ [u'2 + v'2 + ( 1 + wO2]2 nclz = {Q.} (2.32)43According tod aT^aT^au dt ( a{q.„} ) a{q,,} + O{q,} = { c2v}we obtain the equations of motion in qu -coordinate:Mp .iip{Vp} + Pidp{O(L)}foL[(m + p.A)ii.{V} + EJO'{0 1 } + IA01+ p.A.U*(1, 1 {V} — i).{1,1) + (To — p.A.U* 2 )v 1 .1171+ K'AG(0 — v 1 )({0} — {V'})].dz[(0pCOS OpSi714COS Op -^OpSin OpCOS Op + /6 9 0p CO30p COS2 OP•9 COS 2 -^sin0 cos0 —11)PP^P^P^2 cos0 sin0 sin 2 0- Op OpsinepsinOpcosOp — 16p OpcosOpsin 2 Op — Op OpsinOp cosOp — OpHO(L)}+ InPn klypcos0P cos0P sin0P — 9P76psin9Pcos0 sin0P — P lib cos0P sin 2 0P • 2- szn^— dp4cosq5p sin Op — 0.p16. p COS OpCOS 2 Op + Op (6psin Op COS Op• 2 COSOpCOS 2 Op + 11;p9pS2nOpCOS Op] { 0 ( }•+ /1&[OpIkpCOSOp - I,bp2 sin0pcos0p1{0(L)}EJ^L^(7.112cosOsin0 {0}dzfo^/L 12 + v'2 + (1 + w 1 ) 2-^^1L E J 01)12cos20 + 012)(v1)0^[u'2 +V'2 + ( 1 + ID 1 ) 2 ] 2 { 171}CIZrL 9'EJ[^u'2 + v'2 + (1 + w')2^0']10}dz1L K i AG.1310 293 +[ P3 cos0.cos0{171 } + P3 (u'sinO.cos0.cos0 — cosq5.sintiv`+ co.50.cosO.cos0(1+ w')){O} + P 2 sin0{11 + P2 (sin11).sin0•+ v icos0 + cos0.sin0(1+ w l )){0}1dz — K' AGO — v')({O} — {V'})dz+^K'AG132 [P3 (sinO.cos0{11 + P3 (sinip.cos0.sin0)74+ 7344- v'sinq.sin8 + cosVi.cose.sinck){0}Pi sinO{V} + Pi (sin? .sine-u' + vicos0 + cos? I) .sint9(1 + w')){O}]dz- {Qv}According tod aT^aT^au(^. ) = {Q.}dt Ofq„,}^0{qu,}^Ofq.}we obtain the equations of motion in qu -coordinate:rLMpii,p{Wp} + J rnin{W}dzL^(012 COS 2 9 + 0 12 )(W 1 )•^E [11 ,2 +^+ (1 + tv,)2]2 { vliv 7 }dzrL K I AG.131^[P3 (sirn k.sin0 cos71).sinO.cos0) — P2 cosO.cosO]flifildzA + AL^DsinO.cos0 cosl,b.sint9.sin0)— P i coszp.cost9RIIT'ldz10^+^113= {Q.}According tod aT^aT^au. ) ^ + ^ = {Qo}dt^{q4,1^fq,/,}^a fq,/,}we obtain the equations of motion in q4,-coordinate:/P ( 20P2 .cos 2 Op .sin2c¢p 2/1p .Op .cos9p .cos24 — 2 0PS27120p ){4)p }1 •^ 1 •inPn ( — -./4.cos 2 Op .sin2Op — Op .Op .cosBp .cos2Op^2-8p2 sin20p ){40 p }▪ 1-1&(Op 11) .SinOp Op.8p.cos9p){4)p}• I(c (0 — 1b.sine — ?)).e.cos0){4q.dzILG.Icc (01 — ?I/sine) ^{4)} ^(11 12 +^+ (1 + w1) 2 );-- dz45foL(—P3). IC' AG .tan -1 (—E-3-) .{u'(—cosl .cos0 — sinlk.sin9.sin0)Pi + Pi— 11 cos° .sin0 + (1 + w i )(sin11).cos0 — si rok .sinO .cos0)} {4)1 .dzjoi,^K'AG.tan-'(—P3) ,+^(—P3).^' .{u (—cosO.siri0 + sin' i I) .sine.cosO)Pi +P34- v i cost) .cos0 + (1 + w 1 )(si 7* .sin0 + cos'. sin 0.cos0)}.{(1)}.dz= {Q.}In the above equations, {Q u }, {Q v }, {C2,,,} and {Q 4,} are general forces in matrix formandT{U} = Ui (z) U2 (z)^. Un(z)T{V} = [ ii(z) 172 (z)^.^. .^17(z){W} = [ Wi(z) W2(z) Wn(z)W1(z) T 2 (z)^. tlin(z){O} [ 01(z) 0 2 (z)^. On(z) T1 4) } = I Cbi(z) 4)2(Z) .^(1),(Z)where "T" indicates the transpose of the matrix.2.4.4 Generalized forceThe in-line wave force on a slender offshore structure is usually computed using theMorison equation[6][12][29]. The Morison equation was developed by Morison, O'Brien,46Johnson, and Shaaf in describing the horizontal wave force acting on a vertical pile whichextends from the bottom through free surface. Morison proposed that the force exertedby unbroken surface waves on a vertical cylindrical pile which extends from the bottomthrough the free surface is composed of two components, inertia and drag:avu,f = Cm A I at^CD A D Ivto lv i,in which f is the force per unit length of the vertical cylinder and vw is the wave velocity.A 1 = p1D 2AD= 2pDvw^Hcosh(kz)= vu,(z,t) = V„,(z)coswt = Tsinh(kL)coswtavu,^2ir2Hcosh(kz) sznwtOt^T2sinh(kL)The velocity of wave is based on linear wave theory. Figure 2-13 shows the generalform of a wave train propagating in the X direction. H is the wave height, 1 is thewavelength, T is the wave period, c is the wave speed (c=//T), k is the wave numberrelated to the wavelength 1 by k = 27r/l, and w the wave frequency related to the waveperiod T by (.4) = 27r/T.The appropriate values of the inertia coefficient CM and the drag coefficient CD needto be used to predicate the wave force accurately. A vast library of experimental dataon CM and CD is available from numerous laboratory and field tests that allows us tochoose appropriate values for the coefficients [12].The Morison equation in its original form applies to a fixed structure. For a structurefree to oscillate in the presence of waves and current, the Morison equation is modifiedas:f = Cm^—^CDADIvw±vc — itl(v.±vc — it)47Wave speed cLXFigure 2-13: Definition sketch for linear wave theorywhere vc is the current velocity. The negative sign in front of v c represents the currentopposing the direction of water propagation. When used in two directions of x, and y,the equations have the following form[61:fx = Cm '^—^CD AD Ivwx +vcx — V I(v„,x ±vcx —fy CMAll)wy — CA Ali; + CD AD Hwy ±Vcy — '11 I (V ivy +Vcy —where vws , vwy , vcx , vcy , u and v are wave particle velocities in the direction of x and y,current velocities in the direction of x and y, and pipe displacements in the direction ofx and y, respectively.If we add the term CA Atii into the total mass of the structure, the Morison equationwill have a simpler form:fx — CMAliAux + CD AD I Vwx ±Vcx 1.11(Vwx+Vcx 1L)48fy = emit /1).y CD AD Iv.y +vc, — v j(v„,y+vsy — i))The Morison equation is used to evaluate the wave force on the small structure suchas the riser. When the length of the structure is of the same order of magnitude as thewave length, the Froude-Kryloy theory can be used to evaluate the wave force on thestructure such as platform[6].Based on Froude-Kryloy theory, the expression of dynamic wave pressure exerted onplatform is given by[6],Hcosh(kz)P Pg^cos(kx — Lot)2cosh(kL)The total force on the platform in a particular direction is obtained by integrating thecomponent of this pressure in the direction over the submerged portion of the structure.The expressions for the horizontal forces in the x and y directions and vertical forcecomponents in the z directions are written as:Fx = Cx IspnsdSFy = Cy I I spn y dSF2 = Cz f spn z dSin which Cx , Cy and C2 are force coefficients, and n s , nu , n, are direction normals in thex, y and z directions, respectively. The integration is carried out on the entire submergedsurface of the structure.The generalized force corresponding to modal coordinates are follows:QUI = j. Ui(z)[CmAii),,s + CDADIv..±vc. — iLl(v ius ±vcs —11,(L)Cs^sp-ns dS —^tit i ( z)0(C, I I spn z dS)dz^(2.33)49i = 1,2„n.LQ.,^Vi(z)rmArli. + CDAD iv„,y +vcv — 1)1(vv„,±vcy — 15)1dzV.(L)Cy fI spn9 dS — 10 0,(z)0(Cz f I spn z dS)dzi = 1,2„n.^Q,,,^Wz (L)(Cz I spn z dS)i = 1, 2, ^ , n.In the expressions of equations 2.33 and 2.34, the last terms:— f L W,(z)0(Cz f f pn z dS)dzand^—^01(z)0(Cz f f pn z dS)dz(2.34)(2.35)are geometric nonlinearities associated with vertical force[27][26]. kli i ( z ) and 0,(z) aremodal functions of angle b(t, z) and 9(t, z) respectively.2.5 Particular Case: Planar Motion of an Eulerian Beam with Tip MassConsider a simplified model of the general system as showed in Fig. 2-14.The system is reduced to a cantilever Euler-type beam with tip spherical mass, whichis subject to the transverse oscillation in the plane under the influence of wave and currentforce.Neglect the shear effects and the inertia moment and let= 0= 0= 0= 050Wavesintemal FlowGlobalCoordinateSystem' v V'Nv 'WV,/ ' 'V' s' WV" V' V" W Vv ‘kXSea FloorCurrentFigure 2-14: Simplified model51Lfoequation 2.32 will reduce to describe a two dimensional motion. For convenience, letU,(z) be replaced by Y(z) to denote the Eulerian beam modal functions. Noting thatT i (z) =17, 1 (z), we obtain:or in another form,Mpiip{Y(L)} I (7n + pA)ii{Y}dz+foL[EJulY"} + (To — pAU* 2 )u'llqdz^+^+^fLp.A.U*(ii'{Y} — ii,.{}"})dzfL^it"+ 1 u"]{}7"}dzlL^ILI un2o EJ [0 + 1]2 {1"}dz = {Qu } (2.36)MpiipY,(L)+ f (m pA)iiY(z)dzfoL[E Ju"Y",(z) + (To — pAU* 2 )72 1 Y, 1 (z)jdz10Lp.A.U*(inC(z) — ii.1":(z))dz u• EJ[^^uo +1^."]}7",(z)dzfoL EJ [u121U; 21]2 Yzi(Z)dZ = Qui (2.37)Substitutingi = 1,2, ^ , n.U(Z, t)^Y,(z).qui(t) =1:1C(z).qi (t)..1^i=1into the equations 2.36 and 2.37 and noting the orthogonal property of the modes:(m pA)ICYjdz AlpYi (L)Y3 (L) = 0 (i j)^(m pA)1C 2 dz MIC 2 (L)^(i j)fLo52[E/Ifi llYi" + (To — pAU* 2 )1/,'}'; ']dz = 0 (ij)[E../Y „2 + (To — pAU* 2 )11, i2]dz K„ = w, 2 M,,^= 3)yieldsE 4; [Mpri (L)Yi (L) + f (m + pA)1;(z)17,(z)dz] = M,,4 in^Lq; fE [EJY";(z)Y",(z) + (To — pAU* 2 )}71(z)Y:(z)]dz = Kiiqti=1 JoThe final two dimensional equations therefore become,M,24, K „,q, + Jo p.A.U*(UT„(z) — U.17: (z))dzufor,^”EJ[u'2 + 1 ulY",(z)dzUlniaEJ (^u/2 + 1)2 Yil(Z)dZ = Quawhere,172 ( z )[CAMP). + CD AD vwx +vcs — uI (vu, x ±vcs — U)1dzYs (L)Cs^sprixdS —^171,(z)ul (Cz Ifspn z dS)dzn(z, t) = E Y;(z). qz (t)i = 1, 2„ n.2.6 Modal Functions of the System2.6.1 Modal functions of an Eulerian beamThe linear conservative differential equations of the system corresponding to Eulerbeam is,(rn pA)ii (pAU* 2 — To )u" EJun" = 0fLJofor,LQui53with boundary conditions:U = U = 0^at z = 0and,(pAU* 2 — TOIL' + EJu m = Mp iiEJun = 0^at z = LThe differential equations and boundary condition of modal function is,En' (pAU• 2 — T0 )Y" — w 2 (m pA)Y = 0^(2.38)with boundary conditions:Y = = 0^at z 0(pAU• 2 — T0 )Y'^= —co 2 A/p YTiJY" = 0^at z LBy solving the differential equations, we get the modal function:Y(z) = AsinA i z BcosA i z CshA 2 z Dch ,\ 2 zwhere,(I k44a 4k4a 2=^2 -1-a22 +542 pAU* 2a = E J^To> 0if,a 2^a4= - -2 + —4 +k 4a2= 1a 2^a42 -+ ^4 + a4or,if,pAU*2a2 Tu^E^ >0 Jk4w2(m pA)E JApplying the boundary conditions to the modal function, we obtain:19 11 C + b i2 Db21 C b22D —where C and D are arbitrary constants, andbll^(pAU* 2 – To )(A 2 ch,A 2 L – ) 2 cos L)EJ (A 2 3 ch,A 2 L A 2 A l 2 cosA l L)M(shA 2 L – sinA l L)b12 = (pAU* 2 – To )(A 2 shA 2L A i sin) i L)+ EJ(7 2 3s1),2L – A l 3 sinA i L)55L10L10w2 Mp (chA 2 L — cos A i L)b21 = E J (A 2 2 shA 2L A 2 A i sin) i L)b22 = EJ() 2 2 chA 2 L A l 2 cosA i L)The natural frequency is determined frombii b22 — biz b21 =0and the modal function is given by,2Y(z) = D[(chA 2 z — cos A i z) — r(shA 2 z — —szn) l z)]b12^b22T^b11^b21A, and b,3 are specified with respect to natural frequency, and D is an arbitrary constant.The orthogonal conditions for this system are defined as:(1) Orthogonality with respect to mass:for mn,L(m p A)YCn Yndz Alp }',„(L)Y,(L) = 0for m = n,(m p A)Y,, 2 dz A/p Y„, 2 (L) = M„.,,(2) Orthogonality with respect to stiffness:for m#L,[E JY„"Ym n + (To — pAU* 2 )Yn l Y„, i]dz = 056for m n,L[EJY,," 2 + (To — pAU* 2 177, 12 )]clz = K„,,^wn. 2 1117-in2.6.2 Modal function for a Timoshenko beamThe linear conservative differential equations of the system corresponding to Timo-shenko's beam are:with boundary coand,(m pA)ii (PAU *2 — TOO K*(0' — u") = 0hi) — EJ7,b u d- K*(0 — u') =0nditions:u = = 0^at z = 0(pAU* 2 — To )u' K*(1,1) — u')^AlpiiE JO' = —^at z = Lwhere,K* K AGAfter eliminating 1/), we have,EJu""(1 PAU *2 — To) + (PAU*2 — To )u"K*I(m + pA) (94 u(m pA)ii ^If*^at'EJ(m + pA) I(pAU* 2 — To ) a4 u[I -I- = 0K*^K*^az25t257Alternatively, after eliminating u, we have,pAU* 2 — To )E JO""(1 ^K*o) + pAU* 2 07 I(m + pA) + (m + p 41) +K*^at4EJ(m + pA) I(pAU* 2 — To ) a4 11)— [I + ^K*K•^I^ 0az 2at 2This shows that u and 7b have the same natural frequencies. Let,u(z, t) = Y(z)sin(wt +0(z, t)^tli(z)sin(wt + 11))we obtain the modal relations of the system:w2(m pA)Y = (pAU* 2 — To )Y" + K*(' , — Y")^(2.39)^= —EPP" + K*(41 — Y 1 ) = 0^ (2.40)with boundary conditions:Y= =O^at z = 0^(pAU• 2 — To )Y 1 + K*(kli — Y') = —111pw 2 Y.^at z = LEJV Ipw 2 W^at z = LThe frequency equation is,d4Y 2d2 ^— k 4 Y" = 0dz 4 ±a dz 258The solution of the differential equations is,Y(z) = AsinA i z BcosA z Csh.) 2z Dch,A 2 zwhere,a 2^a42 + 4 + k4C14 k44a 2 =(pAU* 2 — T0) + (.0 2( i KJ.(m pA) gpAu•2-21) )>0K• Ef(i.^PAU" -71')K•k4 W2 (M pA)(1 —EJ(i^pAU• 2 — To \K • )or,a 2— 2 + a 4—4 +a 2^a42 1-^4 +k4+ pA )^.1(pAU ' 2 -TO)(pAU* 2 — TO+ w 2(I Eja2 =^ (7'11pAU• 2 — To^K•^1 EJ(l^COK•In order to determine the frequency we need explicit relationship between^and Y(z).From eq. 2.39 and eq. 2.40 one gets,EK* 2 EJw 2 (m pA)j y ,K•(K* — w 2 I)EJ(K* — pAU* 2^y„,K•(K• — w 2I)+ f21""a22if,—if,59from which we have,(z) = Ai(fi — f2 ) 1 2 )(AcosA i z — BsinA l z)+ A 2 (f1 + f2 A 2 2 ) ( Cch A 2 z + D s \ 2 z )Recall,Y(z) = AsinA i z + BcosA i z + Csh.\ 2 z + DchA 2 zApplying boundary condition's to these two equations we obtain:[ b1 ^bizb21 b22CD00The natural frequency is determined from,b11 b22 — b12 b21 = 0where,b11 2(PAU" — TO — K s )[A2ChA2L — A2 (fl f2A2)cosA L]+ K * A2(fi f2)2 2 )(chA2L — cosA i L](fl^f2Al2)+ Mpw2(shA2L A2(fi + f2A22) sinAiL)— f2Al 2 )b12 = (pAU* 2 — To — K*)(A 2 sli)t 2 L + )t i sinA i L)K*[.ki(f1 — f2)1 2 )sinA1L A2(f1^f2)2 2 )3hA2L]+ Mpw 2 (ch)) 2 L — cosA i L)b21 = )2(f1 f2A2 2 )[EJ(A2shA2L A i sin)1 i L)— /pw 2 (citA 2 L — cosA i L)]60b22 = Ei[Al 2 CA — f2 A l 2 )cosA 1 L + A2 2(fi f2A2 2 )chA2L]— /pw 2 [Ai(fi — i2)1 2 )sinA i L A2(fi i2A2 2 )shA2L1Finally, the modal functions of Timoshenko's beam are given by,Y(z) = D[(ch)1 2 z — c Aos2 Aui 1.2 + f2A22)sinA z)]) —L' 21 (shA 2 z^b11^Ai(A — f2A1-,)^i11/(z) = D[Ai(fi — f2)'1 2 )sznA1z^A2(A^.i2A2 2 )810t2zbllbiz A 2 (fi^f2A2 2 )(chA2z — casA i z)1It can be shown that for Timoshenko's Beam, there exist the following orthogonal con-ditions.(1) Orthogonal condition w.r.t. mass:for mn,[(m pA)Ym Yn ./xlim tli n ]dz + AlYm (L)1,(L)+ IpT,(L)Tn (L) = 0.for m =fo [(m pA)Yn 2 /tPn2 ]c/z A/p Yn 2 (L) /p ilin2 (L) = Ainn(2) Orthogonal condition w.r.t. stiffness:when m^n,fLo[K: (W m —^—^EJVmtliIn (T0 — pAU* 2 )1,17,1,:jdz = Uwhen m = n,o [IC*(xli n —YD 2 EJTC,2 + (To — pAU* 2 )V,i jdz^K,,,,61where,K„,, = Ainnwn 262Chapter 3PARAMETRIC NUMERICAL ANALYSIS3.1 Computational ConsiderationsConsider a set of first order differential equations,f(q,t)with the numerical solution at t=nA denoted by q„, where A is the step size used inthe numerical integration. Basically, there are two methods to obtain qn . They areone value and multi-value procedures (prediction-correction). The former uses only q n_ iwhile the latter uses k previous values, a ••••) qn-k) to determine q„. For complicateddifferential equations, the multi-value method gives more accurate and reliable resultsthan the one-value approach.A multi-value method consists of two processes: predictor and corrector. In theprediction process, an approximation to q„ denoted by q„,0 is obtained by an interpolationmethod. The predicted value is then refined by the corrector formula. If the qn , i valuesatisfies the error test with a user specified tolerance, then the calculation will terminateat qn =q„, i . Otherwise, the correction process is repeated until qn ,„„ satisfies the error testand hence qn=q„,,,. The order of a multi-value method refers to the complexity of thecorrection formula. In general, the higher the order, the more accurate is the solution atthe expense of a higher computing cost.Since the system's equations of motion are complicated, a multi-value method was63preferred. The IMSL:DGEAR subroutine was chosen. It is adapted from Gear's subrou-tine DIFSUB and has the following advantages:1. Accuracy. The subroutine is based on the implicit Adam's multi-value methodwith a variable order of up to twelve. This high order ensures the numerical solution tobe accurate.2. Ease of Programming. In the correction process, the Jacobian of the differentialequations is required; hence, it should be coded in by the user. In DEGEAR subroutine,there is an option of evaluating the Jacobian numerically.3. Economy. Based on the result of the error test at each integration step, thesubroutine automatically changes the order of the correction formula or the step size,if necessary. This feature makes the subroutine efficient and economical in terms ofcomputing cost.4. Additional feature. By changing the parameter METH of the subroutine, it canbe used to solve stiff equations based on Gear's method.In order to obtain a numerical solution, the equations of motion were first written asa set of first order differential equations. In the case of a two dimensional marine risersystem, the original equations are:Alia4,+ Kii q, +rL p.A.U*(iNi (z)— iL.IC /(z))dzIL^u„EJ[11'2 + 1^7"1 z ) dzL^U171112Jo EJ (u 12 + 1 )2/711 (z)dz = Q uzwhere,= f Yi(z) [ CmAli).. + CD ADIvu,s ±v.. —^— it)*0+ }CMGs J spnrdS —^11.11(z)1Z(C2 J f spn z dS)dzJQUI64=^Yi (z). qi (t)^= 1, 2, ^ ,n.If we consider only first two modes, the equations may be written as:+ C1242rL 1+ J Ej[ fqin(z)-E 42Y21(z)} 2 + 1^lliqiYill (z) -+^z)11' " 1 (z)dz- foLEJ [glY11(Z)^q2Y21(Z)Ilqini(Z)^q21-2"( )12 Yil (z)dz(1071'(z)i- q 2 12(z)] 2 + 1) 2foLYi (z)[Cm il i i)„,x• CDADIvws+vc. —^z^42 1 2( z )1{v..±vc. — iYi(z) — 4 2 Y2 (z)ydz+ Yi (L)Cx f spnidSfo L i(z)[41Y11 (z) + q2Y21 (z)j(C z f I s pr z z cIS)dzM2242 + K22q2 + C2141L 1o^E J {gi n (z) + 42}721 (z)} 2 + 1^11[41Y,' (z) q2 Y2"(z)11'" 2 (z)dzfoL E J[q i Yii(z) + q2n(z)li41Yi"(z) + q 2 17(z)1 2 ([07(z) + q2y2w 2 + 1)2 Y21 ( z)dz=^y2 (z+ CDADIvwx+vcs — 411"1(z) — 42172(z)l{vwx+vcs 4072(z) — 42172(z)}idz+ Y2 (L)C.^sprix dSL 1712(z)i4iYi(z) 42Y21 (z)1(C. fispnzdS)dzwhere,L=C1242^fo p. All* (1 .11'1 (z) — ii.1'7(z))dz.1.0L=^p. A IP [{4 1 }/11 (z) + 42}721 (z)Wi(z) — {4 1 } (z) + 42Y2(z)Wii(z)lciz65= q2 I 19 - 11- 111n(Z)Y1(Z) 172(Z)YAZ)IdZ0C12 = — C21for a sphere[5], we get,Cx f fspnz dS = Czpl' gHk.sin(wt)CZ f fspn z dS = Cz pVgHk.tanh(kL)cos(wt)in which V is the volume of the sphere. Letting ql = q l , q; = q2 , q3 = 4 1 and q: = 42)the original second order differential equations can be transformed to a set of four firstorder differential equations as follows(3.1)q2^q2 (3.2)m114;^—1(114i C124:L 1 fo EJ[fqlyii(z)+ qn'21 (z)} 2 f 1^l][q:Yn(z) + q2Y:(z)11/" 1 (z)dz—+ foL 7 [qIY11(Z)^qIY2(Z)][qIY111(Z)^q.2'17211(Z)12 1711 (Z)d.Zairl yii( z )^17 ,21, y2/( z )12 + 1) 2foLYi(ZYCMAIiitusCDADIVwx+Vcx q; 171(Z) (4 1 2( 2 )1{Vwx±Vcx — qPri(Z) — 17:172(2)YdZYi(L)Cx ff s lYnxdSY 1 1( 2 )[4TYAZ) q;IA.2)}(Cz f isPnzdS)d.Z(3.3)M2244 = — K224; — C21 1/3*66Lfo^1 EJ[{qTyli(z) + (4Y2'(z)} 2 + 1^1][ql}1 z) + q;Y2"(z)]1/11 2 (z)dz.1)L EJ [q1171(z) 4. q;172/(z)][qP71"(z) + qP72"(z)12 Y;(z)dz+ ([q11711 (z) + q1172'(z)1 2 + 1) 2JOL+ Y2(z)[Cm Ari),..+ CDADIvu,±vcx — Ui (z) — q:Y2 (z)l{v„,z +vcx — qP 2 (z) — q:}72 (z)Edz+ Y2 (L)Cz I I spnzdSJLO Y 12(z)[qIYAz) + q;n(z)](Cz 1 1 spnz dS)dzEquations 3.1 to 3.4 was solved by the IMSL: DGEAR subroutine with the toleranceset at 10 -8 .3.2 Effect of System ParametersIn designing an offshore structures, the dynamic behaviour of the riser system isof increasing importance. A parametric study is carried out to quantify the effects ofthe riser length, flexural rigidity, tip mass, the exciting frequency of the wave, and themagnitude of wave and current, using a computer program developed for the analysis ofthe riser system. The results of this investigation would be important for the designersin making a reliable design. Table 3.1 presents the input data used for parametric study,which are adopted from [29][6][12]. When one of the parameters, say flexural rigidity,is changed, other parameters are kept constant at the same values as shown Table 3.1during the sensitivity study.3.2.1 Flexural rigidityThe first and the second natural frequencies of the system for varying flexural rigidityare comparatively showed in Figure 3-1 With the increase of rigidity EJ, the natural(3.4)67Table 3.1 Data for parametric studyLength of the riser^ L=158mMass of the sphere Afp=87000kgModulus of elasticity of riser pipe^E=2.07X1011N/m2Riser pipe outer diameter^ Do=0.4064mRiser pipe inner diameter D,=0.3747mSphere diameter^ D,=3mDensity of pipe p1 = 7995kg/m 3Density of sea water^ P2 = 1025.18kg/rn 3Density of internal flow p3 = 919.5kg/rn 3Velocity of the internal flow^U*=7.32m/sDrag coefficient^ C D= 0.7Mass coefficient Cm=1.5Sphere horizontal force coefficient^Cx=1.5Sphere vertical force coefficient Cz=1.1Wave height^ H=5mWave Period T=263.8sVelocity of current at surface^1,c=0.126m/sInitial displacement at z 1 and z2^u(z1,0) = u(z 2 ,0) = 0.Initial velocity at z 1 and z2 = it(z2 ,0) = 0.68Figure 3-1: The trend of the first and the second natural frequencywhen the flexural rigidity is changedfrequencies increase greatly.The flexural rigidity has an important influence on the dynamic response of the sys-tem. Figures 3-2 and 3-4 describe the trends of first resonance and second resonancewith the change of the flexural rigidity. It is noticed that the large flexural rigidity EJincreases the stiffness of the riser system, which significantly decreases the displacementof the tip mass in the first resonance and second resonance cases.An FFT analysis is carried out to analyze the frequency response components inthe resonance case. Fig. 3-5 verifies that even in the second resonance case, the firstfrequency response component plays an important role at the initial phase.3.2.2 Length of the riser691.00.50.0-0.5-1.0 Time(s)00.1lobo^2obo^3u1:10-0.25 11 1 11 11111II1111111250 SOO 70 100 12'50 3 10 liTime(s)00U p0.250.00C1111 11111 11111 iu ill11111a b cEJ=0.8 x 10 8 Nm 2f1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzEJ = 0.8 x 109 Arm 2f1=0.012Hzf2=0.16Hzf=f1=0.012HzEJ = 0.2 x 10 11 Nm 2f1=0.06Hzf2=0.85Hzf=f1=0.06HzFigure 3-2: Change flexural rigidity—the first resonance case70a b cEJ=0.8 x 10 8 Nm 2f1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzEJ = 0.8 x 10 9 Nm 2f1=-0.012Hzf2=0.16Hzf=f1=0.012HzEJ = 0.2 x 10 11 Nm 2f1=-0.06Hzf2=0.85Hzf=f1=0.06HzFigure 3-3: Change flexural rigidity—FFT analysis of the first resonance case71a b cEJ = 0.8x 108 Nmf1=0.00379Hzf2=0.0542Hzf=f2=0.054HzEJ = 0.8 x109 Nm 2f1=0.012Hzf2=0.16Hzf=f2=0.16HzEJ = 0.2 x 10 11 Nm 2f1=0.06Hzf2=0.85Hzf=f2=0.85Hz0.000 up a-0.050-0.100-0.150-0.200-0.250-0.300II^lii^dii! 11 ,11^II^' IITime(s)586 10b0 15b0 20b0 2500 3000Figure 3-4: Change flexural rigidity—the second resonance case72UP0.1500.1000.050aUP0.03000.02500.02000.01500.01000.0050Freq.(Hz)00Ulth■—e s IS 1.1 es^•.S se IS-Se^5,5 es 1 . 5.11^s ea b cEJ = 0.8 x 108 Nm 2f1=0.00379Hzf2=0.0542Hzf=f2=0.054HzEJ = 0.8 x 109 Nm 2f1=0.012Hzf2=0.16Hzf=f2=0.16HzEJ = 0.2 x 10" Nm 2f1=0.06Hzf2=0.85Hzf=f2=0.85HzFigure 3-5: Change flexural rigidity—FFT analysis of the second resonance case73The first natural freqency flFigure 3-6: The trend of the first and the second natural frequencywhen the length of the riser is changedAnother important structural parameter is the riser length. In the case of a shortriser length, the amplitude of the displacement of platform grows at the beginning andthen declines to a small amount. With the length of riser increasing, the amplitude ofthe displacement becomes very large.The first resonance response (see Fig. 3-7 b.) show an interesting "beat" phenomenon.The shorter the length of riser, the stronger of the stiffness of the system. This makesthe displacement and the velocity of the system reduce to a significant degree and thehydrodynamic damping play a more important role. The nonlinear hydrodynamic dragterm contains varieties of frequency components. The slowly changing lower frequencycomponents combined with the first natural frequency component show the "beat" phe-nomenon.The FFT analysis in the second resonance case (see figure 3-10) indicates that when74the exciting frequency equals to the second natural frequency the first frequency responsecomponent may have a larger contribution than the second one.3.2.3 Tip massIt is observed from Figures 3-11, the tip mass has an influence in the first naturalfrequency while the second natural frequency is not sensitive to the increasing of the tipmass.Figure 3-12 shows that the heavier the mass, the larger the displacement of the tipmass in the first resonance case, but the influence of tip mass on the first resonance ofthe system is relatively small compared to that of the length and the flexural rigidity ofthe riser system. When the tip mass increases by 100 times, the maximum displacementincreases only 4 times as shown in Fig. 3-12. Similar to the natural frequency, the secondresonance of the system is not sensitive to the tip mass (see Fig. 3-14). This is becausebigger tip mass has a larger inertia which is hardly excited by the high frequency excitingforce.3.2.4 Ocean currentThere exist two kinds of interactions of wave and current, where the current may bein the direction of wave propagation or opposite to the direction of wave propagation.Figures 3-16 to 3-19 present the case in which the current is opposite to the directionof wave propagation. The behaviour of the current loading looks like a static force.When the amplitude of current velocity is much larger than the one of the wave, thedisplacement of the structure increases in the negative direction with small oscillationsaround a negative position. That is why there is only a small difference between themaximum displacement of the first and second resonance responses when the velocity ofcurrent is large.751.:0uP 0 1 11 1 1 11 1 1 1 11II I IIIIIIIIIIiillio:101 , 11111 11110■1 1 11 1il1111111111 111111 it'llitTql""i 101060E3 7.5060E131111111;.poolo- Te (s )0.0-0.5-1.0-1.511CA1 1151050 r-5-10-15 Times(s)^O0E3L^1.606bE4^1.5Qtlo^2.00i)0E4a b cL---10mfl =0.245Hzf2=12.88Hzf=f1=0.245HzL=101mfl =0.0075Hzf2=0.131Hzf=-f1=0.0075HzL=380mfl =0.0008Hzf2=0.0083Hzf=f1=0.0008HzFigure 3-7: Change length of the riser—the first resonance case76U r,0.0350 '0.03000.02500.02000.01500.01000.0050a b cL=10mfl =0.245Hzf2=12.88Hzf=f1=0.245HzL=101mf1=0.0075Hzf2=0.131Hzf=f1=0.0075HzL=380mf1=0.0008Hzf2=0.0083Hzf=f1=0.0008HzFigure 3-8: Change length of the riser—FFT analysis of the first resonance case77506 Jobb^15b0 20b0^25b0^30b0a b cL=1Ornf1=0.245Hzf2=12.88Hzf=f2=12.88HzL=101mf1=0.0075Hzf2=0.131Hzf=f2=0.131HzL=380mf1=0.0008Hzf2=0.0083Hzf=f1=0.0083HzFigure 3-9: Change length of the riser—the second resonance case78a b cL=10mfl =0.245Hzf2=12.88Hzf=f2=12.88HzL=101mfl =0 .0075Hzf2=0.131 Hzf=f2=0.131HzL=380mfl =0.0008Hzf2=0.0083Hzf=f1=0.0083HzFigure 3-10: Change length of the riser—FFT analysis of the second resonance case79(Hz)^The second natural frequency f2Fre0.0600.0500.0400.0300.0200.0100.000Figure 3-11: The trend of the first and the second natural frequencywhen the tip mass is changedThere are similar solutions when the current is in the direction of wave propagation asshown in Figures 3-20 to 3-23. What has happened in this case is that when the velocityof the current is large, the system oscillates around the zero equilibrium line after severalcycles of unsteady oscillations. Again the first and second resonant responses have thesame pattern for the large current velocity. The FFT analysis shows there exist manyfrequency components in the resonance cases.3.2.5 Ocean wavesOcean wave is the main exciting force for the marine riser system. The wave heightrepresents the magnitude of the loading. Comparisons of different amplitudes of thewave height and different wave propagation frequencies are presented in Figures 3-24 to3-28. An increase in amplitude of wave height results in a large displacement of tip mass80IIa^a I—2.5u p2.50.0bTime(s)4doolobo' .200. 430 6UP^a1.51.00.50.0—0.5—1.0—1.5 alla I 1 210 I^a. 1 a 13 Time(s)o08011 11 1111 0 no liffiilill11111111 . 11161111111 loTEime(s)up5.02.50.0—2.5—5.0a b cMp =8700kgfl =0.0088Hzf2=0.062Hzf=f1=0.0088HzMp=87000kgfl =0.00379Hzf2=0.0542Hzf=f1=0.00379HzMp=870000kgf1=0.00128Hzf2=0.053Hzf=f1-=0.00128HzFigure 3-12: Change tip mass—the first resonance case81a b cMp =8700kgfl =0.0088Hzf2=0.062Hzf=f1=0.0088HzMp=87000kgf1=0.00379Hzf2=0.0542Hzf=f1 =0.00379HzMp=870000kgfl =0.00128Hzf2=0.053Hzf=f1 =0.00128HzFigure 3-13: Change tip mass—FFT analysis of the first resonance case82a b cMp=8700kgfl =0.0088Hzf2=0.062Hzf=f2=0.062HzMp =87000kgfl =0.00379Hzf2=0.0542Hzf=f2=0.0542HzMp =870000kgfl =0.00128Hzf2=0.053Hzf=f2=0.053Hz0.000uP b-0.050-0.100-0.150-0.200-0.250-0.300 Time(s)250 500 7t0 1000 f2'50 15)30 17'50 2000u-0.000-0.050-0.100-0.150-0.200-0.250-0.300-0.350Figure 3-14: Change tip mass ^-the second resonance case83a b cMp =8700kgfl =0.0088Hzf2=0.062Hzf=f2=0.062HzMp=87000kgf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzMp =870000kgfl =0.00128Hzf2=0.053Hzf=f2=0.053HzFigure 3-15: Change tip mass—FFT analysis of the second resonance case84UP2.50.0560 lobo 15b0 2060 '1-25.60 dot b56bTime(s)000-2.50.0-2.5UP^A^ATime(s)oto fdoo 15b0 2D'60 '2560 Obb 3566 4obo-5.0a b clic =-0.126m/sfl =0.00379Hzf2=0.0542Hzf=f1=0.00379HzVc=-1.26m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzVc=-5m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzFigure 3-16: Change the amplitude of current(negative)—the first resonance case85a b cvc =-0.126m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzV=-1.26m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzK=-5m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzFigure 3-17: Change the amplitude of current(negative)^FFT analysisof the first resonance case860.000 up-0.050-0.100-0.150-0.200-0.250-0.300a b cVV =-0.126m/s 17,=-1.26m/s 14=-5m/sf1=0.00379Hz f1=0.00379Hz f1=0.00379Hzf2=0.0542Hz f2=0.0542Hz f2=0.0542Hzf=f2=0.0542Hz f=f2=0.0542Hz f=f2=0.0542HzFigure 3-18: Change the amplitude of current(negative)—the second resonance case871^0 . 5^1^II^5.^5.^I^SUP0.1500.1000.050aFreooq.(Hz)L4L.s is i . e is 5 . I is I.5•55^5 . 5 55 51.55^51a b cV=-0.126m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzVc =-1.26m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzV=-5m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzFigure 3-19: Change the amplitude of current(negative)—FFT analysisof the second resonance case88u P10 -CTime(s)250 500 7t0 1000 12'50 1500 1 700 20'00au p560 ^ obo 1500 2bbo 5^obb 350^obl-oi me (s)bUP1050-5-1050o^•2.50.0-2.52Ida) 1500 2600 500 30bb 5500Time(s)4000a b cVc =0.126m/s 17,=1.26m/s 17,=5m/sf1=0.00379Hz f1=0.00379Hz f1=0.00379Hzf2=0.0542Hz f2=0.0542Hz f2=0.0542Hzf=f1=0.00379Hz f=f1=0.00379Hz f=f1=0.00379HzFigure 3-20: Change the amplitude of current(positive)—the first resonance case89a b cVc =0.126m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzV=1.26m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzVc=5m/sf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzFigure 3-21: Change the amplitude of current(positive)—FFT analysisof the first resonance case90Cup10-0.3000.2500.2000.1500.1000.0500.000u p^b5.02.50.0-2.5-5.0Time(s)lobo^2obo^30b0^etobo^Sob()^sob00-10-20Time (S)2t0 560 7t0 1060 1200 1500 1700 2000-30a b c17c =0.126m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzVc=1.26m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzVc=5m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzFigure 3-22: Change the amplitude of current(positive)—the second resonance case91UP a0.1500.1000.050LAb...^ 1—55711-- i—s e 15 to •5 , 55 1 . 1 SI 55.51 5 IFreq.(Hz)00a b cVc =0.126m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzVc=1.26m/sf1=0.00379Hzf2-=-0.0542Hzf=f2=0.0542HzVc=5m/sf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzFigure 3-23: Change the amplitude of current(positive)—FFT analysisof the second resonance case92Or-1yi ( zo 172( zi)Y1 ( z 2 ) Y2 ( z2)(3.6)u(z i , 0)u(z 2 , 0)Y2 (z i )(3.7)Y2 ( z2)and a large nonlinearity. The pattern of the first resonant response is therefore highlyirregular. The results also show that the higher frequency response is not very important.The designers should pay more attention to the basic frequency response.3.3 Initial ConditionsIn order to get the ordinary differential equations of the system, a modal discretizationmethod was used. The direct numerical solutions by IMSL-DEGEAR is under the modalcoordinate. The relationships between natural coordinates and modal coordinates are:^u(z, t) =^yi (z). qz (t)2=1^t) =^yz (z) 4,(t)i=iIf we use two modes to approximate the system, the initial conditions necessary forthe integration procedure can be obtained by solving the following equations:u(z i , 0) = Yi(zi )qi(u) + Y2(zi)q2(0)u(z 2 , 0) = Yi(z2)q1(0) + Y2(z2)q2(0)iL(z i , 0) = 11(z041(0) + Y2(21 ) 42(0)il(z2 , 0) = li(z2)41(0) + Y2( z2)42(0) (3.5)93a b cH=0.53mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzH=5.3mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzH=53mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzFigure 3-24: Change the wave height—the first resonance case94a b cH=0.53mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzH=5.3mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzH=53mf1=0.00379Hzf2=0.0542Hzf=f1=0.00379HzFigure 3-25: Change the wave height—FFT analysis of the first resonance case95-0.000uP b-0.050-0.100-0.150-0.200-0.250-0.300-0.350Time(s)250 560 7t0 1000 1200 1560 170 2000UP C0.50.0-0.511-1.0 Time(s)2to 560 750 1060 1200 1560 1750 20000.000 uP-0.010-0.020-0.030-0.040-0.050a b cH=0.53mf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzH=5.3mf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzH=53mf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzFigure 3-26: Change the wave height—the second resonance case96a b cH=0.53mf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzH=5.3mf1=0.00379Hzf2=0.0542Hzf=f2=0.0542HzH=53mfl =0.00379Hzf2=0.0542Hzf=f2=0.0542Hz0.200 b0.1500.1000.050 Freq.(Hz)00Lab.—I Ciiris so 1 . 1 55 es•os 5 . 5 IS 54.55 IFigure 3-27: Change wave height—FFT analysis of the second resonance case97a b c df=0.000379Hznonresonancef=f1=0.00379Hzfirst resonancef=f2=0.0542Hzsecond resonancef=0.54HznonresonanceU .to a0.50.0ilk-0.5-1.0- 1.5 Time(s)I,^III^se^III^IS^III^IS^•000UP , • b A 02.50.0-2.54 4 1-.,^4 Time(s)5 0^1000-^1 bnOU^20Su0^25,oU^,4011.)^,ibbri^21(1000.000u C-0.050-0.100 --0.150 --0.200. Jiff-0.250 --0.300 :^ Time(s)10^200^300^400^SOU^Ei0U^/0U^800-0.000uP d-0.050-0.100-0.150-0.200-0.250-0.300 Time(s)Is^is^Is^• is^is^.ss^is^:soFigure 3-28: Change the wave frequency98Natural coordinates z 1 and z 2 must satisfy:171 (z i ) Y2(zi.)0 (3.8)Y1 (z2 ) Y2(z2)where u(z i , 0) and u(z2 , 0) are corresponding initial displacements. The influences of theinitial conditions to the response of the system are investigated. Figure 3-29 shows thecase where the exciting frequency is far below the first natural frequency while Figures3-30, Figure 3-31 and Figure 3-32 show the cases where the exciting frequencies are eitherequal to one of the natural frequencies or much larger than the second natural frequency,respectively.Each figure includes four different initial conditions:a. at rest, i.e. both initial displacement and initial velocity are zero;b. initial displacements are none zero but initial velocities are zero;c. initial displacements are zero but one of the initial velocity is riot zero;d. initial displacements and one of the initial velocity are not zero.All four cases considered show that the initial conditions have an influence in thetransient behaviour. After several circles of unsteady oscillations the system enter into itssteady state. Between the initial displacement and initial velocity, the latter has a moresignificant effect on the dynamic response, which makes the amplitude of displacementof the mass go to the steady state rapidly in the first and second resonance cases (seeFig. 3-30 and 3-31 c and d).99• f1 = 0.00379Hz • f2 = 0.054Hz of = 0.000379Hza b c du(z 1 ,0) = 0u(L,O) = 0U(z i , 0) = 0ih(L,O) = 0u(z1 ,0) = 0.17mu(L,O) = 0.4mU(zi, 0) = 0ii,(L,O) = 0u(z 1 ,0) = 0u(L,O) = 0U(zi , 0) = 0iL(L, 0) = lm/su(z 1 ,0) = 0.17mu(L, 0) = 0.4miL(Zi , 0) = 0U(L, 0) = lm/sU .1.0 a0.50.0ilk-0.5-1.0-1.5Time(s)III^III III III Ill III II 'GOOu1.0 b0.50.0-0.5ilk-1.0 Time(s)-1.5 III^III SI GOO IS III SI • 400U C3.02.01.00.0-1.0 Time(s)is^III Is III se III se • 00u d3.02.01.00.0-1.0 Time(s)is^see is iii is see le • 00Figure 3-29: Change the initial conditions—the wave frequencyis far below the first natural frequency100U P2.50.0-2.52.50.0-2.55bo ^ 1Cou '2b130 .3(.VUTime(s)".3511^4000U Pa2.50.050o^1 ubu^151)U^2UUU^.3U^.351\fl \WO-2.5\'251SU Time(s)J J z^T(b00 u00 lbOu ubu'ebOu Juttu 3brit) u00 ime s)djTime(s)25 U .300)^duo°2.50.0-2.550o 1doo ^ 1 bb0 \J20 U• = 0.00379Hz • f2 = 0 054Hz of = fi = 0.00379Hza b c du(zi , 0) = 0u(L,O) = 0U(z i , 0) = 0ii(L, 0) = 0u(zi , 0) = 0.17mu(L, 0) = 0.4mii.(zi, 0) = 0U(L, 0) = 0u(z i , 0) = 0u(L, 0) = 0U(zi, 0) = 0U(L, 0) = lm/su(z i , 0) = 0.17mu(L, 0) = 0.4mU(zi, 0) = 07..( L,0) = lin I sFigure 3-30: Change the initial conditions^ the first resonance case101•L 0.00379Hz of2 = 0.054Hz of = f2 = 0.054Hza b c du(z i , 0) = 0u(L , 0) = 0U(z i , 0) = 0U(L, 0) = 0u(z i , 0) = 0.17mu(L, 0) = 0.4miL(zi, 0) = 0ii(L, 0) = 0u(z i , 0) = 0u(L, 0) = 0U(zi, 0) = 0U(L, 0) = lrn/ su(z i , 0) = 0.17mu(L , 0) = 0.4mii(zi, 0) = 014,0) = 1m/ so.000u • a-0.050-0.100-0.150-0.200 --0.250-0.300 Time(s)es III IS III Is III II •soou. b0.25 '0.00-0.25-0.50 Time(s)is 545 Is III II los It , S00U. C2.0 ,1.00.0-1.0Time(s)II III It III II 111 III •WOU .d2.01.00.0-1.0Time(s)SI III II III II III II •100Figure 3-31: Change the initial conditions—the second resonance case102= 0.00379Hz sf2 = 0.054Hz of = 0.54Hza b c du(z 1 ,0) = 0u(L, 0) = 0iL(z i3 O) = 0i.L(L, 0) = 0u(z 1 ,0) = 0.17mu(L,O) = 0.4mii(zi3O) = 0ii(L,O) = 0u(z 1 ,0) = 0u(L,O) = 0ii(z1,0) = 0ii(L , 0) = lm I su(z 1 ,0) = 0.17mu(L, 0) = 0.4mii(zi, 0) = 0U(L, 0) = lm I s-o.000u-0.050-0.100-0.150-0.200-0.250• a-0.300 Time(s)is is is ,is is .III •s :110u,0.25 ih, b0.00-0.25-0.50 :s0Time(s)so is is • is se .41 esU. C2.01.00.0-1.0Time(s)III II is • Ill II .ii is :IIU• d2.01.00.0-1.0ILLTime(s)is is II • is so .1i se :10Figure 3-32: Change the initial conditions—the wave frequencyis larger than the second natural frequency103Chapter 4CLOSED-FORM SOLUTION USING THE BUTENIN METHODThe governing nonlinear and coupled equations of motion do not admit any exactclosed-form solution. Several approximation procedures can be used to investigate thelarge amplitude motion of the system. The variation of parameters method, which wassuggested by Butenin[30], is intended for differential equations with small nonlinearities.4.1 Simplified Equations of MotionThe two dimensional differential equations in modal coordinates are:+ Knqi + c1242oL^1+ EJ[{41 1711(z) + 42 1721(z)} 2 + 1foLEJ [4X(z) + 421721(4[0'in (z)([40711(z) + 42Y2(4 21] [07 ( z ) + q217(z)11711 1(z)dz_A_ q2 Y2" ( z )1 2 y,i( z )dz1)2=- Yi(z)[CmAri)„,.CDADIvwx+vcs — 4iYi(z) — 42Y2(z)1{v..±ve. — 4iYi(z) — 421/2(z)}idzYi (L)Czpl g Hk.sin(wt)fo L i (z)[(07.11 (z) + q 2 Y;(z)](C z pi' g H k.tanh(kL)cos(wi))dzM2242 + K22q2 C2141rL 1+ ifqiYii(z) + q2Y2(z)} 2 + 1^1][07(z) + q2n(z)1Y"2(z)dzY2fL ,[q13711 (z) + q2Y21(z)][01 "(z) + q2Y211(z)l2,2 z)dz, (o rid ([41}'(z)q2y2,(z)12 + 1 )2104Lfoj.:y2(z)[CmAri)..+ CDADIvavr+vc. — 401(z) — 42Y2(2)1{vwx+vc. — 41Y2(2) — 42}2(z)Ndz+ Y2 (L)CxpV gHk.sin(wt)1112(z ) [qiY1(z)+ q2 YAz)1(Cz pVgHk.tanh(kL)cos(cot ) )dz^(4.2)Some of the nonlinear terms in the equations 4.1 and 4.2 may be approximated byapplying the binomial theorem. Retaining nonlinearities up to and including the thirdorder in coordinate q1 , we get:1^1][qiInz) + q21 72"(z)0 7"1(z)dzfo EJ[{q i YAz) + 92}721 (z)} 2 + 1L EJ[q i YAz) + q2 Y2 (z)llq 1 1","(z) + q2}T(z)]2Yi1(z)dzo^ ([40711(2)^92}2(z)]2 + 1) 23^2 ^3.— ai2402 — a1391q 2 — a i4 q 2 ,fL^1Jo EJ[ {07(z)+q2n(z)}2 +1foL EJ [giv(z)+q2Y2'(z)i[qiYi"(z)+q2172"(z)121.72,(z)dzugly? (z)^92}2(z)]24_ 1)23^2 ^3.— a 22q02 — anqiq 2 — a24q 2 ,where,La n = I 2E J(17) 2 (1;11 ) 2 dzLa12 =^3EJ1/11 1'7 1(I2}TLa13 =^EJ[2171 /1"1 721 32"^(1-7;17" +1711 1/2") 2 1dzLa 14^EJY21211(12T2" l2T2")dzl][q i iin(z ) + 42 172"(z) i Y 11 2(z )dz1053a21 — a12a22 — a13a23 = 3 a14rLa 24 =^2EJ(1:21 ) 2 (}P 2 dzThe differential equations thus can be simplified as,M114 .1 + c1242^K11(/1an 3 ^2^3+ a124 1 42 + a i 3q02 + a i4 q 2• C DA Dly^— 401(2) — 4212(z)1{vwx+vc. — 402(z) - 4 2 12(2)Wiz+ Y2(L)C x p1.7 g H k.sin(wt)Y i i(z)1/211711 (z)-1- q2 1 (z)}(Cz pl 7gHk.tankkL)cos(wMdz11,1224 .2 + C2141 + K11423 +^2 +^ 2^,.3= an qi^a224192 a23viv2^a241/2iLjo Y2(Z)[CMA/ 1)wx▪ CDADIV„x ±Vcx — 41 Yi (z) — 4212(z)l{v..±vc. — 4072(z) — 42} 2 (z)fidz+ Y2(L)C xpV g H k.sin(wt)L2(z)EqiY1(2) + q2Y2I (z)1(C z p1 7 g H k.tanh(kL)cos(bA))dzIntroducing the time scaleWt =T106We have the dimensionless form of differential equations:c1d2qi^ A idq2^ + nl 2 ql7- 2^r1^,^3^2^2^3— L^l11111w 2W 2allqi -F- adiv2:42 + ai39142 + ai442+ lo Yi (z)[C,A iw^drd+ CDADIvt" — v`x — 1'1(4° ddrqi Y2(z)w drq2(v. — v. — Yi(z)w -dq idT — Y2(z)a) dTdq2 )1dz+ Yi (L)Cx pV g H k.sin(r)r i,Jo 1"2 (z)[q 1 Y1(z) +q 2 Y 1 (z)]C z pi ,/ g H k.tanh(k L)cos(r)dz}d2 q2^dqi ,^2^ + A2 ^ n2 42d7 21141 22W 2 r^3^2^2^3^ la2lqi^a22q1q2^a23q1q2^a24q2foL dv,,Y2 (z)[C,,A iw^dTq2dC DADIvwx — vcs — li(z)u) dr^Y2( 4') dr I.(vwx — vex — Yi (z)Lodq,^Y2 (z)wdq2 )]dzdr drY2 (L)C x pV g H k.sin(r)fLJo Y i2(z)(41 }711(z) 92Y2 (z)]C .07 g H k.tanh(k L)cos(T)} dz(4.3)(4.4)where,A2—C12= ^1 (Ai—C1211122Ww l \ 2n 2 = ( _ )(.4.)107W 22LetM11K22AI2 2f0 nLy fz) • CD •AD.Vw2s(z)dZA111w 2foL YI(Z).CD. AD.( srlicoshkz \ 2 dzT sinhkz 1= w2when current velocity is much larger than wave one. u is a dimensionless coefficient (smallparameter) which measures the closeness of the approximation of the given system to alinear one. The differential equations can be expressed as:— A142 + nig]. — 14(41, 41, 42,42,T) + DoinT^(4.5)+ A241 + 71,42 =^41, 42, 42, T) + D 2 sinT^(4.6)where f(41, 41, q2, 42, 7 ) and g(qi, 41, q2, 42,r) are nonlinear functions of 41,41, 42, 42 andT.4.2 Butenin's Approach-A Variation of Parameters Method4.2.1 The eigenvalues and the solutions corresponding to linear systemNeglecting the nonlinear parts and the hydrodynamic force of equations 4.3 and 4.4,we have:^— Al42 + T441 = 0^ (4.7)42 + A241 + 74q2^0 (4.8)Milw 2when wave velocity is much larger than current one, or letfoL (Z).CD. AD.V c2x (Z)dZ108Solution of the problem 4.8 can be written as:ql = Qle (3T +13)q2 = Ci2e ("1-13)where Q i and Q2 are complex eigenvectors, while s is a complexform of eigenvalue problem becomes:(4.9)eigenvalue. The matrix(3 2 +^e (sr+o) —^0(s2 + n22)^Q2( S 2 + n?)A2S— A i s2 + n3)0 = 0let,or,^4^2^■ n^ in e^unS /2 + n 2 + A1A2)3 2 n 2 = It follows fromA l A 2^0and,2^ 1 —( 74 + 74^ 2+ A 2 )±0n i2 -i- n 2 A 1 A^22)2 — 4nn 22S 12 = < 02The eigenvalues are:s l = ik l =[(n i2 + n 22 + A 1 A 2 ) —^/(n2\^1 + ,,,2^,,2,„21'`2^/11 /,2 \2^A^''212[(n? + Al A 2 ) + ^(724 + A1A2)2 — 4n7 n2]3 2 — ik2 2109The complex modes are determined by:[(—q + nD --A i ik;A 2 ik; (—k3? 1- nD _ Q200(j = 1, 2) Normalizing w.r.t. Q i = 1, the complex modal matrix is:[ A =1^1Q (21) Q (22)where,i(k2_ ,2/1(1) = •-^1 i A 2kik22-1^'''2,;(1, 2^02k2QV) = 2a2 ''2 ^k2^1c3^—If = 0, the differential equations will be,— A1q2 + n;. q i = D l sinTq2 — )2qi + n 2 q2 = D 2 sinTIn order to get the solutions more conveniently, let the harmonic function be representedin complex form as F3 (t)=D zr (j=1,2) so that the equations of motion become:^— A142 + T4q1 = D l eIT^ (4.10)^q2 — A2ql + n 2 q 2 = D2e'T^(4.11)Since the actual excitation is given only by the imaginary part of D3 Or, the response willalso be given only by the the imaginary part of q,(t) where (11 (0 (j=1,2) is a complex110quantity satisfying the differential equations 4.10 and 4.11. By coverting to complexmodal coordinate,{q} = [MOwe obtain the steady state solution w.r.t. 6 coordinate as:{6} = H 1 De'where,H-1[A2 + a 2 (74 — 1)]^ij(7q^1) + A1c1 21— [A2 +^— 1)]^— 1) + A l a i ] AA = (a 2 — cx i )[(74 — 1)(r4 — 1) — A1A2]D2}Thus, the steady state solutions w.r.t. 6 coordinate are:D1[A2 + ct2(n — 1)]eir+ D2A [(n 2 — 1) + A 1ck2] et(r+2)2 =^D1^ [A2 a i (r4 — O]ezr2D 2A^— 1) + A i adet ( T+'21)The complete solutions w.r.t. 6 coordinate are:= aeq k1 r+r31) + Q [A 2 + cf 2 (n. — 1)]eirD2 A [( 712 - 1 ) I Aia2]e(T+D1 116 = be (k2 T -L82) —[A 2 +^— 1 )]e irD2— [(n 12 — 1) + A i o l ]el (r+Dwhich can be written w.r.t. q coordinate as:q1 = 6 4-2beo2T+02)D1 2^D2^itr+,(n 2 — 1)e" ^Ale ■ 2/LAI 01q2^iC116 2:a26c aei (k i r +131+^a2bet(k2T+132+1)A2 arr.+ 21)+ D2D 2-r A (n^1)e' ^ e^2"i 1(4.12 )(4.13)The expressions of the actual response are given by the imaginary parts of the equationsof 4.12 and 4.13:where,= asin(k i T 01) + bsin(k 2 T 02)D1^D2 + ^ (n22^1)sinT ^AicosT1 A1q2 = a i acos(k i + 01) a 2 bcos(k 2 7- + /32)D2 2^D1 A2— 1)sinT ^COSTAl= (r4 — 1)(74 — 1) — A1A2(4.14)(4.15)4.2.2 Solutions corresponding to nonlinear systemFirst, we study the case when the natural linear frequencies are different from thefrequency of the external force.112Let /.40 and we want to find solutions of the system 4.5 and 4.6 in the form as 4.14and 4.15, but where a, b, 0 1 and 02 are slowly varying functions of time.By performing the averaging technique, we get the approximate equations for a, b,0 1 and /32 defined as:^da^ A2^ (A1A3^ A1)T^2(k? — 14)^db A2T^2(k? – k3)( AiB3 —B1)^d01^P^A^,T^2a(k? — k3)(A1 A4 + 2—A2)c102 A2^dT^2b(k? —Ic3) ( B4 ^+ (7,-2B2)where,pw 1.o27o27,-47r 3 J ° J^= ^f *cos0 Orielr^ (4.16)1 f2w f2w f2w^A2 = ^*.szn0Oncly^ (4.17)47r Jo Jo Jo471. 3^ /027- 1027r 1 02-ff g * s z n0071 drA3 =1 j(27, j( 27,1 2 ,„gA 44 4^7 3 0 0 0 cosdedndT1^[ 271- f27- f 271-f *cosnckdridTB 1 = ^473 / 0 Jo Ju(4.18)(4.19)(4.20)1f27, pi,B2 = ^47r 3 J 0 Jo Jo f *szraidOwly (4.21)1131^r2x f2ir /27rB3 = ^ * S Zni47r 3 ,10 J o Jo gA thidT1^/27r^27r27r47r 3 J ofoB4 — g* cosriddr dr(4.22)(4.23)Here1c 1 7-^/31k 2 T + 02In order to integrate the nonlinear force, we need to know whenT) = vwx(r) — vc. — Yi(z)4g, 7- ) — 12(z)42(, TpoIn fact, the moment at which x(,7),^= 0 may be arbitrary during the period of 27r.Because Butenin's method itself is an averaging method and the chance for x(,7i,T)>0is the same with one for^77 , T)<0 during the period of 27r, we can simply suppose:X( > T) o^when(0^< 7, 0 < g < 7, 0 < -r < 7)X(4. , n, T) < o^when(R-^. <277-,7r <71<27,7 < T<27r)The integrations of x(,7), T) term in the equations 4.16 to 4.23 thus have the followingform:1^ frfriw{fLYI(z).CD.ADX2(4.,7/, r)dzIcoserkdric/T47 3 0001^[27r [ 27,- .)( 2.7,4^ 73 J R. J R.^{^(z)•CD•ADX2(, 71, T)d.z}cosddridT1^fr r f ir {I L (z).cD.ADx 2 (,n, 7- )dz}sin000-7-47 3 00001 pirf2R- 1271-47r3f xYi(Z).CD.ADX 2(, 7- )dz}sin0OrrIT7,^7r^0AiA;114A 93' =^1 j r ri" rff Y2 (z).0 D. A Dx 2 (6,71, T)dz} sin6d6dr dr47r 3 o 0001 pir[2w f47r3 w^foL Y2(Z).CD. ADXV )71 7)d21.52n ed6c1T dTA: = 147-3 oofirl ow IL Y2 (z).cD .ADx 2(6, n , 7- )dzlcos6c100T1 [2w [2w f 2w f L471-3 J Y2(Z)•CD • ADX 2 , 71 T)d.Z1CO306C171dTw^w^oB: =^1 17r r r { ILY(Z).CD.ADX2(6)77) r)d.Z}COSTId6d/idT47r 3 J oJ fo o Jo1^ w^2w ILY4rr3 f x fi(Z).0 D. ADX 2 (e 71 T)c12} CO371d6dOTw f w 0B2* = 471r 3 1:1:107r{foLYI(Z).CD.ADX2(6, T)d.z} siny dOri dr1 f i2w 2wf 2w IL411-3 /^7} 1(Z)•CD-ADX 2 (6) 71) r)d.Z}Sirind6d7idTww w 0=--^1471-3J^fo or { fL Y2(z).CD•ADX 2(, 7- )dz} simc107-id6ro 1 f i2w rwi2w iL47r 3 i^{^2(Z)•CD.ADX2(6, 71, T)d21.5i7171d6d71dTww w= 4710 iwo fwofow{ L 172(Z).0 D • ADX 2 (6 71, T)C1Z1CO371d6d71dT 01 j-2wf2wf^y2w47_ 3^w 2(Z)•CD•ADX2(6, rhi)d,Z}COSTiddildTBy calculation, we find,da „2^D 1,2 P=^7- 12 ,,^nab^- 16dr^Pl4 Pi5b Pdb D 2 i p j„21D^D24 a + P25 b P2621a^22 u^I 23 a U^I-dr431a^ = B11a 3 + B12ab 2 + B13a 2 /314 b2 + B l ob+ B16 a + /3 17 1) + B18dr115A*(1^c i e2A.P11r)a = (1 — cle2A.P11T)^2P11(4.25)432^ B21b3 + B22ba 2 + B23a 2 + B24b2 B25 ab B26a B27b + B28dTThe first two of these equations are coupled nonlinear equations which have no exactclosed-form solutions. However, from the numerical analysis we found the second pa-rameter b played less important role than the first one a. We can therefore assumethatdb 0dTThe solutions under this consideration have the following form:b= bo^(4.24)2B 11 A*Pii(1 — c i e 2 A*Picr)—^+ lnicie2A.P11T —P121(Bi4b2 B17b MOP* (2* 1n1Q* — P* ce 2A.13" 1- 1p. Q.1311 (41-11 ,4 *2 -^A* B + B 2 )4P121^ Bi3A*/314 b2 Bi7b + B18 Bi2b2^(2*^Bi5b 1316IT+ Oio^ (4.26)PAb+(B25 b + B26) (B21b2 B24b -+ B27 B28/b)(B23 + B22b)(4/3121A •2 —^A* B + B2))7-4/3121 b+ 020^ (4.27)2)1* (B23 + B22 b)2 Iiib(1 _ce2A.PuT)B(B23+?22b)^ 2A.P11,_ ce^116where,A* =^-P-1119 + P14)^(Pi2b2 + Pisb + Pm) ^2P11 P11B = P13b + P14Pi, and B13 are constants determined by structure and nonlinear force coefficients, andc, b0 , Pio and 020 are constants determined by initial conditions.The final solutions given by Butenin's method are:^u(z,t)^Yi (z)q i (t)^172 (z)q2 (t)asin(k i T^)^bs in( k 2 T + 132)+ A1 2(n^D2^ 2 1 )sinT^A1 Al COSTq2^ckiacos(kiT + 131)^Ce2 bCOS( k2 T^Q2)^D2 2^DI A 2^(n i^1)sinr ^ COST^Al Aiwhere a, b, 13 1 and 02 are determined by equations 4.24 to 4.27.Let us now look at the second case when one of the basic frequencies of the system isequal to the frequency of the external force. We suppose that the amplitude of the linearpart and the nonlinear part of the external force have the same order in a neighborhoodof zero as p.In the resonance case, say, k 1 =1, according to Butenin's method, the solutions shouldhave the following form:= asin(7- + 0 1 ) bsin(k2 T + 132)Q2 = oc i acos(7- + 0 1 ) a 2 bcos(k 2 T + Q2)117Compared with the form of the solution of the nonresonance case, the solutions ofthe resonance case have no steady state part. In the same way as before we find thefollowing approximate equations for a, b, 0 1 and 02:da= ^ (Al243 —A1)dT^2(k? – 14)dT^2( k? kdrbII3Q2^ `1"^2B1)a 2d01^ 1A4 + ^A2)dT^2a(k — k)dfr32dr = 2b(k? — kn(A1B4 a2 B2)where1 f J27r^2-nAl = ^f* cosddridr47 3 J o o1/2o 7r /o2712o 7rA2 = 47r3^f* sinckdridrJ JA3 –1 1 27r/ 27r/ 27rg*.sin0007-47r 3 o^o^o1^/27r 274-127rA4= g*cosdOrid-r47r 3 J o Jo o1^ 2r ,2r 27r471- 3 / 0 J 0 J 0 rcos7ickchidr1/2w i 27,1 2.,47r 3 J 0 Jo^f*sinrAdndrB1 =B2 =1181^[ 27r 27r /27rB3 = ^47r 3 i 0 j o^0 .9 * sirrilddridr1 /27, 2ir f airB4 = 4^7 3 0 0 0 g comickchichDifferent from the nonresonance case, the functions of f* and g* not only contain thenonlinear terms but also include the linear force.The approximate equations for a, b, 0 1 and 02 are of the same form as equations 4.24to 4.27 but have different coefficients.4.3 Typical Response Results: Comparison Between Analytical and Numer-ical DataThe validity of the analytical solution was assessed over a range of different excitingfrequencies. The initial conditions for analytical solutions are arbitrarily appointed. Inthe a portion of the Figure 4-1, the initial conditions are at rest'. In the b, c and dportion of the Figure 4-1, the initial displacements and the initial velocities are not zero.It is significant to note that the analytical procedure is able to predict the maximumamplitude as well as the phase of the oscillations with an acceptable degree of accuracy.The error between the analytical solutions and the numerical solutions is within 15%,which may be due to the simplification of the equation of the analytical procedure, ordue to the approximation of the analytical solution, or due to the error accumulation ofboth methods.However, the accuracy drops off for the large nonlinear force. The Fig. 4-2 a indicatesthat when the wave height increases to a large number, the hydrodynamic force termsbecomes highly nonlinear. The dimensionless coefficient p, which measures the closenessof the approximation of the system to a linear one is not a small parameter (in the caseof Fig. 4-2, y=0.89).119a b c du(z i , 0) = 0u(L, 0) = 0ii(z i , 0) = 0ti(L, 0) = 0u(zi , 0) = 0.17mu(L,O) = 0.4mit(zi, 0) = 0ii,(L, 0) = lm/su(z i , 0) = 0.17mu(L,O) = 0.4mii(zi, 0) = 0it(L, 0) = lm/su(z i , 0) = 0.17mu(L,O) = 0.4m14.zi , 0) = 011(L, 0) = lm/su•to0.50.0-0.5-1.0-1.5a(f=0.000379Hz)'s^ Time(s)411114Is^III^II^III^II^III^IS^, esoup2.50.0-2.5b(f=f1=0.00379Hz)kilLA ^Time(s)ss^I^gooU.2.01.00.0-1.0c(f=f2=0.0542Hz),........ ... ^- -Time(s)se^Is^55^•Is^Si^.11^is^:sou2.01.00.0-1.0Figure 4 -d(f=0.54Hz)__:-.-•:"1--"7---------"^ 011:^• ?'5 • ar - 'V^"Tale ''''e.^i: 1.1 ion: III^•^. : t^:'•.I Tsial i120It should be noticed that even in the case of a small parameter it, the closed-form so-lution sometimes cannot predicate the profile of the amplitude of the oscillations exactly.Fig. 4-2 b. and c. represents the case in which the analytical solution can only predicatethe maximum amplitude of the oscillation. The reason is possibly because of setting theparameter b as a constant which reduced the accuracy of the analytical solutions.121up^1.0^101111111111111^0.5 001 I^ill0.0-0.5-1.0 11111111 bo 11111111 1111111300 40 111111 1501)0 Time (s)6000iiiII ip5o0,I - c 101 titU P1.00.50.0-0.5-1.0 Timeo(s)11 1 1 1 1 1 11 1 1 1i iil6i SUdig Ia b cH=53mDashed line-analytical solutionSolid line - numerical solutionf=f1=0.00379HzH=5mNumerical solutionf1=0.012Hzf=f1=0.012HzH=5mAnalytical solutionf1=0.012Hzf=f1=0.012HzFigure 4-2: Comparison of numerical solutions and analytical solutions122Chapter 5CONCLUDING REMARKS5.1 Summary of ResultsThe equations of motion for a marine riser undergoing large three dimensional de-flections and rotations are derived using modal discretization approach. Specifically, thethree dimensional shear effects have been addressed.A simple model based on two dimensional motion is used to carry out the parametricstudies of marine riser system due to excitations of the ocean wave and current force.The parametric study represents an important step forward in understanding the dy-namic behaviour of a marine riser. It is particularly relevant to the next generation ofmarine riser systems which will be located in ever increasing deep depths and harsherenvironmental conditions. Based on previous parametric analysis the following generalconclusions can be made.(i) The most important factor in the list of parameters studied which affect the marineriser system's deflections is the length of riser. In the case of short length of riser theamplitude of the displacement of platform grows at the beginning and then declines toa small amount. With the length of riser increasing, the amplitude of the displacementbecomes very large. Measures, such as enhancing the flexural rigidity or being equippedwith nutation damper, should be taken to reduce the deflection of the riser system.(ii) As expected, with the flexural rigidity increasing, the natural frequency of theriser system increases and the displacement of the platform decreases greatly.123(iii) The mass of the platform has an influence on the first resonance of the system.The heavier the mass, the larger the displacement of the platform. The second naturalfrequency and the second resonance are not sensitive to the mass at all. This is becausethe heavier mass has a bigger inertia which is hardly excited by higher frequency forces.(iv) Perhaps one of the most interesting phenomenon is the 'beat' happened duringfirst resonance. When the stiffness of the system increase, the displacements of thesystem and the derivatives of the deflections are reduced to a relatively small amount. Inthe drag force, the wave and current velocity terms which contain varieties of frequencycomponents play more important roles. The slowly changing lower frequency componentscombined with the first natural frequency component show the 'beat' phenomenon.(v) The behaviour of the current loading looks like a static force. When the amplitudeof velocity of current is much larger than that of the ocean waves and its direction isopposite the direction of wave propagation, the displacement of the structure tends toa negative position with small oscillations around the equilibrium center. It should benoticed that when the velocity of the current is in the direction of wave propagation andits amplitude is large enough, instead of going to a positive equilibrium position, thesystem oscillates around the zero equilibrium position after several cycles of unsteadyoscillations.(vi) Ocean wave force is the main exciting resource of the marine riser system. Thelarge amplitude of wave height produces a large displacement of the marine riser system.The comparisons of the response under the different exciting force frequencies show thatthe first resonance always gives rise to largest displacement.(vii) The initial velocity has a large influence on the dynamic response of riser system,which makes the amplitude of displacement of mass reach to steady state rapidly.1245.2 Recommendation for Future StudyThis thesis represents only a beginning in the exploration of complicated three dimen-sional motion of a marine riser system. Based on the studies of a simplified marine risermodel, the following suggestions are made for future work:(i) Because of the interactions of different direction motions of marine riser, thethree dimensional dynamic behaviour cannot be represented by a two dimensional one.The parametric studies based on three dimensional dynamic analysis and the studies ofinteractions between three translational and rotational motion are needed.(ii) As mentioned at the beginning of the thesis, all the previous work done is basedon Eulerian beam theory. No one up till now has explored the three dimensional sheareffects of the marine riser. The reason is mainly because in the linear case the shear effecthas only small contribution to the low frequency response. The consideration of threedimensional nonlinear shear effects makes the equations of motion extremely complicated.The contribution of shear effect to the nonlinear dynamic behaviour of marine riser, basedon interaction between bending effect and shear effect, should be studied, especially inthe high frequency response case.(iii) In the case of a long marine riser and large amplitude of wave height, the risersystem gives rise to large displacements. Man-made damping such as nutation dampingmay be introduced into the riser system in order to reduce the amount of displacementof the marine riser system. However, the studies based on nonlinear elastic beam theoryof such a system will be very complicated and thus need to be carried out intensively.(iv) The diameter of the marine riser in this thesis is kept as constant along the lengthof the riser. When the geometry, say diameter, is changed along the length or severalmasses and nutation damping structures are placed along the riser, the modal functionsand the nonlinear differential equations of the system in terms of modal coordinates may125be easily handled by means of transfer matrices.126BIBLIOGRAPHY[1] Chateau, G.M. "Oil and Gas Production Facilities for Very Deep Water,"Proceedingof the Third International Behaviour of Offshore Structures Conference, BOSS'82,Boston, Massachusetts, 1982, Vol. I, pp. 50-70.[2] de Oliveira, J.G., Goto, Y., and Okamoto, T., "Theoretical and MethodologicalApproaches to Flexible Pipe Design and Application,"Proceedings of the Seven-teenth Offshore Technology Conference, Houston, Texas, 1985, OTC 5021, Vol. 3,pp. 517-526.[3] Skomedal, N.G., "Dynamic Analysis of Flexible Risers: A Comparison with ARandom Vortex Method and the Morison's Equation" Proceedings of the Ninth In-ternational Conference on Offshore Mechanics and Arctic Engineering, New York,N.Y., 1990.[4] Panicker, N.N, and Yancey, I.R., "Deep water Production Riser,"Proceedings of theFifteenth Offshore Technology Conference, Houston, Texas, 1985, OTC 5021, Vol.3, pp. 517-526.[5] Chakrabarti, S.K., and Frampton, R.E., "Review of Riser Analysis Tech-niques," Applied Ocean Research, 1982, Vol. 4, No. 2, pp. 73-90.[6] Chakrabarti, S.K., "Hydrodynamics of offshore structures, "Computational Me-chanics Publications, Southampton Boston 1987.[7] Sarkaya T., and Isaacson, M., "Mechanics of Wave Forces on Offshore Structures,"Van Nostrand Reinhold, New York, 1981.[8] Wang, E., "Analysis of Two 13, 200 ft. Riser Systems Using a Three DimensionalProgram," Proceedings of the Fifteenth Offshore Technology Conference, Houston,Texas, 1983, OTC 4563, Vol. 3, pp. 435-444.[9] Chucheepsakul, S., "Large Displacement Analysis of a Marine Riser," Ph.D. Thesis,1983, University of Texas at Arlington, Arlington, U.S.A.[10] McNamara, J.F., O'Brien, P.J., and Gilroy, S.C., "Nonlinear Analysis of FlexibleRisers Using Hybrid Finite Elements," ASME Journal of Offshore Mechanics andArctic Engineering, 1988, Vol. 110, No. 3, pp. 197-204.127[11] Vogel, H., and Natvig, B.J., "Dynamics of Flexible Riser Systems," Proceedingsof the Fifth International Offshore Mechanics and Arctic Engineering Symposium,Tokyo, Japan, 1986, Vol. III, pp. 363-370.[12] Isaacson, M., "Wave and Current Forces on Fixed Offshore Structures," CanadianJournal of Civil Engineering, 1988, Vol. 15, No. 6, pp. 937-947.[13] Ismail, N.M., "Effects of Wave-Current Interaction on the Design of Marine Struc-tures," Proceedings of the Fifteenth Offshore Technology Conference, Houston,Texas, 1983, OTC 4615, Vol. 3, pp. 307-316.[14] Brouwers, J.J.H., "Analytical Methods for Predicting the Response of Marine Ris-ers,"Proceedings of the Royal Netherland Academy of Arts and Sciences, 1982, Se-ries B, Vol. 85, pp. 381-400.[15] Love, A.E.H., "A Treatise on the Mathematical Theory of Elasticity," Fourth Edi-tion, Dover Publication, New York, 1944.[16] Novozhilov, V.V., "Foundations of the Nonlinear Theory of Elasticity," GraylockPress, Rochester, New York.[17] Finlayson, B.A., "The Method of Weighted Residuals and Variational Principles,"Academic Press, New York, 1972.[18] O,Brien, P.J., McNamara, J.F., and Dunne, F.P.E., "Three-Dimensional NonlinearMotions of Risers and Offshore Loading Towers," Proceedings of the Sixth Inter-national Offshore Mechanics and Arctic Engineering Symposium, Houston, Texas,1987, Vol. I, pp. 171-176.[19] Nordgren, R. P., "On Computation of the Motion of Elastic Rods", Jam, Sept.1974, pp. 777-780.[20] Garrett, D.L., "Dynamic Analysis of slender Rods", JERT, Vol. 104, Dec. 1982,pp. 302-306.[21] Verner, E.A. et al., "Predicting Motions of Long Towed Pipe Strings", Paper No.OTC 4666, Offshore Technology Conference, Houston, Teaxs, 1984, pp. 159-169.[22] Chung, J.S., Whitney, A.K. and Loden, W.A., "Nonlinear Transient Motion ofDeep Ocean Mining Pipe", JERT, Vol. 103, March 1981, pp. 2-10.[23] Nordgren, R. P., "Dynamic Analysis of Marine Riser with Vortex Excitation",JERT, Vol. 104, March 1982, pp. 14-19.128[24] Finn, L.D., "Dynamic Stress in Offshore Pipelines", Paper No. 1671, ASCE Na-tional Structural Meeting, Cleveland, Ohio, April, 1972.[25] Hall, J.E. and Healey, A.J., "Dynamics of Suspended Marine Pipelines", JERT,Vol. 102, June 1980, pp. 112-119.[26] Hurty, W.C., and Rubinstein, M.F., "Dynamics of Structures," Prentice-Hall, NewJersey, 1964.[27] Clough, R.W., and Penzien, J., "Dynamics of Structures," McGraw-Hill, New York,1975.[28] Key, S.W., "Transient Response by Time Integration: Review of Implicit and Ex-plicit Operators," Advanced Structural Dynamics, Edited by J. Donea, 1978, pp.71-95.[29] Irani, M. B., "Some Aspects of Marine Riser Analysis," Ph.D. Thesis, 1989, De-partment of Mechanical Engineering, Department of Mechanical Engineering, Uni-versity of British Columbia, B. C., Canada.[30] Butenin, N.V., "Elements of The Theory of Nonlinear Oscillations," Blaisdell Pub-lishing Company, 1965, New York. Toronto. London.129Appendix ACOEFFICIENTS OF ANALYTICAL SOLUTIONFor the nonresonance case::P11 =CDADk i2.^f \1y 2 Y2 + 2A 1 1'23 c4^2A2112 Y2idz67 3 (k? — lenw 2 io^M11^M22A i CIDA DI4[Y2Y + Y 3a2 ]dzP12 27 2 (k? — 1c3)w 2 Mii o^1 2^2 2CDADIC1k2 ^L( ^Yl - Y2 )dZA1a1a21/23^A2^.2P13 = 27 2 (k? — k3)w 2 fo^CDADki ^rAil'23a1(e2 d2)^A2 ^3P14 =^ id21C1Z2 71-2(k? knw2 —^ Mu^aliv/22Pis^4CDADA1k2a2 i yL15^71-3(k?^lenw 2 A/_ 11 0 22 ( 1 2e2 + Yi d2 )dzP16CDADA1 47 3 (k? — 1c3)w 2 M 11 Jo 11[2v,2 + 71-vu ,2^1 ,12 7r(d2i^d22 )^y22 7, (e 2,^e22 )^271.),,(Yd+ 2Y1 Y2 7r(d i e i + d2 e 2 )]dzCDA D A,k?P21 — ^ c4I 2 2 )dz47r 2 ( k? — Icpw 2 .....^o Y2(Y12CDAD14^L^ L 2A2 Y2Y2dZ}P22 — 67r 3 (k? — lenw 2^11111Y2(1'12 + 2Y22 aDdz + ^fo 11122 1130CD AD kl k2^.^A2a1 fLP23 = 27^22(k2^k2)(.02[M11 0 Y22 la 1 2^1 d' 4- a2 M22 JoY2Y2dz]4C( vD ADAiklal I LP24^7r 3 ( k? — k 3 )w2^0 ;22 2e2 + Y1d2 )dzCD AD k2 r^ L A 2P25 =^ 1 2 / 2 a 2 e 2 Yi a 2 d2 )dz271.2(k? — kDi.„2 ^Y3 d2 d z]Mil 0 a2 M22^4,7,3(^w2M11^1,2[2vc2^,1,712(d?kC?D — D3 A) P26^7rY22(4^— 27rvu,(Y1 + Y2e2) + 2Y1 Y2 71- (d i + d2e2)]dzBn =- 77111 + a 21 m 2B12 = 2(rnn + a 2 rni2)B13 = ^1CDADIc?^Dial IL1' 172 dz^A22 (}713 2Yi Y 2 c4)dzi67r 3 ( k? — knw 2 Mn o^M22 a loB14 =CDADA2k22 -2^f2 2+ 1 2 a 2 )dz471-2 (k? — kDa i w 2 Al22 JoB15CDADa227r2(14ki k2( Ai /Ak111 A2/M22)10L YiY22dz— — nw 2B16^217111(4^d2 )^277t 12 (4^4rn13(diei^d2 e 2 )CD ADk1 A2 IL I , .1 ,^dia2^/ 2 tz.ZAf22 0271-2 (14— knw 2 M11^o/ 2 2e2^d2 )(Li4CD ADA2a2k2d2^a i 71-3 (1c? — knw2M22^Yi -2(Yid2 + Y2 e2)dzB17 = ^ 131Cf)AD )12B184cc17r3(k? — 1c3)(.4) 2 11122 0^[2vc2^7rvw2 + yi2(d2i ^d2)^7r1,722 (e32^e 2 )- 27rv„,(Yi di + Y2e2) + 2Yi Y27(d 1 e 1 d2 e 2 )]dzB21 = 2 (m21 a 2i m22)B22 = m21 a2 m22B23 =CDADA2k? Yi ( 712 + Y22 c4)dz471-2 (kf — kfla 2 w 2 /1/22 foCD ADIC3 r 2.\1x2 A2B24 —^7 2673(k? knw2 [ M110 1 1 1'2 dz a2M22^(Y13 + 21'1 l' 2 c4 )dz]B25 —CDADCl1k1k2(A2 / Al22 — a1/ A711) IL Y22dz271-2 (14— k2 )w 24CDADA 2 a,B26 — (1273( k?^k2) c.,,2 m22 JoLYi 2(1'e2e2 + d 2 )dzB27 2m21(d2i + d2) + 2m22(4 ++e ^4rn23 (d 1^d2e2 )CDADA2k2 Yi 172 (Y2 e 2 Yi d2 )dz▪ 27 2 (k? — 1c3)w 2 A-12210LCDAD.A2B28[2vc2^7rvw2^71,2 (d23.^d22 )^71/22(612 + „.2)I '2 /4a 2 7 3 (k? — kflw 2 /1122 fo— 27rv i,,(Yi dl + Y2e2) + 2S i Y27r(d1 ei + d2e2)jdzwhere,m11 = — EP(a1a22 + 3a11 a 1a 13m12 = —EP(3ct 1 a 24132a12m13 = — EP(a1a23 + ^ )a l3anm21 = EP(a2a22 + ^ )a2m22 EP(3a2a24 + an)a 2a12m23 — EP(a 2 a 23 ^ )a2EP =8(14 —^Ain Al22(.0 2Di 2^(n 2 — 1)d2 = D2 AlA1el^ 02 (n 2= ^ i — 1)e 2 = —D1 A2AlFor the first resonance case:CDAD ^fLAiY12Y2 + 2A1 Y23 (4. 2A2Y12Y2^jdzPn = 67r 3 (1 — k2 )w 2 Jo M11 M22AiCDADk3^L [y2y2 y23(122P12^] dz= 27 2 (1 — knw 2 Mii o 1P13CDADk2 ^iL(Alaia2)2 +-r A2—^Y22Y )dz27 2 (1 — knw2 J0^16111^a1M2 C12133P14P15=0CDADA1 P16 =^ 171(2Vc2^)dz47 3 (1 — k3)w 2 A/JoCDADA LP21 — ^ 1-20,12^(112 1 -22)dz^41-2 (1 — k3)cv 2 Af^0- CDADI4 fL -^r2P22 =^671-3 (1 — knw 2 0 Mn372(}^21722a)dz foL A2IA:2112Y2dziC D AD k2P23^f L 3^A2a127r 2 (1 — k2),,2 [m^Y2 al ce2dz^v2x72 —^11 01 I 24a2M22 0 LP24 0P25 = 0^CDADA1 ^LP26 =^ 12(2v2 )dz47r 3 (1 — knt.v 2 Mii foB11 = mn ai2 m12B12 = 2(m11 a 2 m12)B13 ^^CDAD^nicei^,2 dz^A2 ^(Y3 + 21 71 Y 2 c4)dz]67-3 (1 — knw2 [ M11 o A/22cti Jo^1CDADA214^yi 712 + Az22 a22JOB14^)aj-= 247 2 (1 — k2 )(1 1 w 2 /1122134CDADa2k2()i^— A20/22)1 . .2B15 —^ dz27r 2 (1 — k2 )w 2B16 = 0B17 = 0B18CDADA 2 ^Yi(2v,2 rv u,2 )dz4a 1 7 3 (1 — k3)w 2 M22B21 2(m21 + 021M22)B22 = m21 a2 m22CDADA2 LB23 = ^ Y( Y-32 + Y22 (4)dz47 2 (1 — 4)c/ 2 4.4) 2 A/22 JoB24 — CDADq^2A102 IL^A2 6730 knw2 ^mil 0 Y22 dza211122^(Y13 + 21;1' 2 2 ciliB25 =CD ADaik2(A21 Al22 - IA111) 10 L 170722 dz27r 2 (1 — knw 2B26 — 0B27 — 0CD AD A2 B28 =^ Y1 (21/c2^7r- vu,2 )dz4a 2 7 3 (1 — k3)(.,, 2 Al22 fowhere,am il = — EP(a a22 + 3 11 )a 1135M12 — — EP( 3a1a241112m13 = — EP(a1a23 + ^ )311 11M21 EP(a2a22 ^ )aza i3 ,M22 = EP( 3a2a24M23 = EP(Ct2a23EP =8(1 — knAliiM22w 2For the second resonance case:P11 -=CDADk 2^L A1Y 2 Y2 + 2A 1 Y23c4 2A2Y12 Y267r 3 (k? — 1)w 2 Jo^M11^M22P12AlCDAD^fL—27r2( k? 1)w2Mii^[Y 2 Y + 173 a 2 k/z1 2^2^ CDADkI ^( 1a1a 2 }?^AZ P13 —270(k? — 1 )(Al 2^11/11 al M22 Y12172 )dzP14 0P15 =0CDADAi P16(22/c2^71- v iu2 )dz47 3 (14— 1)w2M11/a13 )a lalz )azC12136P21^CDADAiki2 ^= 12(1',2^0 21 1'22 )dz47 2 (k? — 0(.4) 2 Mu f)L 2A2Y 2CDAD^r^Al A,^+^ Y dz]P22 = 673(k? — 1)w2 2(Yi -r 412 cz.2)u," fo ^Al22 1 20 M11P23 =CDADk111-23ctia2dz^A2Aa;^l'i2 Y2 dz]271-2 (k? — 1)w 2 [ ._^O^a211/122 0P24= 0P25 = 0CDADA1P26 7-- Y2(2vc2^7rv„,2 )dz471-3 (k? — 1)4.4) 2 /1/11 JuB11 =^ct2im12B12^2(M11^ci2 m12)^CDADIc?^Dial I L2^A2 B13 — (Y3 + 2Y1 Y 2 a)dzi671-3 (1c? — 1)w 2 [^YiY dzM11^2 M22a1 fiL 1LB14 =CDADA2 Y (1 2 Y22 aDdz47 2 (k? — 1)a 1 w 2 11/22CDADa2k1(A1 /M11 — A2/M22 )fLyiY22dzB15 —272(k? — 1)w 2B16= 0B17 = 0137CDADA2 ^Y(2vc ^7rvu,2 )dz4a01-3 (k? – 1 )w 2 M221B18 =^0B21 = 2 (m21 a i2 m22)B22 — m21 a2 77122CDAD A214B23 —( 1712 + Y22 )dz47r 2 (k? – 1 )a 2 w 2 /1/22 foCDAD^2A1a2 f IL "1 1722 dz^A2^(Y13 + 2/ Y 2 a 2 )dz]LB24 =67 3 (k? – 1 )w 2^Mil o^a2M2210B25 =^CDADa i kl(A 2 1 Al22 – 1M11 ) I^2L^.2Yi Y dz272(k? – 0(.4)2^0B26 — 0B27 = UfCD AD A2 B28^ (22/c2^7v,„2 )dz4a 2 7 3 (k? – 1 )w 2 Al22 Jowhere,+ 9 an \M11^EP( al a22 -r^)a 1a13m12 – EP(3a1a24 + ^ )a la12m13 = – EP(ai 6/23 + ^ )a 1m21 EP(a2a22 3a11a2M22an \EP(3a 2 a 24 —)a2m23 = EP(a2a23 a12 )a2C12EP =8(k? – 1)Al11Al224.0 2138
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Dynamic analysis of a marine riser
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Dynamic analysis of a marine riser Guo, Yongzhi 1992
pdf
Page Metadata
Item Metadata
Title | Dynamic analysis of a marine riser |
Creator |
Guo, Yongzhi |
Date Issued | 1992 |
Description | The problem of the dynamic analysis of a marine riser is considered in this thesis. The equations of motion for a marine riser undergoing large three dimensional deflections and rotations are derived using modal discretization approach. The nonlinear expressions for the bending and torsional deformations are obtained for a beam undergoing finite rotations. The elastic strain energy is expressed in terms of these deformations and the Lagrange equations are used to derive the equations of motion valid for large deformations. The three dimensional shear effect based on nonlinear elastic theory is included in the formulation. The dynamic response of a typical marine riser system due to excitations of the ocean waves, current and surface vessel motion is studied. The hydrodynamic loading on the riser is represented by a general form of the Morison equation. The nonlinearity of structure and drag force are completely retained and a numerical procedure in the time domain is developed. The effects of structural and hydrodynamic parameters on motion amplitude in two dimensional case are evaluated. Alternatively, the dynamic response of the simplified riser system is studied by an analytical approach based on Butenin's method. The Butenin's method involves less computational effort and the solution based on this method matches the numerical solutions well when the nonlinearities are small. |
Extent | 4948031 bytes |
Genre |
Thesis/Dissertation |
Type |
Text |
FileFormat | application/pdf |
Language | eng |
Date Available | 2008-08-13 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0080911 |
URI | http://hdl.handle.net/2429/1371 |
Degree |
Master of Applied Science - MASc |
Program |
Mechanical Engineering |
Affiliation |
Applied Science, Faculty of Mechanical Engineering, Department of |
Degree Grantor | University of British Columbia |
GraduationDate | 1993-05 |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-ubc_1993_spring_guo_yongzhi.pdf [ 4.72MB ]
- Metadata
- JSON: 831-1.0080911.json
- JSON-LD: 831-1.0080911-ld.json
- RDF/XML (Pretty): 831-1.0080911-rdf.xml
- RDF/JSON: 831-1.0080911-rdf.json
- Turtle: 831-1.0080911-turtle.txt
- N-Triples: 831-1.0080911-rdf-ntriples.txt
- Original Record: 831-1.0080911-source.json
- Full Text
- 831-1.0080911-fulltext.txt
- Citation
- 831-1.0080911.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0080911/manifest