A New Proposal for the Closed-Loop Orientation Control of a Windﬂoat Turbine System

: It is believed that offshore wind power will play an important role in fulﬁlling the environmental commitments made in the Paris Agreement signed in 2015. The use of ﬂoating solutions in offshore wind power is gaining prominence, since it allows the expansion of these technologies to a large number of regions in which this has not been technically feasible until now. One of the current issues that needs to be studied in detail as regards offshore wind structures are the oscillations produced by the waves and the wind that could cause important problems in the turbine such as the reduction of the efﬁciency and performance or the appearance of structural damages owing to undesired loads on the blades, among others. This signiﬁes that the inclination of the structure must not exceed a certain small value in the operating phase. This paper, therefore, presents a simple dynamic modeling and a nonlinear multivariable control model-based system for a Windﬂoat turbine whose main function is to stabilize the system and to keep it vertical in order to improve its operation and, therefore, maintain the offshore wind turbine in the optimal conditions required to harness energy from the wind. The effectiveness of both the dynamic model and the proposed control algorithm has been demonstrated by developing various simulations in a MATLAB-Orcaﬂex environment, in which the Windﬂoat structure proved to have very small vertical deviations when wind impinges over the rotor, despite the undesirable effects caused by parametric uncertainties (e.g., the effect of the anchoring cables) in the system and unmodeled unmatched perturbations (such as the effects of the wind and the waves on the Windﬂoat turbine).


