Variational Control Approach to Energy Extraction from a Fluid Flow

: Energy harvesting from the environment is an important aspect of many technologies. The scale of energy capturing and storage can involve the power range from mWatt up to MWatt, depending on the used devices and the considered environments (from ambient acoustic and vibration to ocean wave motion, or wind). In this paper, the wind turbine energy harvesting problem is approached as an optimal control problem, where the objective function is the absorption of an amount of energy in a given time interval by a ﬂuid-ﬂow environment, that should be maximized. The interest relies on outlining general control models of ﬂuid-ﬂow-based extraction plants and identifying an optimum strategy for the regulation of an electrical machine to obtain a maximum-e ﬃ ciency process for the related energy storage. The mathematical tools are found in the light of optimal control theory, where solutions to the fundamental equations are in the frame of Variational Control (the basis of the Pontryagin optimal control theory). A special problem, named Optimally Controlled Betz’s Machine OCBM-optimal control steady wind turbine, is solved in closed form, and it is shown that, in the simpler steady case, it reproduces the maximum e ﬃ ciency machine developed in Betz’s theory.


Introduction
The problem of energy extraction from an environmental source is certainly one of the key problems in the history of technology. This issue has set the basis of thermal engines and the related thermodynamic theory since the 17th century and recalls the notion of engine efficiency that reaches its highest point in the statement of the Carnot's theorem. In this frame, one recognizes that, from a heat source, which is a bath of vibrating particles in disordered motion, one can extract only a given fraction of that energy in a usable form. The efficiency of engines has been investigated for the last two centuries and still represents an important challenge in engineering. Engine efficiency does not only attract the interest of engineers, but is spread in other scientific fields. Indeed, recent is the work of the physicists Curzon and Ahlborn [1] about the notion of efficiency at maximum power. This concept moves the investigation from the point of view of the overall Carnot maximum efficiency, reached in reversible conditions, i.e., at zero power, towards the analysis of the efficiency at maximum power. This point of view meets the real interest of building up engines that produce the highest possible power in irreversible conditions, and [1] provides a very simple and elegant efficiency expression. The concept of efficiency at maximum power can be applied to any energy source. Indeed, in the theory of wind turbines, developed by Betz [2] in the 1930s, the maximum efficiency of a windmill that extracts energy from a steady wind flow has been investigated, unveiling the existence of a limit to the maximum extractable power, analogously with the thermodynamic efficiency approached by Carnot.

