Skip to main content

ORIGINAL RESEARCH article

Front. Mech. Eng., 28 August 2023
Sec. Solid and Structural Mechanics
Volume 9 - 2023 | https://doi.org/10.3389/fmech.2023.1156721

Dynamic analysis and validation of a multi-body floating wind turbine using the moving frame method

  • 1Institute of Mechanical and Marine Engineering, Western Norway University of Applied Sciences, Bergen, Norway
  • 2Department of Mechanical and Structural Engineering and Materials Science, University of Stavanger, Stavanger, Norway

This research applies the moving frame method (MFM) to the multi-body dynamic analysis of an OC3 phase IV spar buoy with the NREL 5MW turbine. Further, it verifies previous results obtained through numerical comparisons with commercial software. The long-term goal is to lay the foundation for leveraging the MFM to create a self-contained software system for future analyses that can incorporate effects that are more sophisticated, when commercial codes fall short. In this first evidentiary phase, this project treats the floating turbine as a three-bodied system consisting of the platform (platform + tower), nacelle and rotor (hub + blades). Then the paper presents the MFM in a tutorial style—in the context of this problem’s resolution. The paper supplements the multi-body dynamic equations of motion obtained through the MFM with simplified and reduced hydrodynamic, aerodynamic and mooring loads to simulate the translational and rotational response of the floating turbine under various load conditions. The results closely approximate those found in previous work and, in the process, demonstrates MFM’s analytical advantage. Current results capture the coupled responses in all degrees of freedom and gyroscopic effects occurring when the platform pitches with the spinning rotor. The project thus provides an accurate model for the dynamics of the turbine and opens the door to inserting correct advanced hydrodynamics to validate the model further. The work presents simulations for the different load cases through a 3D web page using WebGL and the ThreeJS library. Users may download all software to verify the results. An undergraduate student conducted the work alone, demonstrating the ease of implementation of the MFM.

1 Introduction

The world’s increasing need for energy makes renewable energy sources attractive. While the demand for renewable energy increases, so does the resistance to onshore wind farms due to the environmental impact, which puts certain species at risk WWF (2022). Offshore wind farms will lessen the environmental impact regarding permanent structures’ installation and yield more favourable wind conditions Burton et al. (2011).

Floating offshore wind turbines (FOWT), however, experience harsher environmental loads and also increased freedom to translate compared to fixed installations. This implies that predicting the movement of the FOWT is a challenge unto itself. To maximise their potential, continually improved analyses for floating offshore wind farms must be implemented.

Currently, there exists reliable commercial software to conduct analyses of such structures. However, as engineers continue to infuse analyses with artificial intelligence, it is advisable to rely on such legacy software alone and, instead, supplement analyses with A.I. and digital twins. A digital twin is a representation of an industrial asset in operation that includes the asset’s condition and relevant historical data. When supplemented with AI, such models can vastly improve the asset’s operation; in this case, includes the modelling of dynamics response. A judicious use of inter-process communication algorithms, in which a digital twin can access a database of both previously solved problems coupled with extant weather conditions; and this can speed the computation for near real-time simulations. In fact, this work is currently underway in a related project of crane motion on ships using the method used and discussed herein.

The Moving Frame Method (MFM) can ease the process of obtaining the equations of motion Storaas et al. (2019) and can substantially reduce simulation time when implemented without the computational overhead of commercial software. The MFM utilizes Lie algebra and Cartan’s notion of Moving frames to extract the equations of motion of both single- and multi-body systems as well as 2D and 3D kinetics in a consistent notational matter Impelluso (2017).

This research constructs the equations of motion using the MFM and implements some of the environmental loads. It will focus on the OC3 Phase IV spar-buoy platform and the NREL 5MW reference turbine due to extensive research already performed on such platforms and turbines. The goal is to replicate previously found results and also to lay the foundation for further use of the MFM when addressing complex multi-body dynamic systems. The authors selected the FOWT as such an implementation due to its nature. It is a complex machine, and experiences equally complex loading. Furthermore, the industry is implementing AI for many of its optimization problems. This makes the MFM method an interesting approach to use. Essentially, this project endeavors to demonstrate advanced computational techniques, such as how the MFM enables easy transitions between scenarios where the tower-nacelle-rotor connection is considered as a rigid body, and, later, coupled with a finite element process, a deformable body.

To understand the method, the reader may review two papers Impelluso (2017); Murakami (2013). In the first half of this paper, we will launch directly into the method in a tutorial style.

2 The moving frame method

The Moving Frame Method (MFM) represents a paradigm shift in the foundation of the discipline of multi-body dynamics.

2.1 MFM background

The MFM exploits the philosophy of Elie Cartan who advocated that all objects be endowed with a moving reference frame Cartan (1986). To relate the moving frames in a multi-body system (one for each link), the MFM utilizes a compact notation, along with the group theory and associated algebra of Sophus Lie but distilled to the simplicity of matrix multiplications Frankel (2011).

In addition, the MFM for multi-bodies adopts a revised approach to the Principle of Virtual work, in which a proper variation of the angular velocity is established. While Wittenburg Wittenburg (2007) advocated the use of virtual power to embed driving forces and motor torques into the variation principle, the MFM adopts a restriction on the variation of the angular velocity. This restriction derives from both the variation of the frame and the variation of the coordinates; and as first advocated by Holm Holm (2008).

The MFM is not merely an alternate method of analysis. Rather, it represents a paradigm shift in a discipline—one which empowers the moving reference frame over the inertial. It eschews the notion of fictious forces that derive from such frames, and predicating itself on the reality that everything moves. The MFM renders a notation that is consistent across the facets of dynamics. It enables 2D, 3D, single and multi-body dynamics to be represented with the same consistent notation. It needs no validation as it produces the same equations of motion as the traditional approach to dynamics, but more expeditiously, and with a suite of equations that more more amenable to coding through its matrix-based formulation.