Introduction
There is currently a great interest in solving environmental problems with the aim of reducing the effects of climate change. The main objective of this is to eliminate the dependence on fossil fuels and the environmental repercussions that their exploitation has caused [1,2]. The reduction in greenhouse gas emissions, which are responsible for climate change, by means of renewable technologies for electricity generation is, therefore, being promoted in order to replace fossil fuels [3][4][5][6]. The European Union (EU) has established the need to achieve a reduction in greenhouse gas emissions of just above 80% with the objective of ensuring that more than 80% of the EU's electricity will be produced by renewable energy sources by 2050 [7,8]. Of the various renewable energy sources available, onshore wind power is now a mature technology and has been consolidated in numerous countries [9]. This renewable energy source has the possibility of being implemented on a large scale in electrical systems, signifying that an offshore installation is an option that will allow the installed wind power to be increased, since there is a larger surface area for the installation of these turbines [10,11].
(iii) an increase in the profitability of the project; (iv) less human intervention or (v) easy operation, among others.
The remainder of the paper is structured as follows: Section 2 presents brief descriptions of both the wind turbine and the proposed active floating structure studied herein. The conditions of hydrostatic equilibrium for the Windfloat turbine are presented in Section 3, after which the proposed dynamic model for the Windfloat turbine is developed in Section 4, resulting in a multivariable, time-dependent and strongly coupled system. In Section 5, this dynamic model is employed as the basis on which to propose the non-linear solution adopted to control the system. The validation setup performed by the research group using the Matlab-OrcaFlex tool is explained in Section 6, while Section 7 deals with the validation results obtained after using the proposed control system in both situations of model parametric uncertainties and unmodeled unmatched disturbances caused by the effects of the wind and the waves on the Windfloat turbine. Finally, Section 8 provides the conclusions of this paper and proposals for future work. Figure 1 shows the computer-aided design (CAD) general view of the platform studied in this work, which was created by keeping in mind the Windfloat concept [22] of a 10 MW DTU turbine [23][24][25]. This is because the platform design was conceived with safety margins and, therefore, had to be capable of supporting the loading forces that can be expected of wind turbines with power ratings of up to 10 MW. From the perspective of our proposal as regards controlling the orientation, the most important elements of the device provide the platform with its floating system, which is achieved by means of three cylindrical columns that contain an inner water ballast system whose main function is to stabilize the turbine and to keep it vertical in order to improve its operation. The water inside the ballast tanks is pumped internally to compensate the significant forces and torques that act on the turbine when the wind interacts with the blades. Each of the three water ballast compartments is located in the upper half of its associated column. The wind turbine itself is located on the top of one of the columns, which provides buoyancy to the whole weight of the turbine. The other two columns must be partially filled with solid ballast in order to compensate the difference in weights of the three columns and obtain the hydrostatic equilibrium that keeps the device floating and vertical in the absence of any kind of external disturbance. Heavy plates are located at the bottom of the columns that increase the weight and the viscous damping due to vortexes that occur at the edge of the plates. The added mass also increases due to a greater volume of surrounding fluid that must move, providing a wider hydrodynamic added inertia to the structure [26,27]. The platform is fixed to the seabed with at least three mooring lines (other works have proposed more lines, such as [22], which proposes six mooring lines). Table 1 summarizes the main dimensions of the elements that comprise the floating platform used in this work, while Table 2 contains the main dimensions of the wind turbine.  Table 2. Main dimensions of the wind turbine (taken from [25] Table 1 summarizes the main dimensions of the elements that comprise the floating platform used in this work, while Table 2 contains the main dimensions of the wind turbine.  In order to model the environmental loads applied to the generator, we chose a location near the coast of Ferrol, in Galicia (whose specific geographical coordinates are 43 • 50 N, 8 • 50 W) [28]. The bathymetry of the seabed and the winds in Galicia allow offshore wind energy exploitation. The sea depth in the neighborhood of the proposed location is approximately 150 m and consequently meets the necessary requirements for wind speed, which are 9.5 m/s at an elevation of 80 m [29].

Hydrostatic Equilibrium
Since the wind turbine is considered as a semisubmersible platform, the first step in the design and conception of the device must be to obtain a nominal hydrostatic equilibrium. In hydrostatic terminology, both the gravitational and the buoyancy forces are denominated as "restoring forces". This section proposes a simple distribution of the gravity forces resulting from the masses of the rigid bodies and the buoyancy forces owing to the ballast tanks located inside the device. This proposal is not intended to coincide with the current disposition of these elements in one of the real devices under study but is used for the sole purpose of this paper to show the proposed closed loop orientation control system. This study can, of course, be easily adapted to other mass arrangements and buoyancy tanks.
The most important assumption in this section is that: only three lumped masses of the same values are considered, which are denoted as m iRB , where subindex i denotes the number of the mass (and the number of cylinders in the device as i = 1, 2, 3) and subindex RB denotes rigid body.
Let us define a fixed reference frame S 0 as an orthogonal reference defined by the orthogonal vector set {x, y, z} whose z-axis is vertically oriented. Its origin is at the geometrical center of an equilateral triangle of side L that comprises the cylinders in the device at zero depth (z = 0), which is computed at the level of the sea (free surface). The horizontal plane (x, y) of S 0 coincide, therefore, coincides with a perfect, undisturbed flat free surface of the sea.
According to Figure 2, the proposed locations of the three lumped masses with respect to S 0 are: This location of masses leads to the local center of gravity (CoG) vector results of where h iR denotes the distance along the z axis for each i-th rigid body mass, in accordance with Figure 2: Let ∇ i0 denote the part of the volume (for i = 1, 2, 3) of fluid displaced by the device in a perfect flat sea, without any kind of mooring system and associated with the i-th mass (and cylinder or column), ρ be the water density and g be the acceleration of gravity, both considered constant. The resulting restoring forces for each mass are balanced and they are: r. Sci. Eng. 2021, 9, x FOR PEER REVIEW 6 of 28 The volume ∇ represents the submerged part of the i-th column (including its flatplate, structural elements as braces, etc.) corresponding to the i-th mass. If the volumes associated with each part of the whole device (including the not submerged parts) are denoted as , where > ∇ , and under the previously mentioned hydrostatic equilib- The volume ∇ i0 represents the submerged part of the i-th column (including its flatplate, structural elements as braces, etc.) corresponding to the i-th mass. If the volumes associated with each part of the whole device (including the not submerged parts) are denoted as V i , where V i > ∇ i0 , and under the previously mentioned hydrostatic equilibrium conditions m iRB , the wind generator floats and remains stable with null orientation.
Owing to the need to obtain a stable nominal hydrostatic equilibrium and good buoyancy conditions, a more extensive study has been carried out by using the MAXSURF Software [30], which makes it possible to obtain the main hydrostatic parameters. We employed the stability conditions and the displacement for the proposed platform based on the 3D CAD model with the main dimensions shown in Table 1, above.
The main hydrostatic parameters were obtained without any kind of mooring lines, in the absence of any kind of environmental disturbances and under a static flat and horizontal sea, and are summarized in Table 3, below. The notation employed corresponds with that usually used in naval architecture. The metacenters, KMs (denominated as KM_T when the platform rotates around the y-axis and KM_L when it rotates around the x-axis) are the vertical lines through the new center of buoyancy for the heel platform with each of the vertical lines through the original center of buoyancy to the lowest part of the platform, denoted as KB.
The platform is considered stable because the metacenter is greater than the height of the center of gravity at the lowest part of the platform, denoted as KG. The resulting torque will give the floated heeled body a rotation in the opposite direction to the tilt, thus allowing the body to return to its initial position.
There is no regulation that is directly applicable to the case of offshore platforms for wind turbines, and we have, therefore, employed the offshore oil and gas standards of Det Norske Veritas DNV [31]. According to this rule, when measuring the intact stability of the platform, the criteria are based on comparing the righting torque curve, at constant displacement, with the heeling torque curve caused by the wind, both as a function of the heel angle. The criteria of this standard are summarized in the following points: i The area under the righting torque curve (ART) at the second intercept must not be less than 30% in excess of the area under the wind heeling torque (AHT) curve at the same limiting angle; ii The first intercept must be less than 15 • ; iii The righting torque curve must be positive for the entire range of angles, from upright to the second intercept.
The stability has been verified only for operating conditions and according to DNV, and the heeling torque is owing to winds with a speed of 50 knots. These curves, which were obtained with the MAXSURF Software, are interesting for the transverse axis of the vessel, but in the case of the platform, the curves for both axes contained in the sea surface plane (x, y), were studied, because they had two degrees of freedom (DoFs), to be controlled. The corresponding graph for the y-axis is shown in Figure 3a, while that for the x-axis is shown in Figure 4a. These figures make it possible to appreciate that criteria ii. and iii. were fulfilled. Point i. of the criteria was verified with the data obtained with the MAXSURF Software, which are summarized in Table 4 (Areas in m 2 ). and the heeling torque is owing to winds with a speed of 50 knots. These curves, which were obtained with the MAXSURF Software, are interesting for the transverse axis of the vessel, but in the case of the platform, the curves for both axes contained in the sea surface plane ( , ), were studied, because they had two degrees of freedom (DoFs), to be controlled. The corresponding graph for the y-axis is shown in Figure 3a, while that for the xaxis is shown in Figure 4a. These figures make it possible to appreciate that criteria ii. and iii. were fulfilled. Point i. of the criteria was verified with the data obtained with the MAXSURF Software, which are summarized in Table 4 (Areas in m 2 ).     Mooring Lines Effect over the Device Figure 5 (below) shows the xy-plane view of the proposed disposition of the mooring lines from a local perspective. They were attached from the bottom of each of the heaveplates to the seabed. A different number and diverse disposition of mooring lines can be used (see [22] a.e.). Only three mooring lines were considered here because this was the minimum number of catenaries required in order to fix the device to the seabed.   Figure 5 (below) shows the xy-plane view of the proposed disposition of the mooring lines from a local perspective. They were attached from the bottom of each of the heaveplates to the seabed. A different number and diverse disposition of mooring lines can be used (see [22] a.e.). Only three mooring lines were considered here because this was the minimum number of catenaries required in order to fix the device to the seabed.  The catenaries are located symmetrically, and the angles of these three catenaries with respect to the z-axis are φ 1 = −30 • , φ 2 = 90 • and φ 3 = −150 • , respectively.

Mooring Lines Effect over the Device
In situations of null orientation, all catenaries will be identically located and the cable tension T of each catenary will produce a force over the base of the column (including its heave plate), which is (notations cα = cos(α) and sα = sin(α) are used for the sake of simplicity): where ξ denotes the angle that forms the i-th catenary with respect to the horizontal plane. The effect of the mooring lines on the device is the following force: It is easily observed that in the absence of any kind of marine current, there will not be any displacement in the x-y plane, but that the mooring lines from the perspective of the hydrostatic analysis will produce only a vertical displacement. By adding this effect to Equation (3), the new hydrostatic equilibrium point, therefore, results in: The hydrostatic force in heave will consequently be different with respect to that provided by Equation (3) in the sense that the heave displacement that the catenaries produce on the device is newly restored by incrementing the displaced volume (∇ i f > ∇ i0 ). The whole difference is: If perfect cylindrical shapes were considered, the increase in the vertical displacement of the device owing to the catenaries could be easily obtained as ∆z = Tsξ , where R C is the extern radius of the columns. Due to the presence of complex shapes, such as structural reinforcements, pontoon or bracing, and without taking into account the effects of disturbances such as biofouling, among others, this vertical displacement was obtained from the exhaustive numerical analysis using MAXSURF. The final coordinates of the CoG with respect to S 0 obtained according to Equations (6) and (7) were provided by Equation (2), with heights of the masses h iR , which were obtained in this final hydrostatic equilibrium and are represented in Figure 2, above.

Dynamic Model
Obtaining a dynamic model of the previously shown device, which is provided with two degrees of freedom (DoFs), was an essential stage when designing the controller. The proposed model falls in the category of a "Control Design Model" [32], and is considered to be a very simple model, if not the simplest. It is based on the following assumptions (in addition to assumption (i) shown in Section 3): (i) The vertical motion of the generator is ignored, because the water ballast is pumped only among the ballast tanks. (ii) The displacements of the center of gravity are negligible when ballast handling moves water among tanks. (iii) The added masses of the submerged part are constant, well known, and considered to be lumped. (iv) Viscous friction is modeled as being constant and is considered to be fully uncoupled.
(v) Restoring effects are computed geometrically. (vi) Only small angle variations are considered, because the main issue of the proposed controller is to maintain the vertical orientation of the generator.
The most important simplification of the dynamic model was stated in the previous section and consists of computing the dynamics for only three lumped masses of identical value. Any number of finite masses with suitable values can, of course, be used.
A new local reference frame S defined by orthogonal vectors {u, v, w} is linked to the turbine and located at its own center of gravity (CoG). Let us define the two DoF (only the rotations around axes x and y are computed) as ϕ x and θ y , respectively, which allow the frame S to rotate with regard to the fixed reference frame S 0 . These variables are grouped as a generalized coordinate vector q = ϕ x θ y T . This reference frame S is, therefore, displaced by only a vertical distance z 0 (|z 0 | > |z G0 |) with respect to S 0 when the turbine is arranged with a null orientation (q = 0 0 T = 0), as shown in Figure 6.
Once the fixed ballast (based on concrete, steel or other) and the fixed buoyancy tanks have both been adjusted in order to obtain good stability performance and null orientation, as shown in Section 3, above, the active hydrostatic actuators are considered (a proposal for active hydrostatic actuators was provided in [33]). These are responsible for creating the restoring forces that move the device, and they are ideally located in the same positions in which each of the three active ballast tanks of water associated with its cylinder are considered to be placed. These actuators will produce only vertical component forces that are obtained as differences between buoyancy forces and forces owing to gravity caused by water handling. For the sake of simplicity, forces from the actuators are considered as only buoyancy forces, which is in accordance with simplification (i), above, which defined masses as constants.
positions in which each of the three active ballast tanks of water associated with its cylinder are considered to be placed. These actuators will produce only vertical component forces that are obtained as differences between buoyancy forces and forces owing to gravity caused by water handling. For the sake of simplicity, forces from the actuators are considered as only buoyancy forces, which is in accordance with simplification (i), above, which defined masses as constants. From an inertial point of view, each of these masses is formed of both the rigid body mass and the added mass _ , which is also considered constant: = + _ (8) Figure 6. Reference frames and main magnitudes.
From an inertial point of view, each of these masses is formed of both the rigid body mass m iRB and the added mass m ADD_i , which is also considered constant: The net force of each of the elements (mainly cylinders) is obtained by computing the following forces: (a) that produced by the actuators (which is the multivariable control variable); (b) the loss of buoyancy, which is computed as the incremental volume that emerges outside the sea because of the motion of the device and (c) the incremental tension of the catenaries. The i-th net force produced by the actuators is, therefore: where ∆V u_i denotes the incremental volume of the i-th hydrostatic actuator and is a time dependent variable: ∆V u_i (t).

Control Torques Produced by Actuators
The forces produced by each of the actuators are applied to the device. A detailed proposal of an actuator system was provided in [33]. The application points of these forces with respect to the local reference frame S are: where h denotes the vertical component of the ballast tanks with respect to the rotating reference frame S. When rotations around fixed axes x and y are produced, the real positions of each point with respect to S 0 are obtained by using the composition of two basic rotation matrices [34], as follows: In order to compute the generalized torques from the hydrostatic actuators, and owing to the fact that the transfer of ballast water occurs internally, the sum of forces is null, that is: f u1 + f u2 + f u3 = 0. The torque vector produced by the actuators is, therefore: Due to the location of the turbine over the 2nd column, as represented in Figure 2, it is considered that f u2 = − f u1 − f u3 , which can be expressed as a function of only forces f u1 and f u3 by employing the following equation: The relationship shown below is then obtained: In Equation (15) it can be observed that the z-component of the torque is null, which justifies the fact that it is not possible to cause rotations along the z axis with only vertical hydrostatic forces in the proposed dynamic model (this result was previously represented in Figure 6 as φ z = 0).
If the torques around the x and y axes are grouped as τ = Γ x Γ y T , and forces f u1 and f u3 are grouped as F = f u1 f u3 T , the final relationship between the control forces and torques results in what it is denominated as the coupling matrix: where matrix D(q), results in a q-dependant well-defined term, and it relates the nondiagonal relationship between the forces and torques involved, which results in a strongly coupled, non-linear and time-dependent system.

Torques Owing to the Loss of Buoyancy
The submerged volume of each of the cylinders of the wind turbine does not remain constant if the orientation does not remain null q = 0. Figure 7 (below) shows the procedure used to obtain the torques owing to the loss of buoyancy. Only a brief schema of the columns of the wind turbine as cylinders with the heave-plates have been represented in this figure, and a rotation of an intentionally exaggerated angle of θ y = 10 • has been drawn.
The submerged volume of each of the cylinders of the wind turbine does not remain constant if the orientation does not remain null ≠ . Figure 7 (below) shows the procedure used to obtain the torques owing to the loss of buoyancy. Only a brief schema of the columns of the wind turbine as cylinders with the heave-plates have been represented in this figure, and a rotation of an intentionally exaggerated angle of = 10° has been drawn. The intersections of each of the axes of their associated cylinders with the sea surface when the device remains vertical are points with notation similar to the points defined in Equation (1) or (10), but with null vertical components, = 0. They are denoted as (the subindexes have been not represented in Figure 7 for the sake of clarity). When the orientation of the device is not null, the straight lines that define the axes of the cylinders can be expressed in a parametrical form as follows: where the points = ( ) • denote the rotated points from , is the parameter of each straight line and denotes the orientation of the rotated -axis from the reference frame that is obtained from the third column of ( ), provided in Equation (12). For each cylinder, the parameters are calculated as the intersection of the line defined by Equation (17) with the free surface ( = 0) (see Figure 4, in which only one dl parameter has been represented in cylinder 3). The intersections of each of the axes of their associated cylinders with the sea surface when the device remains vertical are points with notation similar to the points defined in Equation (1) or (10), but with null vertical components, z = 0. They are denoted as P 0i (the subindexes have been not represented in Figure 7 for the sake of clarity). When the orientation of the device is not null, the straight lines that define the axes of the cylinders can be expressed in a parametrical form as follows: where the points P 0Ri = R(q)·P 0i denote the rotated points from P 0i , d i is the parameter of each straight line and w denotes the orientation of the rotated z-axis from the reference frame S that is obtained from the third column of R(q), provided in Equation (12). For each cylinder, the parameters dl i are calculated as the intersection of the line defined by Equation (17) with the free surface (z = 0) (see Figure 4, in which only one dl parameter has been represented in cylinder 3).
After obtaining the three parameters dl i from the previous equation, the equivalent volumes of the loss (or gain) of buoyancy and the associated vertical forces are: where R CI is the inner radius of the water ballast tanks. The points of application of these vertical forces are: and the torque owing to this loss of buoyancy is eventually obtained from the first two components (components x and y, only, because component z is null) of the following equation:

Torques from Wind
Wind speed changes continuously, so it is necessary to describe it statistically. Several functions can be used to describe the frequency of the wind speed distribution, such as the Weibull distribution [35], W(β, η), where β represents the shape parameter (also known as the Weibull slope) and η denotes the scale parameter. The Weibull distribution of the wind speed is characteristic at each location. For a given wind velocity, the value of the forces on the hub of the wind turbine can then be estimated by means of a complex analysis, being of the form: where coefficient C T is a function of the tip speed ratio (TSR) and the shape of the blades, A is the area covered by the turbine blades, ρ Air is the air density and n denotes the speed of the wind. For any given rotor, the coefficient C T is well known and is tabulated as a function of the TSR. Although this force must be considered as a distributed set of forces over the rotor, because the hub and the low-speed shaft of the generator concentrate all the rotor efforts. These forces can be considered concentrated and their absolute value is denoted as F w . In a matrix notation, where α is the angle of incidence of the force over the hub of the gondola with respect to the horizontal plane x y , this results in: The generalized torque from the wind effects on the hub is easily computed: where distances d and H determine the position of the hub with respect to system S (these were intentionally omitted from Figure 2 for the sake of clarity).

Rotation Dynamics
In this subsection, the Lagrange formulation allows us to obtain the equations of motion for the wind turbine in the absence of any kind of external disturbances and frictions. The generalized coordinates are taken as previously: q = ϕ x θ y T . For the computation of the Lagrange function, the kinetic and potential energies are computed beforehand.

Local Inertia Matrix
The inertia matrix with respect to reference system S is obtained from the proposed lumped mass distribution. It can be written as a symmetric and positive definite matrix: J =   I xx P xy P xz P xy I yy P yz P xz P yz I zz   The values of the coefficients of the inertia matrix are shown in Appendix A. Under the hypothesis of perfect symmetry, m 1 = m 2 = m 3 = m, and m 1 + m 2 + m 3 = 3m, matrix J becomes simpler:

Computation of the Rotation Kinetic Energy
The kinetic energy of the wind turbine when rotations on first the x-axis and then the y-axis are applied is given by: where Ω x and Ω y denote the angular velocities over axes x and y, respectively. By substituting Equation (26) with P xy = P xz = 0 and Equation (12) in Equation (28), we obtain:

Computation of the Rotation Potential Energy
The above-mentioned absence of any kind of disturbance, and the hydrostatic equilibrium depicted in Section 3, signifies that the total potential energy has a zero value under these conditions.

Lagrange Formulation
Upon employing Lagrange's equations, expressed by: where L is the Lagrange function L = K − U = K, j = x, y and τ EXT is the sum of external torques, the final model is obtained (see the coefficients of matrices M(q) and C q, . q in Appendix A.
where the torque owing to the friction terms is approximated by employing two uncoupled quadratic functions with a constant friction coefficient according to assumption (v) and is a . q dependent function.
The remaining control variables on the right-hand side of Equation (31) have already been defined: the control torques τ = D(q)·F from the hydrostatic actuators were provided in Equation (16), the torques owing to the loss of buoyancy τ Buoy were defined as a function of q in Equation (22), the torque owing to the thrust horizontal force over the rotor τ w was defined as a function of q in Equation (25) and, finally, τ wwm denotes the unmodeled torques of the proposed model and the effects of the wind on the rest of the structure, from the waves and from the mooring lines, which are also considered as unmodeled external disturbances.

Proposed Dynamic Model
By grouping all of the torques that have been modeled, with the exception of the control torque, a function V q, . q is defined: Equation (32) can, therefore, be rewritten in a simplified form as:

M(q)
.. Figure 8 illustrates the proposed control system for the two DoF wind float studied herein. Let us define a desired reference for both axes x and y in the form of time functions: ϕ * x (t) and θ * y (t), respectively. A nonlinear feedback controller is used to accomplish the tracking of the desired controlled variables, which is intended to maintain their null value,

Proposed Control
= 0, and which is given by (time dependence is omitted): In which the auxiliary control signal Γ u is synthesized as follows: K P , K I and K D are diagonal positive definite matrices that denote the parameters of a vector-value classical uncoupled multivariable proportional-integral-derivative (PID) controller.
unmodeled torques of the proposed model and the effects of the wind on the rest of the structure, from the waves and from the mooring lines, which are also considered as unmodeled external disturbances.

Proposed Dynamic Model
By grouping all of the torques that have been modeled, with the exception of the control torque, a function ( , ) is defined: Equation (32) can, therefore, be rewritten in a simplified form as: Figure 8 illustrates the proposed control system for the two DoF wind float studied herein. Let us define a desired reference for both axes x and y in the form of time functions: * (t) and * (t), respectively. A nonlinear feedback controller is used to accomplish the tracking of the desired controlled variables, which is intended to maintain their null value, * (t) = ( * (t) * (t)) = , and which is given by (time dependence is omitted):

Proposed Control
In which the auxiliary control signal is synthesized as follows: , and are diagonal positive definite matrices that denote the parameters of a vector-value classical uncoupled multivariable proportional-integral-derivative (PID) controller. As will be noted in Figure 8, the control scheme is composed of an inner non-linear loop based on the proposed dynamic model and the outer loop based on classical PID controllers, which are responsible for obtaining the desired closed loop dynamics of the whole controlled system. The main objective of the inner loop is to obtain a linear and As will be noted in Figure 8, the control scheme is composed of an inner non-linear loop based on the proposed dynamic model and the outer loop based on classical PID controllers, which are responsible for obtaining the desired closed loop dynamics of the whole controlled system. The main objective of the inner loop is to obtain a linear and decoupled input/output relationship in the absence of any kind of external disturbance, which were included as τ wwm .
The sensorial system required in order to compute these non-linear terms has to be located in the device, and this is summarized as follows: (i) An inertial measurement unit (IMU) is required in order to obtain the orientation variables q(t) = ϕ x (t) θ y (t) T and their time derivatives . θ y (t) ) T . All the terms of functions M(q) and V q, . q can be computed from these signals, with the exception of the term τ w (q) and; (ii) The strain gauge load cells located over the bearings of the rotor make it possible to compute τ w (q) by measuring the real thrust over the hub of the rotor and applying Equation (25) or; (iii) By using the information captured by the vane and anemometer located on the gondola and estimating the force of the wind from tables of the known C T coefficients [24] and Equations (23)- (25).
The inner loop also requires the online computation of the inverse of the matrix that relates the control forces also provided by actuators F and control torques τ, D −1 (q) too. This matrix D(q) is invertible, with the exception of values of cϕ x cθ y = 0, which are unattainable. Furthermore, the main objective of the outer loop is to stabilize the overall system. Once the matrix D −1 (q) has been coupled to the dynamics of the system, this part of the controller design is highly simplified, since it operates an uncoupled 2 × 2 multivariable time invariant system.
The dynamics of the closed loop tracking error vector e q = (q * − q) = ϕ x θ y T is obtained after substituting expressions (36) and (37) in (35). In the absence of any kind of external disturbance, the following expression yields: ..
The closed loop error, therefore, evolves governed by the following third order 2 × 2 multivariable differential equation: whose time evolution of the error converges in an asymptotically-dominated manner onto a small as desired neighborhood of the origin of the phase space of the tracking error vector, depending on the design matrices of the outer loop of the controller {K D , K P , K I } designed. These controller matrices must be designed so as to render the following 2 × 2 complex valued diagonal matrix M(s), defined as: The desired closed loop dynamics is then proposed from [36,37]: M des (s) = I 2×2 s 3 + 5as 2 + 8a 2 s + 4a 3 (41) which can be defined as a function of only one design parameter a (in rad/s). A linear system with a monovariable characteristic polynomial provided by Equation (41) has a similar time response to the well-known critically dumped second order system with natural frequency ω n = a and a damping coefficient ξ = 1. The controller design matrices {K D , K P , K I } can then be directly obtained by equaling Equations (40) and (41).