A Generalized Model for an Optimal Unsteady Wind Machine
The formulation of the problem leading to an optimally controlled unsteady wind machine is outlined.
First, Betz's theory is resumed to highlight its structure to develop a more general model that removes some of the hypotheses used in the original optimization process. The key point is the Betz's formulation does not determine the actual speed of the propeller, since the force balance of the propeller is not considered, and does not involve the resistance force and the inertia effects on the rotor. This implies that the analysis leaves unknown one of the physical key quantities of the problem, the propeller speed, since, as explained in the next section, the original Betz's theory uses more variables than equations. In fact, the optimality conditions express considering the generated power as a function of the speed ratio (the one between the speed of the inflow and the speed of the flow leaving the propeller). However, to state an optimal control problem for the wind machine, we need to remove the steady state hypothesis, and we need to introduce explicitly the force balance of the propeller, emphasizing here the role of the resistant torque due to the coupled electrical machine.
To make the approach clearer, two different schemes are used, the ones in Figures 1 and 2, respectively. Figure 1 shows the scheme of the physical device in which are represented the inflow section (area A), the propeller section (area A S ), the final section (area A'). The propeller is nothing Energies 2020, 13, 4913 3 of 20 more than a moving surface that is in contact with the fluid current through the blades and connected by its shaft to an electrical device that continuously absorbs the power generated by the fluid flow. Figure 2 is a more essential and schematic representation of the phenomenon investigated in Figure 1 in which, independently of the real configuration of the turbine, only the basics physical facts are represented. In this scheme, the air flow is impinging the moving mass m, that is the analogous of the inertia effect of the rotor, and the damper c is the dissipator that is a very schematic and elemental representation of the electrical machinery. The velocity of the mass is the equivalent of the angular speed of the turbine and the mass m an equivalent of the rotor moment of inertia.
Energies 2020, 13, x FOR PEER REVIEW 3 of 21 Figure 2 is a more essential and schematic representation of the phenomenon investigated in Figure  1 in which, independently of the real configuration of the turbine, only the basics physical facts are represented. In this scheme, the air flow is impinging the moving mass m, that is the analogous of the inertia effect of the rotor, and the damper c is the dissipator that is a very schematic and elemental representation of the electrical machinery. The velocity of the mass is the equivalent of the angular speed of the turbine and the mass m an equivalent of the rotor moment of inertia.  where V is the inflow velocity and the velocity of the propeller, using the propeller disk area and a pressure coefficient C, the force applied to the propeller is as: For the turbine propeller, the dynamic equation reads: where is the propeller inertia force and is the electrical machine resistant force (assumed simply proportional to the propeller speed), and c is the control parameter. In fact, the force has a twofold role. On one hand, it introduces a force opposed to the wind motor force, the presence of which, together with the inertia force of the propeller, produces an additional differential equation that closes the balance between the number of variables and the number of equations used in the   Figure 2 is a more essential and schematic representation of the phenomenon investigated in Figure  1 in which, independently of the real configuration of the turbine, only the basics physical facts are represented. In this scheme, the air flow is impinging the moving mass m, that is the analogous of the inertia effect of the rotor, and the damper c is the dissipator that is a very schematic and elemental representation of the electrical machinery. The velocity of the mass is the equivalent of the angular speed of the turbine and the mass m an equivalent of the rotor moment of inertia.  where V is the inflow velocity and the velocity of the propeller, using the propeller disk area and a pressure coefficient C, the force applied to the propeller is as: For the turbine propeller, the dynamic equation reads: where is the propeller inertia force and is the electrical machine resistant force (assumed simply proportional to the propeller speed), and c is the control parameter. In fact, the force has a twofold role. On one hand, it introduces a force opposed to the wind motor force, the presence of which, together with the inertia force of the propeller, produces an additional differential equation that closes the balance between the number of variables and the number of equations used in the For example, with reference to Figure 1, expressing the stagnation pressure as p = 1 2 ρ(V − V s ) 2 , where V is the inflow velocity and V s the velocity of the propeller, using the propeller disk area A s and a pressure coefficient C, the force applied to the propeller is as: For the turbine propeller, the dynamic equation reads: V s is the propeller inertia force and c . V s is the electrical machine resistant force (assumed simply proportional to the propeller speed), and c is the control parameter. In fact, the force c . V s has Energies 2020, 13, 4913 4 of 20 a twofold role. On one hand, it introduces a force opposed to the wind motor force, the presence of which, together with the inertia force of the propeller, produces an additional differential equation that closes the balance between the number of variables and the number of equations used in the model (see next section). On the other hand, the modulation of the resistance force, through the parameter c, due to the presence of the electrical machine, is responsible for the actual speed of the propeller. For this reason, this force directly controls the instant power the device is extracting from the inflow. The optimal control problem is formulated in terms of the regulation of this propeller resistance force in a way the total energy extracted by the wind unsteady machine is maximum, during a suitable interval of observation of the phenomenon.

Resume and Comparison with the Betz's Theory
A fluid flow of speed V(t) and density ρ transports the power 0.5ρV(t) 3 A through the section of area A, as in Figure 1. The problem posed by the classical Betz's theory relies on the identification of the conditions that maximize the power energy extraction from that flow.
It is instructive to resume the standard Betz's theory, since we are interested to put it in the simplest form to be treated as an optimal control theory problem.
The Betz's theory introduces three fundamental velocities: the inflow velocity V, passing through a section area A, the velocity V s of the propeller solid surface, passing through the propeller disk area A s , and the velocity V leaving the propeller disk and passing through a section area A The velocity V is assigned, the others two speeds are instead unknowns of the problem. Similarly, the section area A s is given, while A and A' are unknowns.
Several simplification hypotheses are introduced in the original Betz's theory, and namely: (i) the flow is incompressible (the flow speed is largely subsonic), (ii) the velocity of the flow is purely axial (no radial flow speed), (iii) absence of internal mechanisms of flow dissipation, (iv) vanishing hub-shaft size, (v) steady state flow conditions. Note that in the Betz's theory, the involved physical quantities are V, V s , V', A, A s , A' and the force F e the flow exerts on the propeller. These seven variables are not independent each other, since in the classical theory of the Betz's wind turbine, four fundamental relationships are used. Namely: (i) two mass-balance equations, VA = V s A s and V s A s = V A expressing the equivalence of the mass rates through the three considered sections; (ii) the momentum balance, to determine the force F e = dm dt (V − V ) = ρA s V s (V − V ) applied to the propeller disk; (iii) the power balance, expressing the equivalence between the shaft power F e V s = ρA s V 2 s (V − V ), and the kinetic energy rate of the flow when passing through the propeller disk, expressed by These are four equations, in terms of the five unknowns V s , V', A, A', F e , since V and A s are assigned a priori, and the quantities of interest can be determined, except one. In the present approach, we introduce the additional force balance of the propeller, closing the equations system for the model, imposing the equilibrium between the wind force, the resistance force of the electrical machine, and the inertia force of the rotor.
The analysis of the best performance of the wind machine is simply determined by the power expression. Under steady conditions, since the force, in the Betz's theory, is F e = ρA s V s (V − V ), and since V s = 1 2 (V + V ), then F e = 2ρA s V s (V − V s ) and the power finally is P = F e V s = 2ρA s V s 2 (V − V s ). Its maximum is obtained for dP dV s = 0, that simply produces V s = 2 3 V, that substituted back into the power expression provides P max = 16 27 1 2 ρA s V s 3 , that is about 60% of the power 1 2 ρA s V s 3 carried by the inflow. The principle for this steady optimization is essentially the same, whatever the approximation used for the propeller force. For example, in the present model using the static pressure, one has F e = 1 2 ρA s C(V − V s ) 2 , the power is P = F e V s = 1 2 ρA s CV s (V − V s ) 2 , and its maximum is for dP dV s = 0, that simply produces V s = 1 3 V and for the maximum power it follows P max = 4 27 C 1 2 ρA s V s 3 , that is a fraction 4 27 C of the available inflow power 1 2 ρA s V s 3 captured by the wind machine. To make the static pressure model of the force and the Betz's theory equivalent, under the profile of the maximum power, it is sufficient to set 4 27 C = 16 27 , i.e., C = 4. This analysis suggests that, independently of the model selected to represent the wind turbine, it exists an optimal working condition, in correspondence of which the produced power is maximum. This analysis, however, refers to steady state conditions and neglects the dynamic of the wind turbine rotor. Thus, one could seek to continuously optimize the wind power modifying the resistant torque exerted by the electrical machine under the action of an unsteady wind.

Unsteady Optimal Control Model
A fluid flow of speed V(t) and density ρ transports the power 0.5ρV(t) 3 A through section A, as in Figure 1. The power device, i.e., the wind turbine blade, is assumed to consist of an active moving surface of characteristic area S and inertia m (for any given shape, e.g., blade, piston, etc. and any motion, i.e., linear or rotational, m can represent either a mass or an inertia moment), as represented in the dummy model in Figure 2. The motion of the power device is partially induced by the action due to the pressure p(t) or the torque, which is expressed by p(t)SC, where C is the pressure coefficient, and by the action of the electrical generator, according to F e (t) = c(t)V s (t): applied to the moving surface by the electrical machine converting mechanical energy into electrical, this action is proportional to the velocity of surface V s (t) through a regulation/control parameter of the electrical machine c(t) function of time.
The mobile surface has the characteristic speed V s (t) (indicating either a linear or rotational speed). Note that the form F e (t) = c(t)V s (t) for the electrical force is reasonable as long as the moving element of the machine is an electrical circuit, in which the current interacts through the Lorentz force with a magnetic field, applied by the static part of the electrical device, the intensity of which can be controlled by c(t). Additionally, V s (t) is a quantity monitored by a sensor, an information that can be eventually used in the effort to regulate the electrical machine through c(t). OCBM-optimally controlled steady wind turbine problem consists in determining how the regulation law of the electrical machine c(t) should be designed to maximize the energy E extracted from the flow. Note that this approach can be considered a universal prototype model to be applied to different energy control applications.
The mobile surface dynamics is described by: where the sign(.) function is used to include the wind force direction. Assuming the flow pressure p against S expressed by: the following: represents the device dynamics. The instant power produced by the electrical machine is: The electrical energy E the device produces, in each time interval T, is: The energy E becomes a functional depending on two unknown functions: the control action c(t), and the mobile surface speed V s (t). The problem of maximizing the extracted energy E is, in this way, translated into the optimal definition of these two functions. The OCBM problem can be stated in a strict mathematical form as: This problem can be conveniently rewritten in a compact form introducing the restraint represented by the device dynamics, by using the time-dependent Lagrange multiplier λ(t):

The Euler-Lagrange Equations and Their Solution
As stated in the previous section, the problem has been formulated using the Pontryagin optimal control technique [30][31][32] and the functional E c, V s , . V s , λ depends on the three unknown functions c, V s , λ. Variation of the functional yields: The Euler-Lagrange equations associated to this problem define a system of three nonlinear equations: is the Delta function. Some considerations can be deduced before tackling the solution to this problem. The first transversality condition δV s (0) = 0 is always satisfied for any assigned initial conditions V s0 . Since no final-time condition can be prescribed on the speed (the differential equation in V s is of first order, and admits only one initial condition), it follows necessarily δV s (T) 0. Eventually, for the solution of the system of Equation (9) to be a stationary point for E, the solution itself must satisfy also the final condition λ(T) = 0.
To find the solution, investigated for V(t) > 0, the second part of Equation (9), λ(t) = −V s (t), is substituted in the other two, returning: A further simplification can be obtained by eliminating the control variable c from the previous equations, obtaining V s (t) = 1 3 V(t). Eventually, the control variable is recovered, and the final optimal solution of the Euler-Lagrange equations, i.e., the desired solution to the OCBM problem, marked by the star exponent, is determined: The control action can be directly piloted in terms of the measured quantity V s leading to the control system: Note that the condition For an arbitrarily prescribed choice of T, this condition in general fails. Let assume the wind turbine works in the time interval t start = 0, when the wind starts blowing, and t stop = T, large enough to ensure that the wind had stopped, producing V(T) = 0. Clearly, the control logic operates in the same interval and the transversality condition λ(T) = 0 is consequently satisfied.
The nonlinear variational control Equation (12) has a clear physical interpretation. In steady condition, the electrical resistant force F e = cV s = 2 9 ρCSV 2 depends on the wind speed: the higher the wind speed, the higher the resistant force, in a way the optimal speed of the moving surface follows the wind speed: V s (t) = 1 3 V(t). When the wind speed suddenly increases, the term − m .

V(t)
V(t) produces a smaller c(t), reducing the electrical resistance force, thus promoting a quick acceleration of the moving surface and a prompt increase of the power. When the wind speed stabilizes at its highest value, and . V(t) drops down, the electrical force progressively rises to a higher final value according to the wind speed increase.
Note this kind of control can be performed by acquiring the control signals V and its derivative from a wind speed meter. However, the electrical machine regulation can be controlled also by acquiring directly the speed signal from the wind turbine shaft, since, in an ideal case, the dynamics of the system produces V = 3V s .

Local Maximum and Directional Derivatives
At first, a Hessian analysis of (7) is carried out so to check whether the obtained optimal trajectories present an absolute maximum or minimum, a local maximum or minimum, or a saddle point. The solution is analyzed for V > V s > 0 and the cost function E can be expressed as function of V S , . V S , V, by substituting the control Equation (12): By calculating the second gradient in terms of the rotor speed and its acceleration into the Hessian matrix and substituting the optimal trajectories of Equation (11) The examination of the determinant demonstrates the obtained trajectory belongs to a saddle point, since DetH = −m 2 < 0. The Hessian analysis does not clarify whether the application control logic of Equation (12) enables the achievement of the maximum extracted power. Therefore, it is more helpful to analyse the optimal trend of the extracted power P e obtained by Equation (13), in terms of the rotor speed V s and the acceleration/deceleration . V s , which are normalized with respect to their design value: where V s max and . V s max are the maximum rotor speed and acceleration respectively and so doing Figure 3 shows how the power trend P e V s , . V s changes for a fixed wind speed V. In this case, V is set to the maximum normalized rotor speed, namely V = 1. It is apparent how for a generic data set related to the turbine parameters m, ρ, S, C, the extracted power is maximum only if the rotor seed reaches one third, as in Equation (11), condition that corresponds to the red dashed line, which identifies the maximum of P e V * s , . V * s . This is further confirmed by the derivative By calculating the second gradient in terms of the rotor speed and its acceleration into the Hessian matrix and substituting the optimal trajectories of Equation (11) still into the Hessian matrix, one obtains: * , * = −12 The examination of the determinant demonstrates the obtained trajectory belongs to a saddle point, since = − 0. The Hessian analysis does not clarify whether the application control logic of Equation (12) enables the achievement of the maximum extracted power. Therefore, it is more helpful to analyse the optimal trend of the extracted power obtained by Equation (13), in terms of the rotor speed and the acceleration/deceleration , which are normalized with respect to their design value: where and are the maximum rotor speed and acceleration respectively and so doing ∈ [0; 1] and ∈ [−1; 1]. Figure 3 shows how the power trend , changes for a fixed wind speed . In this case, is set to the maximum normalized rotor speed, namely = 1. It is apparent how for a generic data set related to the turbine parameters , , , , the extracted power is maximum only if the rotor seed reaches one third, as in Equation (11,

( , )
The best * path Similar considerations can be deduced by varying the wind speed and selecting a specific dimensionless rotor acceleration . V s , as shown in Figure 4. It is apparent how, regardless of the value of the rotor acceleration, for each value of the wind speed, the maximum values of P e , M, M 1 , M 2 , M 3 occur for a rotor velocity such that the optimal condition V * s = V 3 is observed. The law describing the trend of the maximum points is defined by the nonlinear control Equation (12). Similar considerations can be deduced by varying the wind speed and selecting a specific dimensionless rotor acceleration , as shown in Figure 4. It is apparent how, regardless of the value of the rotor acceleration, for each value of the wind speed, the maximum values of , , , , occur for a rotor velocity such that the optimal condition * = is observed. The law describing the trend of the maximum points is defined by the nonlinear control Equation (12).

Betz and OCBM Power and Efficiency Analysis
This section portrays a comparison between the performance of the new control system and that of a classical Betz's machine, not equipped with a regulation system. In both cases, the actual power is evaluated starting from the assumptions of non-stationary wind conditions. The case of a randomly blowing wind, whose speed V(t) is a random variable is considered. Of this flow, the expected value and the standard deviation { }, where {•} is the expected value operator, are assumed to be known. Accordingly, the wind time history is assumed to be decomposed into an average term and small fluctuations ( ) about it: For the random perturbation , the expected value and the standard deviation are: Equation (16), ( ) = + ( ) + 2 ( ), using Equation (17), produces: In the following paragraphs, the power of both methods is evaluated in the case of nonstationary wind, according to Equation (4), with an appropriate choice for the control action parameter . This coefficient has already been defined in Equation (12) for the OCBM machine but has not yet been defined for the steady wind turbine.

Betz and OCBM Power and Efficiency Analysis
This section portrays a comparison between the performance of the new control system and that of a classical Betz's machine, not equipped with a regulation system. In both cases, the actual power is evaluated starting from the assumptions of non-stationary wind conditions. The case of a randomly blowing wind, whose speed V(t) is a random variable is considered. Of this flow, the expected value V and the standard deviation E V 2 , where E{·} is the expected value operator, are assumed to be known.
Accordingly, the wind time history is assumed to be decomposed into an average term V and small fluctuations εV 1 (t) about it: For the random perturbation ε, the expected value and the standard deviation σ ε are: Equation (16), (17), produces: In the following paragraphs, the power of both methods is evaluated in the case of non-stationary wind, according to Equation (4), with an appropriate choice for the control action parameter c. This coefficient has already been defined in Equation (12) for the OCBM machine but has not yet been defined for the steady wind turbine.

The Steady Wind Turbine Power
The steady wind turbine theory considers a steady problem [2], i.e., a constant wind speed V, and produces the maximum power for In stationary conditions, namely constant rotor speed V s and constant control action parameter c Betz , the power is simply:P As previously mentioned, to calculate the actual power in non-stationary conditions, it is first necessary to evaluate c Betz . It is sufficient to recall the dynamics of the system as in Equation (3), for V > V s and applying the stationary condition for which . V s = 0, and the condition of maximum Betz power (19), leading to: At this point, the unsteady conditions can be introduced. Accordingly, the expected value of the classical Steady wind turbine power can be calculated, based on Equation (4): The only unknown term is E V 2 s , which can be deduced from the dynamic system Equation (3), including the wind speed perturbation Equation (16). This implies the solution V s must be not only function of V, but it should also account for its own fluctuation εV 1 (t), namely: Expanding the rotor speed V s (t, ε) in Taylor series up to the first order in terms of ε: and substituting the expansions for both velocities, i.e., Equations (16) and (24), into the Equation of motion (3), still in the case of V(t) > V s (t), the rotor dynamic becomes: By expressing the dependence on ε, and ordering with respect to its first two increasing powers, the cascade of equations is: For t → T , the first of Equations (26) tends to the stationary solution; indeed, excluding nonphysical solutions for V s > V and for a constant regulation c = c Betz , V s0 = V s = V 3 . Note that this first result already highlights an important feature of the expected Betz power: substituting in Equation (22), the series perturbation Equation (24) yields to: Eventually, the time average expressed by · symbol returns: Energies 2020, 13, 4913 11 of 20 Equation (28) unveils the double nature of the power P Betz , which is composed of a stationary contribution, already defined in Equation (20), and a term linked instead to the fluctuations, according to the assumptions related to the velocities of the wind and the rotor. To investigate how this fluctuating term varies, a frequency analysis is applied so to evaluate V s1 .
Substituting V s0 into the second of Equation (26), the fluctuation relation is obtained: By Fourier transforming Equation (29), V s1 (ω) is easily deduced: To return into the time domain, a harmonic wind perturbation V 1 (t) = V 1 (ω)e iωt is considered. As the rotor speed must be of the same nature, V s1 (t) = V s1 (ω)e iωt . The problem of the evaluation of the time average of the expected power E{P Betz } is solved by simply considering only the real part V s1 (t) = R V s1 (ω)e iωt and substituting Equation (30) into Equation (28):

The OCBM Power
As for the Betz case, the aim of this section is the evaluation of the time average of the expected power E{P OCBM } . The procedure results simplified with respect to the steady wind turbine case, as the expected OCBM power has already been calculated and the relation between the wind speed and the rotor speed defined. Indeed, the calculus of the maximum power P * e to be extracted from the system, solution of the general problem posed in Equation (6), returns the optimal solution V * s = V 3 and the optimal regulation c * , as in Equation (12). It follows: Applying the decomposition Equation (16) and time averaging, as for the Betz case, the last term E{P OCBM } is obtained: where the E . VV = 0 .

Steady wind Turbine and OCBM Comparison
The two analytical expressions in Equation (31) and (33) define the characteristic power extraction of the two methods, underlying how the intrinsic difference stands only in the fluctuation term, whilst the average term follows the stationary behavior. Based on this consideration, the following analysis focuses on the fluctuating component V 1 (ω) of the wind speed, in order to not only provide the conditions in which it is better to apply the OCBM controller rather than the steady wind turbine, but also an estimation of the amount of increased energy storage.
The optimal control problem solved for the OCBM machine automatically states greater energy storage than through the Steady wind turbine: Energies 2020, 13, 4913 12 of 20 where the ∼ apex is introduced to represent the respective power fluctuation E{P} . Recalling The solution of Equation (34) gives the following inequality: and it is easy to notice how the relation Equation (36) is always satisfied since ω > 0 and V > 0, namely the combined contribution of the fluctuations and of the optimal active control produces performances of the OCBM control always superior with respect to the Steady wind turbine.
To evaluate these increased performances, the relative produced energy E rel is considered: By substituting into Equation (37) the expression of the two energy sources as in Equation (31) and (33), the nondimensional form of E rel is obtained: with: The nondimensional parameter η provides a measure of the machine technical properties; indeed, it is the ratio between inertial forces and the inertia due to the wind action. On the other hand, the nondimensional parameters χ and ξ are purely related to the wind properties: they compare the frequency of oscillations and their amplitude with respect to the average wind speed, respectively. The two plots in Figure 5 present the variation of the trend of the relative energy gain with ξ, in the presence of harmonic, stationary wind disturbances: Figure 5a shows the different curves varying with χ, Figure 5b with η. In both plots, it is apparent how the higher the amplitude of the wind oscillations V 1 (ω) V, i.e., the lower the ξ, the higher the energy gain. This is explained by the intrinsic feature of the OCBM machine, intentionally designed so to follow the wind fluctuations. Moreover, in Figure 5a to the increase of the fluctuation frequency corresponds the increase of energy gain, instead in Figure 5b the reduction of η, namely when a large dimensions rotor is chosen, engenders a higher energy gain. A slender rotor is intrinsically more reactive and would autonomously adjust its dynamics according to the fluctuations, with no need for any active control and consequently no sensible energy gain. On the other hand, a heavier rotor, because of its higher inertia, needs the activation of the control to efficiently chase the wind oscillations.
The results confirm what was analytically observed, as the OCBM always generates a larger amount of energy with respect to the optimal steady wind turbine, especially with the increase of the wind fluctuations. The control logic driving the OCBM depends on the optimal parameter c * , which, modulated in time, varies assuming both positive and negative values. This implies the OCBM to either absorb or require energy, respectively. This variation, according to Equation (12), ensures the overall produced energy is still much larger than the contribution provided by the steady wind turbine. The results confirm what was analytically observed, as the OCBM always generates a larger amount of energy with respect to the optimal steady wind turbine, especially with the increase of the wind fluctuations. The control logic driving the OCBM depends on the optimal parameter * , which, modulated in time, varies assuming both positive and negative values. This implies the OCBM to either absorb or require energy, respectively. This variation, according to Equation (12), ensures the overall produced energy is still much larger than the contribution provided by the steady wind turbine.

Numerical Solutions
In this section, numerical simulations are performed for different types of wind excitations, confirming the reliability of the analytical results and the better achieved performances of the OCBM control with respect to the Betz one. The results are evaluated in terms of the nondimensional parameter , for given machine properties (see Table 1), i.e., fixed . The considered wind excitations are: (i) stationary wind flow with harmonic fluctuation, (ii) non-stationary wind in the shape of gust, and (iii) random wind flow.

Numerical Solutions
In this section, numerical simulations are performed for different types of wind excitations, confirming the reliability of the analytical results and the better achieved performances of the OCBM control with respect to the Betz one. The results are evaluated in terms of the nondimensional parameter ξ, for given machine properties (see Table 1), i.e., fixed η. The considered wind excitations are: (i) stationary wind flow with harmonic fluctuation, (ii) non-stationary wind in the shape of gust, and (iii) random wind flow. Note the expression for ξ as in Equation (39) is legitimate for harmonic fluctuations only, i.e., for a stationary wind flow, because of the impossibility to define the exact value for V 1 (ω). To overcome this hypothesis, also accounting for non-stationary wind conditions, a new expression for ξ needs to be introduced: where σ V 1 is the integral of the wind power spectral density: and G V 1 ( f ) represents the one-side power spectral density of V 1 (t) in a frequency bandwidth f up to the maximum frequency f max of the wind content. Given Equation (37), high values of ξ represent stationary wind conditions, whereas low values correspond to high gust wind. The selection of f max sets, not only the value of ξ, but also of χ.
In performing the numerical simulations, there is a further aspect not to be neglected: Equation (12), solution of the optimal control problem, does not respect the causality law, since . V s is not observable; numerically, this issue is solved by considering . V s recorded at previous time instants. This implies V s does not exactly follow at each time instant the optimal condition V s = V 3 , analytically obtained. Moreover, to prevent the rotor to turn in the opposite direction, a saturation condition over the control action c is introduced: Table 2 reports the applied conditions for each simulated case and summarizes how the increase of gust, that is decreasing ξ, produces the increase of relative energy storage evaluated with Equation (37), using the values obtained by the numerical simulations. Note that the wind speed average is evaluated over the entire simulation time, that amounts to about 5 h for each simulation, except for the random test which lasts almost 60 h. Each case of Table 2 is discussed in the figures from Figures 6-8, where the optimal control label represents the OCBM numerical solution of the system Equation (3) applying the c(t) control variable calculated in Equation (21). Instead, the Betz label indicates the passive solution obtained by using the parameter c Betz (12) evaluated through the Betz theory.
The first analyzed case, for high values of ξ, applies an input wind of V = 15 m s , harmonic fluctuations V 1 = ±8 m s , with 0.033 Hz frequency, as shown in Figure 6a. Figure 6b reports the wind signal Fourier transform, highlighting the sinusoidal component of the wind fluctuations. Figure 6c depicts the comparison between the instant power, which is constant in the Betz case, being the passive control design, and the OCBM case. It appears the power is fluctuating, taking both positive and negative values. In fact, the system we are proposing of an active type and the electrical machine can work as a generator or as a motor; however, the net energy is always positive and larger than the one produced by the classical Betz machine. Note that V s (t) follows the wind fluctuations and finds its optimum value for V V s 3, as in Figure 6d. Note how the discussed control system is of an active type. Therefore, the negative values of the instant power are physically meaningful, disclosing situations in which the rotor needs to be accelerated, as long as the overall energy stored remains positive, as demonstrated in Figure 6e. It shows the energy stored during the simulations unveiling a sensible, even though still limited, increase of energy storage of the OCBM with respect to the Betz case: The relative energy storage surplus, as reported in Table 2, amounts to 12%. its optimum value for ≅ 3, as in Figure 6d. Note how the discussed control system is of an active type. Therefore, the negative values of the instant power are physically meaningful, disclosing situations in which the rotor needs to be accelerated, as long as the overall energy stored remains positive, as demonstrated in Figure 6e. It shows the energy stored during the simulations unveiling a sensible, even though still limited, increase of energy storage of the OCBM with respect to the Betz case: The relative energy storage surplus, as reported in Table 2, amounts to 12%.  Figure 7 is related to the second study case of Table 2 in which the wind parameters are tuned so to perform a wind shape closer to a gust: a wind speed ratio = 0.5 is obtained with an average wind speed of 2.4 over the 5 h of simulation time; the fluctuations 1 = 3.4 are evaluated by the integral of the power spectral density, in Figure 7b, up to = 3 . As shown in Figure 7a, the gust is confined within a time window of about 2 h: starting from a rest condition, the wind abruptly increases to settle at an average wind speed of 10 , to drop back down at rest again. In the presence of gust conditions, the OCBM produces noticeably enhanced performances with respect to the steady wind turbine, as confirmed by the following figures. Indeed, Figure 7c shows the variation of the power in time: the quick variations of the OCBM case disclose the improved attitude to follow the wind fluctuations; accordingly, also the outcome average value of the power is much larger. This is further underlined by the comparison between the wind speeds in Figure 7d: after a sudden increase due to the inertia of the rotor, which starts from a null velocity, the OCBM promptly returns to the optimal condition that the steady wind turbine is not able to achieve. Figure 7e concludes with the comparison of the energy storage surplus and the evident overall gain of the OCBM.  Figure 7 is related to the second study case of Table 2 in which the wind parameters are tuned so to perform a wind shape closer to a gust: a wind speed ratio ξ = 0.5 is obtained with an average wind speed of 2.4 m s over the 5 h of simulation time; the fluctuations σ V 1 = 3.4 m s are evaluated by the integral of the power spectral density, in Figure 7b, up to f max = 3 Hz.
As shown in Figure 7a, the gust is confined within a time window of about 2 h: starting from a rest condition, the wind abruptly increases to settle at an average wind speed of 10 m s , to drop back down at rest again. In the presence of gust conditions, the OCBM produces noticeably enhanced performances with respect to the steady wind turbine, as confirmed by the following figures. Indeed, Figure 7c shows the variation of the power in time: the quick variations of the OCBM case disclose the improved attitude to follow the wind fluctuations; accordingly, also the outcome average value of the power is much larger. This is further underlined by the comparison between the wind speeds in Figure 7d: after a sudden increase due to the inertia of the rotor, which starts from a null velocity, the OCBM promptly returns to the optimal condition that the steady wind turbine is not able to achieve.   Table 2 in which the wind parameters are tuned so to perform a wind shape closer to a gust: a wind speed ratio = 0.5 is obtained with an average wind speed of 2.4 over the 5 h of simulation time; the fluctuations 1 = 3.4 are evaluated by the integral of the power spectral density, in Figure 7b, up to = 3 . As shown in Figure 7a, the gust is confined within a time window of about 2 h: starting from a rest condition, the wind abruptly increases to settle at an average wind speed of 10 , to drop back down at rest again. In the presence of gust conditions, the OCBM produces noticeably enhanced performances with respect to the steady wind turbine, as confirmed by the following figures. Indeed, Figure 7c shows the variation of the power in time: the quick variations of the OCBM case disclose the improved attitude to follow the wind fluctuations; accordingly, also the outcome average value of the power is much larger. This is further underlined by the comparison between the wind speeds in Figure 7d: after a sudden increase due to the inertia of the rotor, which starts from a null velocity, the OCBM promptly returns to the optimal condition that the steady wind turbine is not able to achieve. Figure 7e concludes with the comparison of the energy storage surplus and the evident overall gain of the OCBM.  The previous results have been produced with no constraints applied to the actuation forces and on velocity and acceleration of the rotor. Here, because of the highly unsteady wind profile, which is characterized by random fluctuations and provides a realistic model for the wind input, reported in Figure 8a, the control parameter is saturated such that the maximum acceleration of the rotor does not overcome 8 −2 force. In Figure 8b, the related power spectral density is displayed once again up to = 3 , as the further frequency content has a negligible contribution, and generates 1 = 1.4 . The energy gain in Figure 8c confirms the more efficient performances of the OCBM, driven by the optimal logic control and the stronger the gust, the larger the amount of absorbed energy: during The previous results have been produced with no constraints applied to the actuation forces and on velocity and acceleration of the rotor. Here, because of the highly unsteady wind profile, which is characterized by random fluctuations and provides a realistic model for the wind input, reported in Figure 8a, the control parameter c is saturated such that the maximum acceleration of the rotor does not overcome 8e −2 g force. In Figure 8b, the related power spectral density is displayed once again up to f max = 3 Hz, as the further frequency content has a negligible contribution, and generates σ V 1 = 1.4 m s . The energy gain in Figure 8c confirms the more efficient performances of the OCBM, driven by the optimal logic control and the stronger the gust, the larger the amount of absorbed energy: during strong variations of the wind, as for instance the one occurring at the 15th hour, the increase of energy storage extraction produced by the OCBM appears much larger than the Betz one.
Energies 2020, 13, x FOR PEER REVIEW 18 of 21 strong variations of the wind, as for instance the one occurring at the 15th hour, the increase of energy storage extraction produced by the OCBM appears much larger than the Betz one. In a nutshell, the presented cases confirm the reliability of the analytical results. A sensible increase of stored energy is observed with the increase of nonstationary wind conditions, characteristic of highly perturbed wind. This is achievable by the OCBM control thanks to its ability to chase fluctuations.

Conclusions
This paper investigates the problem of energy harvesting from a wind turbine. The introduced optimum control strategy has the intent to maximize the amount of energy absorption, identifying the optimal law for the regulation of an electrical machine so to obtain a maximum-efficiency process for the related energy storage. The proposed control strategy, based on the Euler-Lagrangian approach and belonging to the variational feedback control class, is defined even with an objective function characterized by the absence of an absolute maximum point and by the presence of a saddle point. Nevertheless, the obtained wind turbine control law, named OCBM, shows how to optimally store the energy, overcoming the well-known method for steady optimization of the wind turbine (by the Betz's theory). The OCBM machine maximizes the energy extraction in unsteady wind conditions by optimally following the wind oscillations, with best performances in the presence of strong and sudden gusts of wind. Numerical simulations compare Betz machine performances and the OCBM ones based on three nondimensional parameters accounting for machine properties and In a nutshell, the presented cases confirm the reliability of the analytical results. A sensible increase of stored energy is observed with the increase of nonstationary wind conditions, characteristic of highly perturbed wind. This is achievable by the OCBM control thanks to its ability to chase fluctuations.

Conclusions
This paper investigates the problem of energy harvesting from a wind turbine. The introduced optimum control strategy has the intent to maximize the amount of energy absorption, identifying the optimal law for the regulation of an electrical machine so to obtain a maximum-efficiency process for the related energy storage. The proposed control strategy, based on the Euler-Lagrangian approach and belonging to the variational feedback control class, is defined even with an objective function characterized by the absence of an absolute maximum point and by the presence of a saddle point. Nevertheless, the obtained wind turbine control law, named OCBM, shows how to optimally store the energy, overcoming the well-known method for steady optimization of the wind turbine (by the Betz's theory). The OCBM machine maximizes the energy extraction in unsteady wind conditions by optimally following the wind oscillations, with best performances in the presence of strong and sudden gusts of wind. Numerical simulations compare Betz machine performances and the OCBM ones based on three nondimensional parameters accounting for machine properties and wind conditions. These simulations, performed for different wind excitations, from harmonic to random fluctuations, confirm the analytical results.