A pedagogical analysis of the method and an assessment of student learning has been presented previously Impelluso (2017) for single body, 3D dynamics. Regarding research, the MFM has been used to model ROV motion Austefjord et al. (2018), friction Rykkje et al. (2021), gyroscopic wave energy converters Korsvik et al. (2019), a knuckle boom crane Eia et al. (2019), and ship stability Flatlandsmo et al. (2019).

In this project, we apply the MFM to study the dynamics of a wind turbines. In this way, we demonstrate the power of the MFM and its ability for easy deployment in the analysis of floating wind turbines. For in this era of smart machines, we must consider recoding legacy software if we are to imbue existing software with artificial intelligence.

2.2 Applying the moving frame method

Figure 1 illustrates how each body of the FOWT will be endowed with its own reference frame at the body’s center of mass. The combined body of the platform and the tower will hold the first frame from which the inertial frame will be deposited. In other words, the first body’s frame is primary and the inertial frame is merely a decal to be peeled off at the onset of the analysis to compute position data. The second moving reference frame is attached to the nacelle, and the third is attached to the combined body of the hub and the blades. All frames will be aligned with the centre of mass of the respective bodies at the start time. In the following section, the kinematics of each moving frame will be described. This will lay the foundation for obtaining the equations of motion. We now use the Special Euclidean Group, SE(3), to obtain each component’s angular velocity and translational velocities systematically. However, to ease understanding, the reader is advised to consult the primary results of each of the following three subsections, as per Table 1, below:

FIGURE 1
www.frontiersin.org

FIGURE 1. Turbine with frames attached to each body.

TABLE 1
www.frontiersin.org

TABLE 1. Kinematic results.

2.3 Kinematics of the platform, body 1

The platform will have all rotational and translational freedoms.