Validation Setup
Given the obvious difficulty of building scale prototypes for the experimental validation of the proposed control algorithms, the use of computer simulation systems is required. The choice of OrcaFlex [38] as modeling and simulation software is justified, as it is a certified code and is widely used in the offshore and oil and gas industry, and provides a good treatment for the implementation of various funding systems. It also has extensive parameter entry possibilities, along the presentation of hydrodynamic results, and numerical and virtual reality graphics.

Orcaflex Model
The dynamic model required by OrcaFlex results in a very complex and complete model, in which the previous calculation of a large number of parameters is required. These can be grouped into the following: (i) environment parameters: state of the sea, waves, currents, depth, type of seabed, etc.; (ii) device shape parameters: dimensions, materials, frictions and added masses, among others and (iii) parameters of the anchoring system, which are mainly the types of cables and types of systems used to anchor devices to the seabed. Keeping in mind the geometric shapes of the elements that make up the Windfloat (see Figure 9), in order to model in OrcaFlex (see Figure 10), both the columns and heave-plates used were cylindrical. As can be seen in Figure 6, the position of the center of gravity of the water ballast and the moment of inertia of the hydrostatic actuator based on [31] have been modeled by taking into account its filling volume and its inclination. The Windfloat behavior model includes the remaining hydrostatic forces: the thrust produced by the hull volume, the external forces dependent on the state of the sea over the whole structure, the hydrodynamic lift of the hull and forces produced by the anchoring system, which can be modeled even with non-constant buoyancy, the variable cross-sectional area, etc. The model includes the effects of waves and the effects of wind on the whole structure, and not only on the rotor, as it was modeled for the purpose of this paper in Section 4.
parameter entry possibilities, along the presentation of hydrodynamic results, and numerical and virtual reality graphics.