• The rotation of the platform—pitch, yaw and roll—from the inertial frame, is found by post-multiplying the inertial frame by a rotation (a member of the special orthogonal group, SO(3).

• The translation of the platform—surge, heave and sway from the inertial frame is found by post-multiplying the inertial frame by translation coordinates.

e1t=eIR1t(1)
rc1t=eIxc1t(2)

Obtaining these rotational and translational coordinates in the inertial reference frame will be the main focus of later parts of the analysis. However, at this time, we proceed beyond SO(3) and adopt the Special Euclidean Group, SE(3), which will enable us to combine both motions in one data structure.

Next, the frame and its position are group in a structure called the frame connection. Below, both the inertial frame connection and the moving frame connection (the former asserts the origin) are displayed:

eI0(3)
e1trc1t(4)

The two frame connections are related through a frame connection matrix, defined implicitly by Eq. 5. Notice that the expression below recapitulates the two equations that commenced this section.

e1trc1t=eI0E1t=eI0R1txc1t0T1(5)

The frame connection matrix, E(1)(t), is a member of the Special Euclidean Group, SE(3), and carries associated group properties and algebra. As a member of SE(3) the inverse can be found analytically:

E1t1=R1TtR1Ttxc1t0T1(6)

As an aside, note that the superscript (1) implies an absolute frame connection matrix to relate the inertial and moving frame connections.

Next, we obtain the rate of change of the frame connection matrix as follows:

ė1tṙc1t=eIt0Ė1t(7)

Taking the derivative of the frame connection matrix Ė(1)(t), yields:

Ė1t=Ṙ1tẋc1t0T0(8)

Since the inverse of the frame connection matrix is analytically known through the algebraic properties of SE(3), the rate of change of the frame connection matrix can be expressed in terms of the same frame connection matrix:

ė1tṙc1t=e1trc1tE1t1Ė1t(9)
ė1tṙc1t=e1trc1tΩ1t(10)

The two equations above implicitly define Ω(1)(t). By multiplying the inverse of the frame connection matrix and its rate, we obtain:

Ω1t=ω1tR1Ttẋc1t0T0(11)

From this the foundation has been laid to extract the translational and rotational information in the following steps:

ω1t=R1TtṘ1t(12)

As previously mentioned the angular velocity matrix in its skew-symmetric form is used to express the time rate of change of the body frame in terms of the body frame:

ė1t=e1tω1t(13)

From Eq. 11 the translational velocity information of the body in the moving frame can also be extracted:

ṙc1t=e1tR1Ttẋc1t(14)

Which is easily referred back to the inertial frame using the inverse frame relation:

ṙc1t=eItẋc1t(15)

This concludes the kinematic of the platform and the next section will continue through the tree structure of the platform, nacelle and blades. From here, the same structures will apply, however, they will be formulated as relative frame connections.

2.4 Kinematics of the nacelle, body 2

The yaw bearing motor enables the rotation of the nacelle about the tower’s vertical axis. As a simple first pass, this study assume the motor and Nacelle as one body. To construct the frame connection matrices, one may progress from the first body to the second as follows:

• Translate to yaw bearing, henceforth known as point J, using the turbine frame. Ultimately, for simplicity, the analysis will assert that the second body is directly above the first: s2(J/1)(t)=s1(J/1)(t)=0. However, for now, the analysis leaves it as general coordinate form:

sJ/1t=e1tsJ/1t=e1ts1J/1ts2J/1ts3J/1t(16)

• Next, to reach the second frame, we allow the second body to rotate from the first, using a member of SO(3) to affect the rotation, about the vertical axis:

R2/1t=cosϕsinϕ0sinϕcosϕ0001(17)

• Finally, we translate to the center of mass of the nacelle. For simplicity it can be assumed that we progress directly outward: s2(2/J)(t)=s3(2/J)(t)=0.

s2/Jt=e2ts2/Jt=e2ts12/Jts22/Jts32/Jt(18)

With the above mentioned rotational and translational coordinates, the relative frame connection matrix - that relates the platform frame connection to the nacelle frame connection can be obtained as the non-commutative product of the following matrices (which can be stated, in words, as translate, rotate, translate):

E2/1t=I3x3sJ/1t01R2/1t001I3x3s2/Jt01(19)

Performing the matrix multiplication yields:

E2/1t=R2/1tR2/1ts2/Jt+sJ/1t01(20)

The relative frame connection relationship is now complete:

e2trc2t=e1trc1tE2/1t(21)

Note that the superscript on E(2/1)(t) implies a relative frame conenction matrix, between two moving frame connections.

However, by taking advantage of the platform’s frame connection matrix, Eq. 5, the absolute frame connection of the nacelle is found using the Group property of closure under matrix multiplications:

e2trc2t=eI0E1tE2/1t(22)

Define:

E2t=E1tE2/1t(23)

Using what was previously found it can be shown that:

E2t=R2tR2ts2/Jt+R1tsJ/1t+xc1t0T1(24)

Where R(2)(t) in the upper left of E(2)(t) is the absolute rotation matrix of the nacelle frame, defined in the following manner:

R2t=R1tR2/1t(25)

The term located in the upper right quadrant of E(2)(t) describes the nacelle’s center of mass location from the inertial frame.

As in the previous section the inverse of the absolute frame connection matrix is obtained:

E2t1=R2TtR2TtR2ts2/Jt+R1tsJ/1t+xc1t01(26)

The time rate of change of the absolute frame connection matrix:

Ė2t=Ṙ2tṘ2ts2/Jt+Ṙ1tsJ/1t+ẋc1t00(27)

Finally, using the data structure of SE(3) the time rate of change of the second frame connection in terms of the second frame connection is:

ė2tṙc2t=e2trc2tE2t1Ė2t(28)

Define:

Ω2t=E2t1Ė2t=Ω112tΩ112t00(29)

Expanding on the terms in the Omega matrix:

Ω112t=R2/1TtR1TtṘ1tR2/1t+R1tṘ2/1t(30)
Ω122t=R2/1TtR1TtṘ2ts2/Jt+Ṙ1tsJ/1t+ẋc1t(31)

From Eq. 30 the absolute angular velocity matrix can be extracted and used to find the angular velocity vector:

ω2t=R2/1Ttω1t+ω2/1t(32)

Using Eq. 31 the translational velocity of the second body expressed in the inertial frame can be found:

ẋc2t=R2tsc2/JtTω2t+R1tscJ/1tTω1t+ẋc1t(33)

For the sake simplicity, the recursive relationships for only rotations, SO(3) and the combined rotation and displacement SE(3) are presented here.

ωαt=Rα/α1Ttωα1tRα/α1t+Rα/α1TtṘα/α1t(34)
Ωαt=Eα/α1TtΩα1tEα/α1t+Eα/α1TtĖα/α1t(35)

As we conclude the nacelle kinematics it must be emphasized that the nacelle will not yaw in the current case studies, and this will be addressed in later parts of this projec. Despite this, we have installed the possibility to insert such a motor-driven yaw. Also the relative angular velocity between the nacelle and the turbine is analytically known to involve only one rotation since it is driven by a motor torque. It is written as:

ω2/1t=ϕ̇001=ϕ̇e3(36)

Recognizing that the column, above, is the uskewed form of the basis of skew symmetric angular velocity matrices for a single rotation:

e3=010100000(37)

Combing Eq. 36 and 32 the absolute angular velocity vector can be stated:

ω2t=R2/1Ttω1t+e3ϕ̇(38)

2.5 Kinematics of the blades, body 3

In this project, the three blades and the hub will be considered one body; however, nothing prevents future research from extending this and modeling each individual blade.

The rotation will be induced by the wind speed, and a motor (the motor will be applied to slow the spinning of the blades down and generate power).

Since the center of mass of the blades lies coincident with the spin axis of the blades there is no need to translate after accounting for the rotation of the blades.

• Translate from the nacelle CM to the joint (which incidentally coincides with the CM of the blades)

sc3/2t=e2tsc3/2tsc13/2tsc23/2tsc33/2t(39)

• Rotation of the third frame from the second frame will happen about the 1 axis and the rotation matrix is shown:

R3/2t=1000cosψsinψ0sinψcosψ(40)

The non-commutative building blocks of the relative frame connection matrix follows:

E3/2t=I3x3sc3/2t01R3/2t001(41)

Through matrix multiplication it is known:

E3/2t=R3/2tsc3/2t01(42)

Which is used in the context of frame connections in the following way:

e3trc3t=e2trc2tE3/2t(43)

The absolute frame connection matrix is found by pre-multiplying this result by the absolute frame connection matrix for the nacelle:

E3t=E2tE3/2t(44)

The result of this matrix multiplication is shown below, however every step is not shown:

E3t=R3tR2tsc3/Jt+R1tscJ/1t+xc1t01(45)

Since the translation for the yaw bearing to center of mass of the nacelle happens in the same frame as the translation from the nacelle to the center of mass of the blades it is known that:

sc3/Jt=sc3/2t+sc2/Jt(46)

The time rate of change of the frame connection of the rotor frame in terms of the rotor frame can also be found:

Ω3t=E3t1Ė3t(47)
Ω3t=Ω113tΩ123t00(48)

Expanding on the terms in the matrix:

Ω113t=R3TtṘ3t(49)

The angular velocity vector can be extracted:

ω3t=R3/1Ttω1t+R3/2Ttω2/1t+ω3/2(50)
Ω123t=R3/2Ttsc3/JtTω2t+R3/1TtscJ/1tTω1t+R3Ttẋc1t(51)

Which leads to the translational velocity of the rotor expressed in the inertial frame:

ẋc3t=R2tsc3/JtTω2t+R1tscJ/1tTω1t+ẋc1t(52)

Before delving deeper into the solution, note for the sake of simplicity, that the relative angular velocity of the rotor from the nacelle is analytically known to involve only one rotation. Thus, this rotation coordinate can be infused into the analysis in a simpler way as follows.

ω3/2t=ψ̇100=ψ̇e1(53)

Similarly to the nacelle kinematics the column is the unskewed form of the basis of skew symmetric angular velocity matrices fora single rotation about the shared first axis.

e1=000001010(54)

And therefore the angular velocity vector can be restated as:

ω3t=R3/1Ttω1t+R3/2Ttϕ̇e3+ψ̇e1(55)

This concludes applying the frames to the turbine, in the next section the generalized essential coordinates will be examined.

In accordance with Table 1, we summarize the results we have obtained. The work performed yielded the following Cartesian result presented below. The first two of which are tautological:

ẋc1t=ẋc1t(56)
ω1t=ω1t(57)
ẋc2t=R2tsc2/jtTω2t+R1tscJ/1tTω1t+ẋc1t(58)
ω2t=R2/1Ttω1t+ω2/1t(59)
ẋc3t=R2tsc3/JtTω2t+R1tscJ/1tTω1t+ẋc1t(60)
ω3t=R3/1Ttω1t+R3/2Ttω2/1t+ω3/2t(61)

2.6 Generalized coordinates

The next step is to reformulate these expressions and extract the Cartesian coordinates as a function of a minimal set of generalized coordinates. Since there are three bodies the Cartesian coordinate rates will be a 6 × 1 column vector holding the translational and rational velocity of each respective body, where each element in the column will be 3 × 1 column. The full Cartesian rate column will be 18 × 1:

Ẋt=ẋc1tω1tẋc2tω2tẋc3tω3t(62)

In the previous equation, the first row consits of surge, sway and heave; the second row consists of roll, pitch and yaw of the turbine body: this is followed by the Cartesian rates for body 2 and 3. Meanwhile the rate of change of the generalized essential coordinates:

q̇t=ẋc1tω1tϕ̇ψ̇(63)

In Eq. 63, the first two rows are as before, but the third row is one coordinate (the turning of the nacelle) and the last row, also one coordinate, is the spin of the blades.

The Cartesian rates are related to the time rate of change of the generalized essential coordinates through the B-matrix as follows. The system is linear in that the generalized rates are readily extracted:

Ẋt=Btq̇t(64)

For the FOWT the B-matrix will be a 6 × 4 matrix in compact notational form:

Bt=I303x303x103x103x3I303x103x1I3B32B3303x103x3R2/1Tte303x1I3B52B5303x103x3R3/1TtR3/2Tte3e1(65)

Where:

B32=R2tsc2/JtTR2/1Tt+R1tscJ/1tT(66)
B33=R2tsc2/JtTe3(67)
B52=R1tscJ/1tT+R2tsc3/JtTR2/1Tt(68)
B53=R2tsc3/JtTe3(69)

2.7 Hamilton’s principle and kinetics

To obtain the equations of motion for the FOWT, Hamilton’s Principle is applied. To this end, the Lagrangian is established in a structured form to account for all masses. The Lagrangian describes the states of a dynamic system as the difference between kinetic and potential energy for conservative forces. By extending Hamilton’s Principle with the Principle of Virtual Work, we can account for non-conservative forces such as the applied motors and dissipative fluid forces.

L̂=K+W(70)

Both translational velocity and angular velocity contribute to the kinetic energy of dynamic multi-body-systems. However, this energy will be expressed using the linear and angular momentum for each body, α = 1, 2, 3, stated here (using the mass mα, and mass matrix of inertia, Jc(α), of each body) as:

Lcα=eItLcα=eImαẋcαt(71)
Hcα=eIHcα=eItJcαωαt(72)

Next, the goal is to group the previous two equations in one structured form, as follows. A generalized mass matrix M is constructed. It consists of the alternating masses and mass moment of inertia of the bodies in the multi-body system: (1)-Platform, (2)-Nacelle, (3)-Rotor. Each entry on the diagonal below is, itself a 3 by 3 matrix. The impact of a sparse matrix will be handled in the Matlab code.

M=m1I3000000Jc1000000m2000000Jc2000000m3000000Jc3(73)

Using the generalized mass matrix, the generalized momenta can be constructed:

H=Mẋc1tω1tẋc2tω2tẋc3tω3t=MẊt=Lc1Hc1Lc2Hc2Lc3Hc3(74)

The kinetic energy of the entire multi-body system can be stated as:

K=12ẊtTH(75)

With this expression for the kinetic energy and an understanding of the applied loads, we obtain the equations of motion through variational methods. However, there will be a difficulty in formulating the variation of the angular velocity and a digression must be made.

The variation of the generalized velocities needs to be obtained to take the variation of the Lagrangian. The variation of any frame connection matrix has been defined by Murakami at UCSD Murakami (2013) and was obtained independently by Holm (2008).

δΠα=δπαtRαTtδxcαt0T0(76)

Where δπα(t) is the variation of the angular displacements defined as:

δπαt=RαTtδRαt(77)

Note that the unvaried form does not exist. This term above, is merely a definition that will be structurally consistent with related terms in this analysis. Essentially, the rotation matrix is varied, δR(α)(t), and pulled back to an inertial frame by post-multiplication of the transpose of the original form: R(α)T(t). Next the virtual generalized displacement matrix column, δX̃(t) is presented:

δX̃t=δxc1tδπc1tδxc2tδπc2tδxc3tδπc3t(78)

The commutativity of the time differentiation and the variation of the frame is essential during the integration of the Lagrangian. The restriction on the variation of the angular velocities is proved below:

δωαt=ddtδπαt+ωαtδπαt(79)

The variation of the translational velocities:

ddtδxcαt=δẋcαt(80)

Next the two equations above are consolidated into one form. This is done by construction another sparse matrix system. D is a skew symmetric sparse matrix used to account for the restriction on the variation of the angular velocities:

D=0000000ω1t0000000000000ω2t0000000000000ω3t(81)

Using the above equation the variation of the velocities can be stated as:

δẊt=ddtδX̃t+DδX̃t(82)

The sparsity of the D matrix leaves the variation of linear velocities in standard forrm, yet handles the restriction on the variation of the angular velocities. In Eq. 65 the relationship between the generalized Cartesian velocities and the generalized essential velocities is presented, and the same can be shown for the variation of the generalized Cartesian displacements and the essential generalized displacements:

δX̃t=Btδqt(83)

The virtual work done by all the forces, conservative and non-conservative can be stated as:

δW=δX̃tTQt(84)

Where Q(t) holds all the applied forces and moments:

Qt=Fc1ItMc1tFc2ItMc2tFc3ItMc3t(85)

Recollecting Eq. 83 and 84 the generalized forces and a new expression for the virtual work can be stated:

F*t=BtTQt(86)
δW=δqtTF*t(87)

Varying the Lagrangian will now give the following expression:

δL̂=δẊtTMẊt+δqtTF*t(88)

The Lagrangian can now be integrated noting that the boundary condition of no variation at the boundaries must be met:

t0t1δẊtTMẊt+δqtTF*tdx(89)

Taking advantage of integration by parts and the boundary conditions it can be shown that the equation above yields the following:

BtTMBtq̈t+BtTMḂt+DtMBtq̇t=F*t(90)

Hamilton’s Principle extended by the Principle of Virtual Work has now provided the matrix form equations of motion for the FOWT.

2.8 Reconstruction of the rotation matrix

The platform is free to rotate in pitch, roll and yaw directions and therefore the rotation matrix will be unknown and must be reconstructed from the calculated angular velocity of the platform. If the angular velocity is constant, it can be shown that:

Ṙ1t=R1tω1t(91)

By assuming a constant angular velocity and exploiting a mid-point integration scheme for each time step in the chosen numerical method to solve the equations of motion the reconstruction formula can be used:

Rt+Δt=RtexpΔtωt+Δt/2(92)

When the Cayley Hamilton Theorem is applied to a skew symmetric matrix, it is known that:

exptω0=I+ω0ω0sintω0+ω0ω021costω0(93)

2.9 Hydrodynamics

Solving the hydrodynamic problem in general requires knowledge about two different aspects of the fluid; velocity and pressure Perez and Fossen (2007). The overarching Navier-Stokes equations are unnecessarily too advanced for this first pass analysis. However, by assuming an inviscid, incompressible, irrotational fluid flow the velocity can be derived by taking the gradient of a potential function Faltinsen (1990). With the formulations presented below the MFM will be used to relate the hydrodynamic loads to the FOWT.

To assess the hydrodynamic loads that act upon the platform several assumptions and simplifications must be made. The most important assumption is that the hydrodynamic problem can be linearised. This assumption is mostly valid for deep water sea states where it is implied that the wave elevation is relatively small compared to the wave length which in turn eliminates breaking of waves in most cases Jonkman (2007).

By means of linearisation, the superposition principle can be applied, and the true linear hydrodynamic equation Ogilvie (1964) is presented below:

FiPl=FiH+FiW0tKijtτq̇jτdτ+FiD+FiLAijq̈j(94)

• Hydrostatics: FiH is the resultant force and moments from the body’s displacement, altering the position of the buoyancy center and the amount of displaced fluid. It is calulated by summing the pressure acting on the submerged body.

• Wave Excitation Loads: FiW is the forces and moments from waves. These forces are obtained using a corresponding wave spectra and the normalized vector of wave excitation loads which, in turn, is calculated from the geometry of platform by means of numerical software.

• Radiation loads: The convolution integral term represents the energy that dissipates as the platform oscillates in the fluid.

• Viscous Drag: FiD represents the loads from viscous drag due to the platform velocity relative to the fluid velocity.

• Mooring lines loads: The combined forces and moments from all three mooring lines are encapsulated in the FiL term. A quasi-static model will be used to approximate the mooring loads for various displacements of the FOWT Grindheim (2022).

Aij() is the added mass component at infinite frequency. The translating and rotating platform will also need to displace the surrounding water which is known as added mass and is calculated in the frequency domain using a numerical panel method with software such as Wamit.

2.10 Implementing the model in a simulation

The numerical integration used in this work does not exploit Matlab’s Runge-Kutta 4th order (RK4) solver. The authors developed their own. The Matlab solver does not come equipped with the necessary call out to a function to establish the rotation matrix for the base frame (in this case, the turbine body). For while this work uses Runge-Kutta-4 to integrate, it must deal with one other issue, in tandem: there is no standard rotation matrix for the turbine body (it may include roll, pitch and yaw). Consequently, the code cannot state an explicit rotation matrix and its time derivative for ease of updating the B-matrix between times steps. Thus, the approach used was to use the Cayley-Hamilton theorem, and the Rodriguez formula (the latter of which derives from Cayley-Hamilton), to update the rotation matrix from the angular velocity matrix, between the times steps of RK4.

The code is kept as general as possible and flags are implemented to handle specific scenarios that are simulated. Through the MFM the equations of motion were extracted, recall Eq. 90, and they can be written in compact form:

M*q̈t+N*q̇t=F*T(95)

Where:

M*=BtTMBtN*=BtTMḂt+DtMBt(96)

The process will start with calculating the forces acting upon the FOWT at the current time step before initializing the Runge-Kutta method. Between every partial step of the Runge-Kutta method the rotation matrix for the platform that holds the pitch, roll and yaw displacements must be reconstructed using 93.

Once the essential generalized velocities are approximated the generalized coordinates are also approximated through a simple mid-point integration scheme explained further in the next section. The Cartesian coordinates are then calculated and saved for the visual 3D simulation.

A list of pertinent data used in the Matlab code is presented below:

• Current time

• Time step

• Generalized essential velocities

• The rotation matrix describing the full 3D rotation of the FOWT

• Angular velocity matrix of the platform

• Centre of mass location of the platform, nacelle and rotor

• Mass of the individual bodies

• Mass moment of inertia matrix for each body

• Forces and moments from external loads

Rearranging Eq. 95 and discretizing it:

q̈n+1=Mn*1Fn*Nn*q̇n(97)

The Runge-Kutta 4th-order method with the above equation reduces to:

k1=Mn*1Fn*Nn*q̇nk2=Mn+0.5Δt*1Fn+0.5Δt*Nn+0.5Δt*q̇n+12k1Δtk3=Mn+0.5Δt*1Fn+0.5Δt*Nn+0.5Δt*q̇n+12k2Δtk4=Mn+Δt*1Fn+Δt*Nn+Δt*q̇n+k3Δtq̇n+1=q̇n+Δtk1+k2+k3+k416(98)

A simple mid-point integration scheme is implemented to calculate the general essential coordinates:

qn+1=qn+q̇n+q̇n+12Δt(99)

For the load cases which where simulated in this paper a fixed step integration was used. The initial conditions used in the various load cases were based upon the work done in El Beshbichi et al. (2021a).

A complete system description of the structural properties of the floating wind turbine can be found in Table 2.

TABLE 2
www.frontiersin.org

TABLE 2. Structural data Jonkman (2010).

2.11 Traditional methods of validation

Currently, tools such as Bladed by DNV, SIMO/RIFLEX by Marintek, HAWC2 by DTU and FAST by NREL are used to analyze and validate FOWT design Robertson et al. (2014); Cordle and Jonkman (2011); dnv (2023); Jonkman et al. (2010). These are fully coupled dynamic tools which account for hydrodynamic and areodynamic loads, structural dynamics, mooring and control systems in their prediction of the system response. These methods are known as the aero-hydro-elastic-servo fully-coupled approach Cruz and Atcheson (2016) and while accurate in their predictions can be time-consuming to use El Beshbichi et al. (2021b). These methods also assume a generic FOWT design which limit their usability in exploration of new design which differs significantly from the traditional FOWT.

3 Results

To facilitate comparison between existing validated work, this section will present some of the results from individual load cases. A brief overview of the load cases can be seen in Table 3 The response of the FOWT for free decay surge and pitch (Load Case 1.4 El Beshbichi et al. (2021a)) scenario as well as irregular wind and wave response [Load Case 5.2 El Beshbichi et al. (2021b)].

TABLE 3
www.frontiersin.org

TABLE 3. Load cases, based on El Beshbichi et al. (2021b).

3.1 Free decay surge

The platform is initially displaced 20 m in the surge direction. It is released and the response is simulated accounting for hydrostatic, mooring, drag and damping loads. The results that are presented shows the response of the platform in the surge (Figure 2) and pitch (Figure 3) direction. The results show good correlation with previous findings but there are discrepancies in the amplitude and natural period due to the simplified loads which are applied in these simulations.

FIGURE 2
www.frontiersin.org

FIGURE 2. Surge response of the FOWT from Free Decay Surge.

FIGURE 3
www.frontiersin.org

FIGURE 3. Pitch response of the FOWT from Free Decay Surge.

3.2 Free decay pitch

The platform is tilted to 10° in the pitch direction. It is also translated such that the sea-water line is not displaced in the surge direction. The results presented show the response in the surge direction (Figure 4) and pitch direction (Figure 5). As in the previous section the results are not an exact replica but correlate in orders of magnitude and frequency. The differences are readily explained by the simplified implementation of loads.

FIGURE 4
www.frontiersin.org

FIGURE 4. Surge response of the FOWT from Free Decay Pitch.

FIGURE 5
www.frontiersin.org

FIGURE 5. Pitch response of the FOWT from Free Decay Pitch.

3.3 Irregular wave

For this load case an irregular wave spectra with significant wave height of 6m and peak spectral period is 10s. A turbulent wind field is also created with a mean wind speed of 10m/s. Below the response in the surge direction (Figure 6), pitch direction (Figure 7) and the yaw direction (Figure 8) which is induced from the gyroscopic effect of the pitching FOWT as the rotor spins.

FIGURE 6
www.frontiersin.org

FIGURE 6. Surge response of the FOWT LC5.3 El Beshbichi et al. (2021a).

FIGURE 7
www.frontiersin.org

FIGURE 7. Pitch response of the FOWT LC5.3 El Beshbichi et al. (2021b).

FIGURE 8
www.frontiersin.org

FIGURE 8. Yaw response of the FOWT LC5.3 El Beshbichi et al. (2021a).

4 Discussion

The animation is presented on this web site: https://olemg.github.io/FinalSim/. The reader is free to download the Matlab code to reproduce the results.

This project deployed the MFM to analyze the dynamic response of a floating wind turbine. This work has demonstrated how equations of motion for the analysis are easily obtained with the MFM, and in a form suitable for numerical solution by standard numerical methods. The focus of this project was not a more advanced and detailed analysis, but a validation of previous work, with a method that is more easily deployed. Using the MFM, analysts need not turn to commercial software and potential “black box” functions, but are able to create the entire software system themselves.

The results were compared to previous work done on the same FOWT El Beshbichi et al. (2021b). Figure 9 shows some of the free decay results which were found in El Beshbichi et al. (2021a). The discrepancies between this analysis and previous work are relatively simple and readily explained. The obstacles were the wave radiation, damping coefficients and related issues on the location of the center of mass. Seeing as these variables were not properly adjusted to fit this simulation results that where similar in order of magnitude and frequency were deemed sufficient. These issues do not detract from the mathematical simplicity of the MFM, nor its power.

FIGURE 9
www.frontiersin.org

FIGURE 9. Free decay response taken from El Beshbichi et al. (2021b). “Present Model” refers to the model in that article. The heave response in the figure was removed. Licensed under CC by 4.0 (https://creativecommons.org/licenses/by/4.0/).

A validation must be conducted on the best way to solve the wave radiation damping state-space model and make it account for all couplings in all DOF’s. It could also be of interest to use the MFM to model the delta connection (crowfoot) of the mooring lines correctly such that there is no need for the added yaw spring stiffness.

The discrepancies are easy to overcome with greater attention paid to damping factors by a specialist in fluid mechanics. The primary goal of this project was to demonstrate the MFM and to create the overarching model so that such more refined analyses can be conducted. By revisiting the hydrodynamic coefficients such as the added mass matrix, damping matrix and this analysis could be more accurate (or, it might turn out that previous analyses could be improved; or both).

The results presented in this project do not precisely replicate what was previously found but this was expected due to the previously mentioned simplifications and assumptions (all of which can be revisited, improved, or, if found to be faulty, corrected). Damping can easily be introduced into this model. Once the form of the damping (say, due to wave viscosity) is assessed, one then adds it to the force vector, prior to multiplying by B-star to ensure that is involved as a generalized force.

However this project demonstrates that the moving frame method is capable of complex multi-body dynamic analysis.

5 Conclusion

There were two overarching facets to this work.

The first goal was to commence the construction of a software system to model the dynamics of floating wind turbines. In this regard, initial results have shown to be promising. The overall results capture the response of the floating turbine under various conditions, and they correspond to results predicted by commercial software. While there is some minimal discrepancy, the authors attribute this to an incomplete accounting of viscous effects and wave radiation. Ongoing work will rectify this in the next round.

The second goal of this work aligns with the philosophy of the American naturalist philosopher, Henry David Thoreau who opined, with concern, that “Man has become a tool of his tools.” In that spirit, this work has demonstrated that, with the Moving Frame Method, dynamical analysis is reduced to one equation that is readily implemented through a combination of a simple Runge-Kutta algorithm in a matrix-based environment. Furthermore, through WebGL and its interface, ThreeJS, the work is visualized through two-dimensional graphs, and with three-dimensional web pages (viewable on cell phones). In fact, the entire coding was conducted as a Bachelor project under the guidance of one professor who provided judicious advice on hydrodynamics, and another professor who taught the student the Moving Frame Method. In this way, future researchers need not rely on large-scale legacy software and readily implement analyses dedicated to the focus of their research.

6 Future work

For one of the authors, the future work consists of the dissemination of the MFM in the undergraduate curriculum. The undergraduate discipline of dynamics has become corrupted. There is, today, too heavy a reliance on the inertial frame, too much of a focus on 2D planar dynamics, and an ill-thought use of energy principles. Thus, students never learn the nuances and beauty of the discipline when focused on planar problems. Some schools split the assorted topics into other courses, further disintegrating the discipline, undermining a complete understanding, and giving rise to the disengagement many students now feel. The MFM brings reason and simplicity to dynamics, while uplifting it to three dimensions and solving advanced 3D problems at the undergraduate level. It leverages energy principles not as back-of-envelope confirmation of solutions in 2D (as is done in most textbooks today), but in its rightful place: Hamilton’s Principle and variational calculus. This coming fall, the method will be taught in three universities. The author invites any reader interested in this method to contact him to discuss access to all pedagogical content.

Data availability statement

The raw data supporting the conclusion of this article will be made available by the authors, without undue reservation.

Author contributions

YX guided on the hydrodynamic, wind and mooring loads. TI guided on the Moving Frame Method. O-MG developed the code. All authors contributed to the article and approved the submitted version.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

References

Austefjord, K., Hestvik, M., Larsen, L.-K., and Impelluso, T. (2018). Modelling subsea rov motion using the moving frame method. Int. J. Dynam. Control 7. V04AT06A010. doi:10.1115/IMECE2018-86191

CrossRef Full Text | Google Scholar

Burton, T., Jenkins, N., Bossanyi, E., Sharpe, D., and Graham, M. (2011). Wind energy handbook. 3. New Jersey, United States: John Wiley & Sons.

Google Scholar

Cartan, E. (1986). “On manifolds with an affine connection and the theory of general relativity,” in Monographs and textbooks in physical science (Bibliopolis).

Google Scholar

Cordle, A., and Jonkman, J. (2011). “State-of-the-art in floating wind turbine design tools,” in Proceedings of the International Offshore and Polar Engineering Conference, Beijing, China, June 20-25, 2010.

Google Scholar

Cruz, J., and Atcheson, M. (2016). Floating offshore wind energy: The next generation of wind energy. Cham: Springer Science & Business Media. doi:10.1007/978-3-319-29398-1

CrossRef Full Text | Google Scholar

DNV G1 (2023). enWind turbine design software | Bladed. Norway: DNV.

Google Scholar

Eia, M. E., Vigre, E. M., and Rykkje, T. R. (2019). “Modeling a Knuckle-Boom crane control to reduce pendulum motion using the moving frame method,” in ASME International Mechanical Engineering Congress and Exposition. Available at: https://asmedigitalcollection.asme.org/IMECE/proceedings-pdf/IMECE2019/59414/V004T05A069/6512768/v004t05a069-imece2019-10436.pdf.

CrossRef Full Text | Google Scholar

El Beshbichi, O., Xing, Y., and Ong, M. C. (2021a). An object-oriented method for fully coupled analysis of floating offshore wind turbines through mapping of aerodynamic coefficients. Mar. Struct. 78, 102979. doi:10.1016/j.marstruc.2021.102979

CrossRef Full Text | Google Scholar

El Beshbichi, O., Xing, Y., and Ong, M. C. (2021b). An object-oriented method for fully coupled analysis of floating offshore wind turbines through mapping of aerodynamic coefficients. Mar. Struct. 78, 102979. doi:10.1016/j.marstruc.2021.102979

CrossRef Full Text | Google Scholar

Faltinsen, O. M. (1990). “Sea loads on ships and offshore structures,” in Cambridge ocean technology series (Cambridge: Cambridge University Press).

Google Scholar

Flatlandsmo, J., Smith, T., Halvorsen, r. O., and Impelluso, T. J. (2019). Modeling stabilization of crane-induced ship motion with gyroscopic control using the moving frame method. J. Comput. Nonlinear Dyn. 14. doi:10.1115/1.4042323

CrossRef Full Text | Google Scholar

Frankel, T. (2011). The geometry of physics: An introduction. Cambridge: Cambridge University Press. Google-Books-ID: gXvHCiUlCgUC.

Google Scholar

Grindheim, O.-M. R. (2022). Dynamic analysis of a multi-body floating wind turbine. Shock Vib. 2022. doi:10.1155/2022/6883663

CrossRef Full Text | Google Scholar

Holm, D. (2008). Geometric mechanics, part ii: Rotating, translating and rolling. Singapore: World Scientific Publishing Company. doi:10.1142/P549

CrossRef Full Text | Google Scholar

Impelluso, T. (2017). The moving frame method in dynamics: reforming a curriculum and assessment. Int. J. Mech. Eng. Educ. 46, 158–191. doi:10.1177/0306419017730633

CrossRef Full Text | Google Scholar

Jonkman, J. (2010). “Definition of the floating system for phase IV of OC3,”. Report No:NREL/TP-500-47535. doi:10.2172/979456

CrossRef Full Text | Google Scholar

Jonkman, J. M. (2007). Dynamics modeling and loads analysis of an offshore floating wind turbine. Report No:NREL/TP-500-41958. doi:10.2172/921803

CrossRef Full Text | Google Scholar

Jonkman, J. M., Larsen, T. J., Hansen, A. D., Nygaard, T. A., Maus, K., Karimirad, M., et al. (2010). “Offshore code comparison collaboration within iea wind task 23: phase iv results regarding floating wind turbine modeling,”. Report No: NREL/TP-5000-48191. doi:10.13140/2.1.3576.5768

CrossRef Full Text | Google Scholar

Korsvik, H. B., Rognsvaag, E. S., Tomren, T. H., Nyland, J. F., and Impelluso, T. J. (2019). Dual gyroscope wave energy converter. Energy ASME Int. Mech. Eng. Congr. Expo. 6. doi:10.1115/IMECE2019-10266.V006T06A008

CrossRef Full Text | Google Scholar

Murakami, H. (2013). “A moving frame method for multi-body dynamics,” in Proceedings of the ASME 2013 International Mechanical Engineering Congress and Exposition, San Diego, California, November 15–21, 2013. In IMECE2013.

CrossRef Full Text | Google Scholar

Ogilvie, T. F. (1964). Recent progress toward the understanding and prediction of ship motions. Delft, Netherlands: The Delft University of Technology.

Google Scholar

Perez, T., and Fossen, T. (2007). “Tutorial on modelling and simulation of marine system dynamics,” in IFAC Conference on Control Applications in Marine Systems (CAMS), Trondheim, Norway, September 13-16, 2016.

Google Scholar

Robertson, A., Jonkman, J., Vorpahl, F., Popko, W., Qvist, J., Frøyd, L., et al. (2014). Offshore code comparison collaboration ContinuationWithin IEAWind task 30: phase ii results regarding a floating semisubmersible wind system. Ocean Renew. Energy Int. Conf. Offshore Mech. Arct. Eng. 9. doi:10.1115/OMAE2014-24040.V09BT09A012

CrossRef Full Text | Google Scholar

Rykkje, T. R., Leinebø, D., Bergaas, E. S., Skjelde, A., and Impelluso, T. J. (2021). Modeling friction in robotic systems using the moving frame method in dynamics. Int. J. Mech. Eng. Educ. 49, 25–59. doi:10.1177/0306419019844160

CrossRef Full Text | Google Scholar

Storaas, T., Virkesdal, K., Brekke, G., Rykkje, T., and Impelluso, T. (2019). “Stabilizing ship motion with a dual system inertial disk,” in ASME 2019 International Mechanical Engineering Congress and Exposition, Salt Lake City, Utah, November 11–14, 2019.

CrossRef Full Text | Google Scholar

Wittenburg, J. (2007). Dynamics of multibody systems. Cham: Springer Berlin Heidelberg.

Google Scholar

WWF (2022). Vindkraft og naturen.

Google Scholar

Keywords: dynamics, moving, frame, method, digital, twin, multi-body

Citation: Grindheim O-M, Xing Y and Impelluso T (2023) Dynamic analysis and validation of a multi-body floating wind turbine using the moving frame method. Front. Mech. Eng 9:1156721. doi: 10.3389/fmech.2023.1156721

Received: 01 February 2023; Accepted: 11 August 2023;
Published: 28 August 2023.

Edited by:

Samir Emam, American University of Sharjah, United Arab Emirates

Reviewed by:

Paweł Kwiatoń, Częstochowa University of Technology, Poland
Mohammed Zikry, North Carolina State University, United States

Copyright © 2023 Grindheim, Xing and Impelluso. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Thomas Impelluso, tjm@hvl.no

Download