Orcaflex Model
The dynamic model required by OrcaFlex results in a very complex and complete model, in which the previous calculation of a large number of parameters is required. These can be grouped into the following: (i) environment parameters: state of the sea, waves, currents, depth, type of seabed, etc.; (ii) device shape parameters: dimensions, materials, frictions and added masses, among others and (iii) parameters of the anchoring system, which are mainly the types of cables and types of systems used to anchor devices to the seabed. Keeping in mind the geometric shapes of the elements that make up the Windfloat (see Figure 9), in order to model in OrcaFlex (see Figure 10), both the columns and heave-plates used were cylindrical. As can be seen in Figure 6, the position of the center of gravity of the water ballast and the moment of inertia of the hydrostatic actuator based on [31] have been modeled by taking into account its filling volume and its inclination. The Windfloat behavior model includes the remaining hydrostatic forces: the thrust produced by the hull volume, the external forces dependent on the state of the sea over the whole structure, the hydrodynamic lift of the hull and forces produced by the anchoring system, which can be modeled even with non-constant buoyancy, the variable crosssectional area, etc. The model includes the effects of waves and the effects of wind on the whole structure, and not only on the rotor, as it was modeled for the purpose of this paper in Section 4.

Matlab-Orcaflex Integration
Given that OrcaFlex does not allow the computation of feedback variables, the Matlab environment and its interface with the OrcaFlex Application Programming Interface (API) was chosen for the implementation of the control algorithms proposed in Sec-

Matlab-Orcaflex Integration
Given that OrcaFlex does not allow the computation of feedback variables, the Matlab environment and its interface with the OrcaFlex Application Programming Interface (API) was chosen for the implementation of the control algorithms proposed in Section 4. The result of this integration was a joint tool that made it possible to validate in simulation the control systems designed, implemented in Matlab [39], on the dynamic model of the floating device previously developed in OrcaFlex. This integration was carried out by means of a Dynamic Link Library (DLL), which acted as an interface between both software packages. The integration scheme used is shown in Figure 11.

Matlab-Orcaflex Integration
Given that OrcaFlex does not allow the computation of feedback variables, the Matlab environment and its interface with the OrcaFlex Application Programming Interface (API) was chosen for the implementation of the control algorithms proposed in Section 4. The result of this integration was a joint tool that made it possible to validate in simulation the control systems designed, implemented in Matlab [39], on the dynamic model of the floating device previously developed in OrcaFlex. This integration was carried out by means of a Dynamic Link Library (DLL), which acted as an interface between both software packages. The integration scheme used is shown in Figure 11.    Figure 11 also shows the data signals flow of the integrated system: OrcaFlex simulates the dynamic behavior of the device using the parameters entered by the user as regards environmental forces, while the RunTime of Matlab computes and communicates the forces generated by the control system according to the proposed control system in Equations (36) and (37) and Figure 8.
The OrcaFlex time manager is responsible for generating the synchronism and the sampling periods for the full simulation. In this case, fixed intervals of sampling periods T s = 10 ms have been programmed, and each instant of double integration is, therefore, provided by t k = k T s , k ∈ N. Both the forces computed by the controller (from Matlab to OrcaFlex) as the variables of the position and orientation of the device (OrcaFlex to Matlab) are communicated through the use of the aforementioned DLL. OrcaFlex makes it possible to obtain temporary responses from the device, producing files of the simulations, videos of them and temporal responses of the variables for their subsequent analysis and treatment. The tool designed allows OrcaFlex to generate a large amount of variables from the simulation, such as the coordinates of centers of gravity, speeds, accelerations, rotation matrix, etc., even in the online mode. Matlab, meanwhile, takes the parameters of the controller and the reference from temporary files, and executes the controller algorithm, thus generating the aforementioned force vector produced by the dynamics of the hydrostatic actuators modeled and applied to the OrcaFlex model. It also allows the storage of the various variables and their subsequent graphic and numerical treatment.

Validation Results
This section shows the results of various simulations carried out to validate the correct behavior of the proposed control system with which to maintain the orientation of the previously proposed platform provided to a ballast system close to the null value.
The dynamic model implemented in OrcaFlex and the control system (37) described above were used to carry out numerical simulations with which to demonstrate the goodness of both in terms of controllability and the ability to achieve the desired orientation with the forces derives from the change in ballast.

Validation of System Decoupling
The set of simulations is formed for a step reference signal in the absence of external disturbances, such as the wind and waves, in order to validate the transient and steady state responses and test the dynamic decoupling between independent movements. Actuator saturation is not considered here, and a first order filter was used to obtain the time derivative of the output signals. The simulations in Figure 12a show the temporal responses obtained for the variable ϕ x (t) before a reference step signal, ϕ x * (t), with an increase of ±1 • and transitions in 250 s, 550 s and 850 s. The regulator gains were obtained from (41) with the design parameter a, which takes the values 0.05 rad/s, 0.1 rad/s and 0.2 rad/s. Figure 12b shows the reference signal θ y * (t), corresponding to a step signal with the same increases and transitions at 100 s, 400 s and 700 s, and the temporal responses for the variable θ y (t), for the same values a. The initial operating point is taken as being  Step responses for the y-axis rotation. Figure 13a shows the zoom of the temporal responses ϕ x (t) in the time interval [350, 450] scan. The response was better while the closed loop dynamics was faster. The steady error was not fully canceled in the neighborhood of t 2θ = 400 s , which was owing to the maximum change in the orientation θ y * (t 2θ = 400), and there was no perfect decoupling at all. Figure 13b similarly shows the zoom of the temporal orientation responses θ y (t) for the time interval [500, 600] s, in which the maximum change in the orientation ϕ x * (t) occurred at t 2ϕ = 550 s . Step responses for the -axis rotation Figure 13a shows the zoom of the temporal responses ( ) in the time interval 350, 450 scan. The response was better while the closed loop dynamics was faster. The steady error was not fully canceled in the neighborhood of = 400 , which was owing to the maximum change in the orientation * ( = 400), and there was no perfect decoupling at all. Figure 13b similarly shows the zoom of the temporal orientation responses ( ) for the time interval 500, 600 s, in which the maximum change in the orientation * ( ) occurred at = 550 .

External Disturbances: Waves and Wind
Various environmental models can be used to achieve realistic simulations of any kind of marine device. In this work, a classical wave spectrum ( , ) was used as a function of both the frequency and the direction , with a Pierson-Moskowitz (PM) spectrum [40]. Due to the closed loop dynamics expected for the device, owing to the presence of the superposed low frequency (LF) model and wave frequency (WF), only the LF effects will be considered as external disturbances. Figure 14 shows the PM spectrum for a significant wave height of 1.5, 2 and 2.5 for the peak periods of 7, 8, 9, 10 and 11 each (the values selected are the most frequent in the chosen location [29]), along with the Bode (a linear rather than a logarithmic plot for easy comparison) magnitude for the proposed closed loop desired dynamics with a different design parameter , whose values were 0.05, 0.075, 0.1, 0.15 and 0.2 / .
Upon observing the graph, the interaction with the waves of the desired closed loop dynamics with = 0.15 / and = 0.2 / , which is greater than that of the rest, can be appreciated. The regulator chosen to carry out the following simulations will be that which has a value for the design parameter of 0.1 / , because the dynamics corresponding to the regulators with = 0.05 / and = 0.075 / are considered too slow.

External Disturbances: Waves and Wind
Various environmental models can be used to achieve realistic simulations of any kind of marine device. In this work, a classical wave spectrum W(ω, ψ) was used as a function of both the frequency ω and the direction ψ, with a Pierson-Moskowitz (PM) spectrum [40]. Due to the closed loop dynamics expected for the device, owing to the presence of the superposed low frequency (LF) model and wave frequency (WF), only the LF effects will be considered as external disturbances. Figure 14 shows the PM spectrum for a significant wave height of 1.5, 2 and 2.5 m for the peak periods of 7, 8, 9, 10 and 11 s each (the values selected are the most frequent in the chosen location [29]), along with the Bode (a linear rather than a logarithmic plot for easy comparison) magnitude for the proposed closed loop desired dynamics with a different design parameter a, whose values were 0.05, 0.075, 0.1, 0.15 and 0.2 rad/s. Upon observing the graph, the interaction with the waves of the desired closed loop dynamics with a = 0.15 rad/s and a = 0.2 rad/s, which is greater than that of the rest, can be appreciated. The regulator chosen to carry out the following simulations will be that which has a value for the design parameter of 0.1 rad/s, because the dynamics corresponding to the regulators with a = 0.05 rad/s and a = 0.075 rad/s are considered too slow.
Upon observing the graph, the interaction with the waves of the desired closed loop dynamics with = 0.15 / and = 0.2 / , which is greater than that of the rest, can be appreciated. The regulator chosen to carry out the following simulations will be that which has a value for the design parameter of 0.1 / , because the dynamics corresponding to the regulators with = 0.05 / and = 0.075 / are considered too slow.  From the point of view of the frequency domain, the controlled platform exhibits a closed-loop dynamics given by Equation (41) and behaves as a low-pass filter with the cut-off frequency set in parameter "a", without any kind of resonance. Then, the wave frequency zone acts on the platform in the rejection zone of the low-pass filter as external disturbances that are all rejected with a minimum attenuation of approximately −130 dB for a value of a = 0.10 rad/s, which results in an imperceptible effect on the time domain responses.
Furthermore, in the first set of simulations in this section, the behavior of the platform to step inputs of both angular coordinates was shown. A simulation with null reference signals and step variations of the wind magnitudes is presented as follows. This simulation does not, of course, represent a realistic case, but this kind of simulation is very useful as regards demonstrating the goodness of the controller [41,42]. Both the orientation responses of the platform under open loop and closed loop dynamics for the controller with the previously chosen design parameter a = 0.1 rad/s and when it has been subjected to wind disturbances with the next step-type signals with transitions in both the speed, n(t) (m/s), and orientation, α(t) ( • ), are depicted in Figure 15.
The first graph shows, from the open loop time response, that the platform cannot maintain a null orientation in ϕ x (t), reaching a heel angle of −23 • after the first transition in the wind speed with an undamped oscillating response of around −14 • . When the second transition in wind speed occurred, the oscillating response of the platform increased its inclination to -18 • . Owing to this response, this orientation hardly varied when the wind decreased its speed or changes its orientation, and it was possible only to appreciate a decrease in the oscillation value. With regard to the orientation θ y (t) for the open loop response, which can be seen in the second graph, the change was perceptible only after the transition in the wind direction, increasing the inclination up to approximately 6 • .
The behavior of the platform after the same disturbances created by the wind according to Equation (42), but in a closed loop, is then plotted on the same graphs, which makes it possible to verify that the proposed control system is able to compensate for the wind effects, but that the non-modeled effects caused by the unmodeled terms owing to the mooring cables are still present. The effects of the sudden changes in wind speed and orientation on the behavior of the platform can also be appreciated in both Figure 15a,b because of the non perfect decoupling term based on the inverse of matrix D(q) defined in Equation (16).
The first graph shows, from the open loop time response, that the platform cannot maintain a null orientation in ( ), reaching a heel angle of −23° after the first transition in the wind speed with an undamped oscillating response of around−14°. When the second transition in wind speed occurred, the oscillating response of the platform increased its inclination to -18°. Owing to this response, this orientation hardly varied when the wind decreased its speed or changes its orientation, and it was possible only to appreciate a decrease in the oscillation value. With regard to the orientation ( ) for the open loop response, which can be seen in the second graph, the change was perceptible only after the transition in the wind direction, increasing the inclination up to approximately 6°. The behavior of the platform after the same disturbances created by the wind according to Equation (42), but in a closed loop, is then plotted on the same graphs, which makes it possible to verify that the proposed control system is able to compensate for the wind effects, but that the non-modeled effects caused by the unmodeled terms owing to the mooring cables are still present. The effects of the sudden changes in wind speed and orientation on the behavior of the platform can also be appreciated in both Figure 15a,b because of the non perfect decoupling term based on the inverse of matrix ( ) defined in Equation (16).
In the last simulation, the wind values were obtained from the data recorded by the La Langosteira II buoy (43.35N, 8.56W) [29], which is located very close to the point chosen for the location of the platform. The registered values of wind speed and orientation are, according to the Weibull distribution, which are (8.34,12.29) and (220.07,10.05) respectively, expected for this location. The compensation of the wind signal was obtained after being processed by a low pass filter, while Gaussian noise from a normal distribution (0,2 • 10 ) was intentionally added. The saturation and non-ideal dynamics of the actuators were programmed. We also modeled a regular wave with a period of 9 and a significant height of 1.5 , which is the most frequent wave according to the data collected by the buoy.  In the last simulation, the wind values were obtained from the data recorded by the La Langosteira II buoy (43. 35 N, 8.56 W) [29], which is located very close to the point chosen for the location of the platform. The registered values of wind speed and orientation are, according to the Weibull distribution, which are W n (8.34, 12.29) and W α (220.07, 10.05) respectively, expected for this location. The compensation of the wind signal was obtained after being processed by a low pass filter, while Gaussian noise from a normal distribution N 0, 2 × 10 −3 was intentionally added. The saturation and non-ideal dynamics of the actuators were programmed. We also modeled a regular wave with a period of 9 s and a significant height of 1.5 m , which is the most frequent wave according to the data collected by the buoy.
The platform was then subjected to the force of the wind, resulting in angle variations of up to approximately ϕ x (t) = 23 • and θ y (t) = 6 • when it operated in open loop. The time intervals in which its oscillations, owing to the effect of waves only, were close to the desired null value for both angles in the first and last section of the simulation interval, correspond to time intervals with an intentionally null speed of the wind but in the presence of wave disturbances.
With regard to the closed loop time responses, which have been plotted in an overlapping manner with the open loop time responses, it can be seen in Figure 16 that the proposed control system is able to maintain the platform with almost null orientation (less than 0.

Conclusions
The use of floating solutions is gaining prominence in offshore wind power owing to the fact that floating wind power multiplies the energy potential of this form of energy several times by obtaining better use of the wind and much higher capacity factors. An important question that needs to be investigated in offshore wind structures are the oscillations produced by the waves and the wind that could cause important reductions in efficiency and performance of the turbine or the appearance of structural damages owing to undesired loads on the blades, signifying that the inclination of the structure must not exceed 10° in the operating phase. This paper proposed a new dynamic model and its corresponding control algorithm for Windfloat turbines with the objective of substantially reducing/eliminating changes in the orientation of the floating structure and obtaining an optimal use of the wind and much higher capacity factors. The dynamic model proposed was based solely on the definition of three lumped masses that permit the forces from the three buoyancy actuators to be converted into three torques that are responsible for the rotations of the Windfloat structure. Furthermore, the definition of the nonlinear control system can be decomposed into the following three states: (i) the definition of an uncoupling matrix that is responsible for decoupling the Windfloat turbine motions for nonsingular orientation values; (ii) a compensation of the centrifugal and Coriolis torques by means of a nonlinear uncoupling-based model and (iii) the definition of an auxiliary proportional-integral-derivative (PID) control action.
Numerical simulations were provided, in which the effectiveness of the proposed control method is assessed with respect to parametric uncertainties (e.g., the effect of the

Conclusions
The use of floating solutions is gaining prominence in offshore wind power owing to the fact that floating wind power multiplies the energy potential of this form of energy several times by obtaining better use of the wind and much higher capacity factors. An important question that needs to be investigated in offshore wind structures are the oscillations produced by the waves and the wind that could cause important reductions in efficiency and performance of the turbine or the appearance of structural damages owing to undesired loads on the blades, signifying that the inclination of the structure must not exceed 10 • in the operating phase. This paper proposed a new dynamic model and its corresponding control algorithm for Windfloat turbines with the objective of substantially reducing/eliminating changes in the orientation of the floating structure and obtaining an optimal use of the wind and much higher capacity factors. The dynamic model proposed was based solely on the definition of three lumped masses that permit the forces from the three buoyancy actuators to be converted into three torques that are responsible for the rotations of the Windfloat structure. Furthermore, the definition of the nonlinear control system can be decomposed into the following three states: (i) the definition of an uncoupling matrix that is responsible for decoupling the Windfloat turbine motions for non-singular orientation values; (ii) a compensation of the centrifugal and Coriolis torques by means of a nonlinear uncoupling-based model and (iii) the definition of an auxiliary proportional-integral-derivative (PID) control action.
Numerical simulations were provided, in which the effectiveness of the proposed control method is assessed with respect to parametric uncertainties (e.g., the effect of the anchoring cables) in the system, and unmodeled unmatched perturbations. The results obtained illustrate that the Windfloat turbine performed extremely well with the desired settling time and with no overshoot and achieved almost null orientation as regards both angles of the platform, in addition to obtaining a small variation in angles of approximately 0.02 • of them. Finally, the study of different control strategies that improve the harnessing of the wind resource and higher capacity factors will be the topics of our future research.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Acknowledgments:
The authors wish to thank the staff of the Grupo de Investigación Tecnológico en Energías Renovables Marinas from the UPM -GIT-ERM Group-, and F.R.V. for programming the aforementioned DLL allowing the realization of this work.

Conflicts of Interest:
The authors declare no conflict of interest.

Appendix A
Matrices M(q) and C q, . q and coefficients of the inertia matrix This appendix shows how the matrices M(q) and C q, . q were obtained. The partial derivatives with regard to angular speeds, which were noted as Ω x and Ω y, , were obtained.