A Model for Performance Estimation of Flapping Foil Operating as Biomimetic Stream Energy Device

: During the recent period intensive research has focused on the advancement of engineering and technology aspects concerning the development and optimization of wave and current energy converters driven by the need to increase the percentage of marine renewable sources in the energy-production mix, particularly from offshore installations. Most stream energy-harvesting devices are based on hydro-turbines, and their performance is dependent on the ratio of the blade-tip speed to incident-ﬂow speed. As the oncoming speed of natural-occurring currents varies randomly, there is a penalty for the latter device’s performance when operating at non-constant tip-speed ratio away from the design value. Unlike conventional turbines that are characterized by a single degree of freedom rotating around an axis, a novel concept is examined concerning hydrokinetic energy converters based on oscillating hydrofoils. The biomimetic device includes a rotating, vertically mounted, biomimetic wing, supported by an arm linked at a pivot point on the mid-chord. Activated by a controllable self-pitching motion the system performs angular oscillations around the vertical axis in incoming ﬂow. In this work, the performance of the above ﬂapping-foil, biomimetic ﬂow energy harvester is investigated by application of a semi-3D model based on unsteady hydrofoil theory and the results are veriﬁed by comparison to experimental data and a 3D boundary element method based on vortex rings. By systematical application of the model the power extraction and efﬁciency of the system is presented for various cases including different geometric, mechanical, and kinematic parameters, and the optimal performance of the system is determined.


Introduction
In nature, many species, such as aquatic animals, insects and birds, exploit energy directly from the fluid around them by controlling and maneuvering their bodies' locomotion via oscillation mechanisms either actively and/or passively controlled. For example, tuna fish, dolphin and sharks exhibit excellent hydrodynamic performance with high cruising speed, high efficiency and low noise by extracting water energy through their tail and/or fin's flapping motion.
Inspired by the above examples, biomimetic energy-harvesting devices have been developed for operation in the marine environment. The latter are based on the oscillatory motions of foils, which typically employ a fish-like profile modelled after the tail fins of marine animals; see [1,2]. Unlike conventional turbines, there are several prominent features of such bio-inspired energy converters: (i) They are environmentally friendly in terms of noise generation due to their relatively low tip speed, thus reducing impact on the communication and navigation of aquatic animals. (ii) Conventional energy harvesters of naturally-occurring water currents that are typically based on propeller-like wind or hydro turbines have scale limitations because they require relatively high speeds for proper operation. Thus, they perform rather poorly in the case of slow currents. 2 of 20 (iii) The impact of marine turbines on the marine ecology can be significant. The tips of turbine blades can reach high speeds, generating noise that will disturb wildlife.
An additional advantage of semi-activated biomimetic stream energy converters is that they perform oscillations with respect to more than one degree of freedom, and one oscillatory motion is controllable. This supports an optimized operation as compared to conventional hydro-turbines that are characterized by a single degree of freedom and their performance is dependent on the ratio of the blade-tip speed to incident-flow speed and the uniformity of the incoming flow. As the speed of natural-occurring currents varies randomly, there is a penalty for the latter device performance when operate at non-constant tip-speed ratio away from the design value.
Several design configurations have been proposed as biomimetic stream energy devices [1], with the main difference in the way the foil motion is constrained. One type represents systems which demand controlling/actuating the pitching motion. The existing flapping type energy harvesters in industry are often based on this design. Thereby, energy input is needed to activate the pitching motion, whereas energy harvesting is achieved through the resulting heaving motion generated by fluid dynamic lifting forces. Positive net energy extraction is possible only if the energy extracted from the heaving motion is higher than the energy expenditure to activate the pitching motion.
Another type represents systems which rely on flow-induced instabilities to generate oscillatory motions and vibrations in the heaving and pitching directions. This configuration simplifies the mechanical design since no activation system is needed. Furthermore, Kinsey and Dumas [3] considered a system with forced pitching and heaving motion. They investigated the power generation capability of a single NACA0015 foil, pivoted about a point one third of the chord downstream of the leading edge. Both heave and pitch motions were considered to be sinusoidal, fully prescribed and their phase difference was fixed at 90 • . Using appropriately the feathering limit notion (see [3], Figure 5) the transition from propulsion to energy extraction is defined. Zhu et al. [4] considered a semi-passive case in which the foil moves linearly in the heaving direction. They studied small amplitudes of motion of a zero-thickness plate using the two-dimensional analytical methodology of Theodorsen [5], and a three-dimensional boundary element model. The maximum power extraction was theoretically predicted and it was shown how to choose a non-dimensional damper strength parameter accordingly. In addition, it was found that the optimal performance of the system was where there was no spring on the plunge motion, and the foil pivot point was placed at the middle of the chord. Under these conditions the best frequency, pitch angle and damper strength parameter combinations yielded an efficiency of 25% based on plunging amplitude (which would be reduced when total swept distance was used). It was also illustrated that three-dimensional effect associated with small span-to-chord ratios would decrease the energy harvesting efficiency, while ground effect could enhance it. Moreover, Shimizu et al. [6] considered the same case as before, and the object of their study was to maximize the performance of the system. They used SOM, a data mining technique which is useful for visual presentation of multidimensional results and a Navier-Stokes code in order to evaluate the power and the efficiency taking viscous effects into account. Low-frequency and large amplitude heaving motions lead to large power extraction, whereas high frequency and small amplitudes result in high efficiency. They also showed that Navier-Stokes solutions gave better efficiency estimation up to 36.6%, in contrast with the potential flow predictions.
In the same direction of research, Huxham et al. [7] conducted a laboratory experiment in a water tunnel in the University of Sydney. The device consists of a vertically aligned aluminum NACA0012 foil with span equal to s = 0.34 m and chord c = 0.1 m, pivoted at a point located at 1/4 th of its chord from the leading edge, and linked to a vertical axis by an arm of 0.3 m length. Additionally, the edges of the foil are placed near the ceiling and the floor of the water tunnel in order to reduce the 3D effects and the induced drag. The heave of the foil is not prescribed but is dependent on the response of the system due to hydrodynamic forcing on the foil. The object of their study was to determine the influence of the pitch amplitude and reduced frequency on the power generation and the efficiency of the system. Many tests were conducted in order to determine the effects of these parameters with flow speed of U = 0.5 m/s and given power take-off characteristics. The above concept has been further exploited by the BioSTREAM device which is developed by BioPower Systems Pty Ltd., Sydney, Australia; see [8]. The above device is based on the principle of oscillating hydrofoils for tidal stream power generation. The foil is attached to an arm, and is therefore restricted to a circular motion. The component of the forces acting on the hydrofoil tangentially to the arm serves to either accelerate or decelerate the foil, while the hydrofoil resembles a fish tail. A prototype device is to be built and installed in a tidal stream off Flinders Island (see Figure 4 in Reference [8]).
In this work, the performance of the above flapping-foil, biomimetic flow energy harvester is investigated by development and application of a semi-3D model based on unsteady hydrofoil theory. The numerical results are verified by comparison to experimental data and a 3D boundary element method based on vortex rings. The considered biomimetic device operates by active pitch control of the self-pitching motion, in order to generate significant lift forces to drive the rotational degree of freedom associated with the transverse motion of the oscillating foil; see Figure 1, where also an indicative 3D view of the system is included. ment in a water tunnel in the University of Sydney. The device consists of a vertically aligned aluminum NACA0012 foil with span equal to s = 0.34 m and chord c = 0.1 m, pivoted at a point located at 1/4 th of its chord from the leading edge, and linked to a vertical axis by an arm of 0.3 m length. Additionally, the edges of the foil are placed near the ceiling and the floor of the water tunnel in order to reduce the 3D effects and the induced drag. The heave of the foil is not prescribed but is dependent on the response of the system due to hydrodynamic forcing on the foil. The object of their study was to determine the influence of the pitch amplitude and reduced frequency on the power generation and the efficiency of the system. Many tests were conducted in order to determine the effects of these parameters with flow speed of 0.5 U = m/s and given power take-off characteristics. The above concept has been further exploited by the BioSTREAM device which is developed by BioPower Systems Pty Ltd; see [8]. The above device is based on the principle of oscillating hydrofoils for tidal stream power generation. The foil is attached to an arm, and is therefore restricted to a circular motion. The component of the forces acting on the hydrofoil tangentially to the arm serves to either accelerate or decelerate the foil, while the hydrofoil resembles a fish tail. A prototype device is to be built and installed in a tidal stream off Flinders Island (see Figure 4 in Reference [8]).
In this work, the performance of the above flapping-foil, biomimetic flow energy harvester is investigated by development and application of a semi-3D model based on unsteady hydrofoil theory. The numerical results are verified by comparison to experimental data and a 3D boundary element method based on vortex rings. The considered biomimetic device operates by active pitch control of the self-pitching motion, in order to generate significant lift forces to drive the rotational degree of freedom associated with the transverse motion of the oscillating foil; see Figure 1, where also an indicative 3D view of the system is included. The combination of pitching and heaving motion of the oscillating wing sets the hydrofoil sections at an angle of attack relative to the incoming stream flow. This results in an unsteady lift force acting approximately perpendicular to the current direction, with strong component to the rotation direction around the vertical axis producing torque and power output. The magnitude of the lift force depends on various parameters such as the current velocity, the density of the medium in which the hydrofoil operates, its planform shape and surface area, the aspect ratio of the hydrofoil and its sectional characteristics, and the lift and drag coefficients which are dependent on the angle of attack. By systematic application of the present model the power extraction and efficiency of the system is presented for various cases including different geometric, mechanical, The combination of pitching and heaving motion of the oscillating wing sets the hydrofoil sections at an angle of attack relative to the incoming stream flow. This results in an unsteady lift force acting approximately perpendicular to the current direction, with strong component to the rotation direction around the vertical axis producing torque and power output. The magnitude of the lift force depends on various parameters such as the current velocity, the density of the medium in which the hydrofoil operates, its planform shape and surface area, the aspect ratio of the hydrofoil and its sectional characteristics, and the lift and drag coefficients which are dependent on the angle of attack. By systematic application of the present model the power extraction and efficiency of the system is presented for various cases including different geometric, mechanical, and kinematic parameters, and the results are verified by comparison to measured data and other methods, and the optimal performance of the system is determined.

Problem Formulation
We consider the unsteady oscillatory motion of the biomimetic wing shown in Figure 1, which operates in a uniform stream of velocity U. The foil performs simultaneously an angular oscillatory motion around the vertical axis Oz due to an enforced prespecified self-pitching angular oscillatory motion around the pivot axis PP. In this way, the foil operates as a semi-activated energy device extracting power from the incident flow. The oscillating hydrofoil is connected with the vertical rotation axis with an arm of length R from to Oz to the pivot point, and has two degrees of freedom defined by the angular motions θ and ϕ; see Figure 2. For small angles ϕ, the combined motion corresponds to enforced pitch θ and an induced heave Rϕ in the transverse y-direction, with phase difference around 90 deg, resembling the motion of a flapping wing. The pitch motion is prescribed, while the heave motion is induced by the interaction of the incidence flow with the hydrofoil at an instantaneous angle of attack α, taking into account the effect of the supporting system modeled as a spring, and the power take-off system modeled as a damper. The flapping wing is assumed to have symmetric sections and its cross-section is presented in Figure 2. A special case of wing is considered with an orthogonal planform shape of span s with sections of constant chord length c, where p is the location of the foil self-pitching axis. However, wings with more general planform area and geometrical characteristics can be considered. The arm R is mounted on a spring-damper base B modelling the Power-Take-Off.

Problem Formulation
We consider the unsteady oscillatory motion of the biomimetic wing shown in Figure 1, which operates in a uniform stream of velocity .
U The foil performs simultaneously an angular oscillatory motion around the vertical axis Oz due to an enforced prespecified self-pitching angular oscillatory motion around the pivot axis PP. In this way, the foil operates as a semi-activated energy device extracting power from the incident flow. The oscillating hydrofoil is connected with the vertical rotation axis with an arm of length R from Oz to the pivot point, and has two degrees of freedom defined by the angular motions θ and φ; see Figure 2. For small angles φ, the combined motion corresponds to enforced pitch θ and an induced heave Rφ in the transverse y-direction, with phase difference around 90 deg, resembling the motion of a flapping wing. The pitch motion is prescribed, while the heave motion is induced by the interaction of the incidence flow with the hydrofoil at an instantaneous angle of attack α, taking into account the effect of the supporting system modeled as a spring, and the power take-off system modeled as a damper. The flapping wing is assumed to have symmetric sections and its cross-section is presented in Figure 2. A special case of wing is considered with an orthogonal planform shape of span s with sections of constant chord length c , where p is the location of the foil self-pitching axis. However, wings with more general planform area and geometrical characteristics can be considered. The arm R is mounted on a spring-damper base B modelling the Power-Take-Off. As discussed, after the self-pitching motion of the foil is enforced, the hydrofoil undergoes an angular motion. For small responses, the system performs periodic oscillations of the angle of attack of the hydrofoil, and as a result, circulation is created, which is fundamental in order to produce lift and alternating torque around the vertical axis.
In order to formulate the dynamical equation of the system the following simplifying assumptions are made (see also Katz and Plotkin [9]): (i) The fluid is considered to be homogeneous (constant density) and the flow is ideal.
An additional assumption is that the flow leaves the trailing edge smoothly. In this sense the flow is considered to be irrotational and vorticity is concentrated on vortex sheets shed from the trailing edge, permitting flow with circulation around the wing sections. (ii) In addition, the wing span is large and the section thickness is relatively small permitting the approximate application of a semi-3D unsteady lifting flow model based As discussed, after the self-pitching motion of the foil is enforced, the hydrofoil undergoes an angular motion. For small responses, the system performs periodic oscillations of the angle of attack of the hydrofoil, and as a result, circulation is created, which is fundamental in order to produce lift and alternating torque around the vertical axis.
In order to formulate the dynamical equation of the system the following simplifying assumptions are made (see also Katz and Plotkin [9]): (i) The fluid is considered to be homogeneous (constant density) and the flow is ideal.
An additional assumption is that the flow leaves the trailing edge smoothly. In this sense the flow is considered to be irrotational and vorticity is concentrated on vortex sheets shed from the trailing edge, permitting flow with circulation around the wing sections. (ii) In addition, the wing span is large and the section thickness is relatively small permitting the approximate application of a semi-3D unsteady lifting flow model based on a stripwise along the span application of 2D sectional results based on Theodorsen unsteady thin hydrofoil theory [9,10]. (iii) Finally, it is considered that extended flow separation does not manifest and viscousflow effects are small and can be approximately considered through sectional drag coefficients expressed in terms of empirical relations based on Reynolds number and instantaneous angle of attack.
The foil is subjected to the prescribed pitching motion θ(t) and the induced transverse motion h(t) = R ϕ(t), that appears as transverse heaving motion; see also the sectional plot of Figure 2. Thus, assuming that the fluid surrounding the foil and the wake is inviscid, incompressible and irrotational and that each section of small thickness wing undergoes small-amplitude pitching and heaving harmonic oscillations, the hydrodynamic lift force per unit span is obtained as the combination of an unsteady lift force and added mass effects; see Katz and Plotkin [9], Section 13.8. The sectional lift coefficient is given by where an over-dot denotes time-differentiation, k = ωc 2U is the reduced frequency, and C(k) is a coefficient determined through the Theodorsen function. The latter under the assumption of harmonic oscillatory motion with an exp(iωt) time dependence is expressed as a ratio of the Hankel functions of the second kind as follows: see Katz and Plotkin [9], Newman [10]. Moreover, the moment per unit span required to pitch the foil with respect to the pitching axis, is approximated by the following equation (see [9], Section 13.8): The corresponding 3D coefficients are obtained by a stripwise integration along the span of the wing where A is the wing planform area and c m the midchord length. In addition, C 3D is a correction factor taking into account 3D effects which is usually estimated from steady lifting line theory of elliptical wings as C 3D = AR/(AR + 2) where AR = s 2 /A denotes the aspect ratio of the wing. Obviously, for a flapping wing with an orthogonal planform shape, with constant spanwise chord length distribution, it holds C L,

The Dynamic System Modelling the Time-Domain Performance of the Biomimetic Device
The hydrofoil is connected through the arm R to the vertical axis and the power-takeoff system modelled by using a rotational spring of constant C and damper with coefficient B; see Figure 2. The dynamics governing the oscillatory motion of the considered device is governed by the balance of angular momentum about the vertical axis Oz, where ϕ(t) represents the angular response of the system (see Figures 2 and 3), I is the moment of inertia of the foil including the connecting arm, B is the damping coefficient representing the energy extraction from the generator (PTO), C is the rotational spring constant, and T represents the generated torque.
the reduced frequency  To proceed we assume that that the response of the system activated by the self-pitching motion (defined by Equation 6) is a function of time and, thus, the generated torque is equal to: where F θ is the projection of the hydrodynamic forces in the tangential direction (see Figure 2).
As already discussed in the above, the present model is based on the unsteady thin hydrofoil theory. In order to take into account the leading-edge suction force effects the additional assumption is made that the instantaneous lift force on the hydrofoil acts in the transverse y-direction as shown in Figure 3, and thus, the first term in the right-hand-side of Equation (13) is given by: and the corresponding component due to the drag force acting in the horizontal x-direction by: Figure 3. Sketch of the biomimetic stream energy converter with force and velocity vectors.
The power extraction is achieved by the generated torque which is essentially provided by the lift force of the dynamic foil and affected by the corresponding drag force. Thus, maximizing the torque can increase the amount of the output power. We next consider the self-pitching oscillation of the foil defined as a simple periodic function of time of angular frequency ω and amplitude θ 0 : where g(t) is a filter function enabling the smooth enforcing of the self-pitching motion, the generation of lifting flow with angle of attack and the angular response of the system. Moreover, the derived transverse heaving motion of the foil is: which in the steady-state condition and small responses |ϕ| << 1 results in harmonic oscillations with amplitude ϕ 0 and induced transverse (heaving) motion of amplitude h 0 ≈ Rϕ 0 ; see Figure 2.
As the foil oscillates in the uniform current, the relative flow velocity due to the rotational motion of the wing about the vertical axis Oz, in conjunction with the incident parallel flow, is: where U f in is the corresponding rotational velocity of the wing: ϕe ϕ (9) acting in the tangential direction with e ϕ denoting the corresponding unit vector (see Figure 2). On the other hand, the self-pitching motion of the foil is taken into account for the formation of the angle of attack, and results in the change of the effective angle of attack of the flow relative to the foil is given by: where the latter component describes the effect of the induced transverse heaving motion due to the rotational along the vertical axis oscillatory motion of the wing: Based on the above, the most important parameters concerning the response and the performance of the examined system are: the dimensionless heave amplitude h 0 /c, -the amplitude of the enforced self-pitching motion θ 0 , -the location of the foil pivot axis with respect to the foil leading edge p/c.
In addition, the most important flow dimensionless parameters, which define the dynamic behavior of the system are: the reduced frequency k = ωc 2U and the Strouhal number the amplitude of the effective angle of attack α 0 , -the Reynolds number Re = Uc υ through which, in conjunction with the instantaneous angle of attack, the drag forces are estimated To proceed we assume that that the response of the system activated by the selfpitching motion (defined by Equation (6)) is a function of time and, thus, the generated torque is equal to: where F θ is the projection of the hydrodynamic forces in the tangential direction (see Figure 2).
As already discussed in the above, the present model is based on the unsteady thin hydrofoil theory. In order to take into account the leading-edge suction force effects the additional assumption is made that the instantaneous lift force on the hydrofoil acts in the transverse y-direction as shown in Figure 3, and thus, the first term in the right-hand-side of Equation (13) is given by: and the corresponding component due to the drag force acting in the horizontal x-direction by: The lift force F L is estimated through the lift coefficient by means of Equations (1) and (4), and the drag force is obtained by a drag coefficient calculated by means of the following empirical formula, where the first term C f (Re) is dependent on the Reynolds number and represents skin friction effects, and the second term models the effect of the dynamic angle of attack α(t) expressed through the empirical coefficient C a ; see also Filippas and Belibassakis [11] and the references cited there. Substituting the above equations into the equation of motion, Equation (5), in the simple case of a wing with orthogonal planform (constant chord distribution along the span) we finally obtain the following dynamic equation modelling the angular response and the performance of the examined system: Since the power take-off system is idealized as a simple damper, the instantaneous power output is calculated by: and the time average power output can be calculated using the results of the system, which in the case of harmonic response at angular frequency ω is: The external moment M required to drive the enforced pitching motion equals the hydrodynamic moment, Equations (3) and (4), and thus, we can calculate the power input required to rotate the hydrofoil. The average power spent in order to pitch the foil is expressed by: Before continuing it should be noted that is useful to choose the appropriate measure for the performance of the examined system. In the relevant bibliography, there are two definitions for the performance of the device, whose difference can be significant, especially if large pitching amplitudes are considered. The difference between the two concerns the cross-sectional area normal to the velocity of the current swept by the foil considered. The first definition uses the area swept by the pivot axis, while the second takes into account the whole movement of the foil, including the pitching motion. The energy extraction efficiency η is chosen in the present study to be defined as the ratio between the cycle-averaged power extracted and the energy flux of the incoming flow, as follows: where d denotes the distance between the extreme points in the transverse direction during the oscillatory motion of the system, and the quantity 1 2 ρU 3 sd is the incoming kinematic power flux of the uniform incident flow through the section swept by the oscillating hydrofoil.

Numerical Results
In order to study the effect of various parameters on the response of the system including the the effect of non-linearities, such as the terms involving the effect of the response on the direction of the tangential forces generating the torque, and the quadratic terms associated with the drag force acting on the foil, the equation of motion is integrated starting from rest. For this purpose, defining the new variables ϕ 1 (t) = ϕ(t), ϕ 2 (t) = . ϕ(t), the dynamic system is put in the following first-order form: where the right-hand side (forcing) of the dynamical system is defined by: The above system is numerically integrated by using a 2-3rd order Runge-Kutta method in Matlab ® , which provides an easy to implement yet powerful tool for numerical solution of the initial value problems governed by the ordinary differential equation (17).
In order to validate the accuracy of the method, the results are first compared against the experimental data provided by Huxham et al. [7]. For this purpose, a hydrofoil with NACA0012 sections was modelled. Following the experimental set-up [7], the foil had constant chord distribution along the span which was equal to s = 0.34 m, with chordlength c = 0.1 m. In the specific arrangement the pivot point was located at p/c = 1/4 of its chord, and was connected to the vertical axis by using an arm of R = 0.3 m length. The aspect ratio of the foil was, therefore, AR = 3.4. However, the foil edges were located very close to the water tunnel ceiling and floor which perform like end-plates in order to minimize the tip vortices. This results in the reduction of the induced drag and the 3D effects on the flow, and consequently the efficiency of the configuration and the performance of device are increased. The lift and drag of the foil create a moment at the axis of rotation of the arm and a sensor is placed at the shaft to measure this moment (see Figure 2 in Reference [7]). Two identical torque transducers were used to measure the torque needed for the self-pitching motion of the foil. A dashpot was placed on the arm shaft to dampen the motion modelling and also the power take-off generator.
In order for the present method to better simulate the experiment, an artificial increased value of the aspect ratio was considered to take into account the effect of the nearby boundaries in the vicinity of wingtips on the reduction of the tip wake vortices. Following Koutsogiannakis et al. [12], an artificial aspect ratio AR = 10 was used in the numerical simulations, and the results obtained were compared with the measured data provided by [7].
In Figure 4, the response of the system was illustrated, as calculated by the present method. Results are obtained by considering no spring constant (C = 0) and a dimensionless damping coefficient B = B/ 1/2ρUsc 2 R = 29.5, as reported in [7]. More specifically, the time history of the angle of attack and the enforced self-pitching angle are shown, for onset flow speed U = 0.5 m/s, and reduced frequency k = 0.1. It is observed in Figure 4 that the amplitude of the response is calculated to be ϕ 0 = 7.85 • , and the corresponding amplitude of the effective angle of attack is found to be α 0 = 37.54    Figure 5, illustrates the calculated ratio of the non-dimensional heaving response normalized to the chord of foil for reduced frequency k = 0.1, and self-pitching amplitude θ 0 = 50 • . Again, the results obtained by the current model after the numerical response has converged to periodic state are found to be very close to the experimental data indicated by using symbols. The maximum heave response, (corresponding to the amplitude of oscillatory motion of the device in the transverse direction) is calculated to be h/c = 0.41, and compares quite well with the experimental value which is h/c = 0.43; see reference [7].  Figure 4. Experimental data from [7] are indicated by using symbols.
Finally, the calculated time-history of the input and output power are presented for the same as before case in Figure 7. It is clearly observed in these figures that the present method results are in good agreement with the experimentally measured data that are available from [7] when the system has reached its steady-state (periodic) response.
It is seen in these figures that the results obtained from the present method indicate an underestimation of the self-pitching torque and consequently also of the input power. This is mainly due to the calculation of the drag force which is based on a simplified empirical model, Equation (16), and does not accurately treat flow separation effects.
More specifically, for the considered configuration, flow and operation conditions, the power required to drive the self-pitching mechanism of the foil is 0.51 W and the output power is 1.42 W, and as a result the calculated efficiency is 26.1%, which is increased by a relatively small amount as compared to the experimental value which is reported to be 23%. Nonetheless, the present simplified method is shown to represent well the dynamic response and performance of the considered system.  Figure 4. Experimental data from [7] are indicated by using symbols.
Moreover, for the same as before case, the output and input torque are depicted in Figure 6 for the case of reduced frequency k = 0.1 and pitching amplitude a 0 = 50 • . The maximum output and input torque are 3.28 Nm/rad/sec and 0.18 Nm/rad/sec, respectively, and they are quite good compared with the experimental data by [6] (indicated by using symbols).
Finally, the calculated time-history of the input and output power are presented for the same as before case in Figure 7. It is clearly observed in these figures that the present method results are in good agreement with the experimentally measured data that are available from [7] when the system has reached its steady-state (periodic) response. Finally, the calculated time-history of the input and output power are presented for the same as before case in Figure 7. It is clearly observed in these figures that the present method results are in good agreement with the experimentally measured data that are available from [7] when the system has reached its steady-state (periodic) response. Figure 6. Input and output torque as a function time as calculated by the present model for the considered system for the same case as in Figure 4. Input torque is shown by using dashed lines and output torque by solid lines, respectively. Experimental data from [7] are indicated by using symbols.  Figure 4. Input power is shown by using dashed lines and output power by solid lines, respectively. Experimental data from [6] are indicated by using symbols: crosses for input power and circles for output power, respectively.

Evaluation of the Dynamic Model through Application of Unsteady Vortex Element Method
The present dynamic model described in detail above, is based on a semi-3D stripwise application of unsteady hydrofoil theory. In order to investigate the significance of 3D flow characteristics associated with finite aspect ratio, induced drag, and nonlinear trailing vortex sheet dynamics and their effects on the simulation of the performance of the device, in this section the results obtained by the time-domain model are compared  Figure 4. Input torque is shown by using dashed lines and output torque by solid lines, respectively. Experimental data from [7] are indicated by using symbols.

Figure 7.
Input and output power as a function time as calculated by the present model for the considered system for the same case as in Figure 4. Input power is shown by using dashed lines and output power by solid lines, respectively. Experimental data from [6] are indicated by using symbols: crosses for input power and circles for output power, respectively.

Evaluation of the Dynamic Model through Application of Unsteady Vortex Element Method
The present dynamic model described in detail above, is based on a semi-3D stripwise application of unsteady hydrofoil theory. In order to investigate the significance of 3D flow characteristics associated with finite aspect ratio, induced drag, and nonlinear trailing vortex sheet dynamics and their effects on the simulation of the performance of the device, in this section the results obtained by the time-domain model are compared  Figure 4. Input power is shown by using dashed lines and output power by solid lines, respectively. Experimental data from [6] are indicated by using symbols: crosses for input power and circles for output power, respectively.
It is seen in these figures that the results obtained from the present method indicate an underestimation of the self-pitching torque and consequently also of the input power. This is mainly due to the calculation of the drag force which is based on a simplified empirical model, Equation (16), and does not accurately treat flow separation effects.
More specifically, for the considered configuration, flow and operation conditions, the power required to drive the self-pitching mechanism of the foil is 0.51 W and the output power is 1.42 W, and as a result the calculated efficiency is 26.1%, which is increased by a relatively small amount as compared to the experimental value which is reported to be 23%. Nonetheless, the present simplified method is shown to represent well the dynamic response and performance of the considered system.

Evaluation of the Dynamic Model through Application of Unsteady Vortex Element Method
The present dynamic model described in detail above, is based on a semi-3D stripwise application of unsteady hydrofoil theory. In order to investigate the significance of 3D flow characteristics associated with finite aspect ratio, induced drag, and nonlinear trailing vortex sheet dynamics and their effects on the simulation of the performance of the device, in this section the results obtained by the time-domain model are compared with numerical predictions obtained by application of the 3D unsteady vortex element method including non-linear free-wake analysis.
We consider the same geometrical configuration as that considered in the experimental study by Huxham et al. [7], as previously shown in the examples of Figures 4 and 5. In order to illustrate the ability of present semi-3D model to represent the effects of the 3D flow we consider oscillating pitching and heaving wing of trapezoidal planform and swept-back wing, with the same dimensions as the experimental configuration, in infinite domain, without any effects from nearby boundaries in the tip region, as shown in Figures 8 and 9, respectively. For the simulation we have used amplitudes of the selfpitching motion θ 0 = 50 • and heaving motion h/c = 0.41 and phase difference between the two oscillatory motions the same as were obtained in the test case examined in the previous section (see also Figure 10), for flow incidence U = 0.5 m/s, reduced frequency k = 0.1 and Strouhal number St = 0.082. The 3D analysis is obtained by application of the unsteady vortex element method (VEM); see Belibassakis et al. [13]. The latter method is based on the assumption that the flow around the oscillating wing(s) is inviscid and incompressible. As the wing is steadily translated, the vorticity created in the boundary layers of its upper and lower surfaces is continuously shed into the wake. Furthermore, this vorticity, which is subsequently convected away from the wing with the local velocity, is assumed to be concentrated into sheets of infinitesimal thickness, constituting the trailing vortex wakes. Under the previous assumptions the problem, consists of the Laplace equation subject to the following conditions (see also Katz and Plotkin [9]): (i) Solid boundary condition. The normal component of the relative velocity (w) on the surface of the wing must vanish. (ii) Trailing vortex sheet boundary conditions. Trailing vortex sheets must be force-free material surfaces and the vorticity field (concentrated on these surfaces) must be divergence-free. Shed vorticity is downstream convected with the local velocities as calculated at the nodal points of each free vortex ring. (iii) Kutta condition. Among the alternative approximate forms of the Kutta condition in unsteady flow, in the present work this condition is implemented by enforcing the vorticity generated at the trailing edge and the tips of the wing to shed with the local velocity. (iv) Quiescence condition at infinity.      The results in one period of oscillation concerning the time history of motion and the calculated 3D lift coefficient are plotted in Figure 10. In particular, the calculated coefficient from the non-linear 3D vortex element method is plotted for the trapezoidal swept-back wing by using a solid line and the corresponding result concerning the rectangular wing by using dashed line, respectively. The latter result is compared with the prediction by the present semi-3D model obtained by means of Equation (4), which is indicated by using a red line, and is shown to produce comparable prediction. We also observe in this figure that the swept back wing, resembling a fish-tail configuration for the studied biomimetic system, provides increased lift force that could be found to be beneficial concerning the performance of the studied system in a realistic environment, and this finding is left to be further studied in future work using a fully 3D method.

A linearized Model with Application to the Performance Characteristics of the Studied Device
Since the present semi-3D method is shown to provide good results (compared to experiments and 3D models) in the case of oscillating foils of simple rectangular planform shape, in this section a linearized version of the model is derived in order to be systematically applied in the frequency domain and study the effect of various parameters concerning the performance optimization of the system.
Considering harmonic time dependence, the hydrofoil is subjected to a prescribed self-pitching motion represented by: where ω is the angular frequency and 0 θ the self-pitching amplitude .
Furthermore, considering the range of responses where the non-linear effects can be approximately omitted, the angular response of the system is described by: The problem is treated by implementing a boundary integral representation for the disturbance velocity field, and is numerically solved by a 3D panel method consisted of quadrilateral vortex elements and discrete source elements. More details are available in references [13,14].
Plots concerning the vortex wake development downstream the oscillating wing, as calculated by the vortex element method using free wake analysis, staring from rest and enforcing the kinematic conditions due to the combined oscillatory motions heave (in the transverse direction) and self-pitching motion, are shown in Figures 8 and 9. In particular, in Figure 8 the wing with rectangular trapezoidal planform shape with aspect ratio AR = 3.40 corresponding to the configuration by Huxham et al. [7] is considered, and in Figure 9 the corresponding wing with trapezoidal planform area and same main dimensions, but with a leading edge swept angle Λ = 70 • , leading to increased aspect ratio AR = 5, and about 42% reduced wing planform area, is shown. The midchord length of both wings is the same c = 0.1 m, as well as the span which is 0.34 m. The intensity of the vortex rings in the wake is illustrated by using color plots. We observe in these figures that the effect of vortex sheet separation from the tips is reduced in the case of the swept back wing which, in conjunction with the increased aspect ratio, leads to increased lift coefficient.
The results in one period of oscillation concerning the time history of motion and the calculated 3D lift coefficient are plotted in Figure 10. In particular, the calculated coefficient from the non-linear 3D vortex element method is plotted for the trapezoidal swept-back wing by using a solid line and the corresponding result concerning the rectangular wing by using dashed line, respectively. The latter result is compared with the prediction by the present semi-3D model obtained by means of Equation (4), which is indicated by using a red line, and is shown to produce comparable prediction. We also observe in this figure that the swept back wing, resembling a fish-tail configuration for the studied biomimetic system, provides increased lift force that could be found to be beneficial concerning the performance of the studied system in a realistic environment, and this finding is left to be further studied in future work using a fully 3D method.

A Linearized Model with Application to the Performance Characteristics of the Studied Device
Since the present semi-3D method is shown to provide good results (compared to experiments and 3D models) in the case of oscillating foils of simple rectangular planform shape, in this section a linearized version of the model is derived in order to be systematically applied in the frequency domain and study the effect of various parameters concerning the performance optimization of the system.
Considering harmonic time dependence, the hydrofoil is subjected to a prescribed self-pitching motion represented by: where ω is the angular frequency and θ 0 the self-pitching amplitude. Furthermore, considering the range of responses where the non-linear effects can be approximately omitted, the angular response of the system is described by: where Re(.) denotes the real part of complex quantity. We consider that the above angular motion is small, |ϕ 0 | << 1, and corresponds to a harmonic oscillatory motion in the transverse direction: and the corresponding harmonic heaving motion is where h 0 = Rϕ 0 is the complex heaving amplitude. Thus, for small oscillatory motion amplitudes the angle of attack is similarly described as: Using the above the hydrodynamic sectional lift force is obtained in the following form: and as a result, the following relation is obtained for the harmonic response of the system from Equation (17) and omitting the effects of non-linear terms: In the sequel the above linearized model derived in the frequency domain is first verified by comparing results with numerical simulations by the time domain model Equation (17) and the experimental data provided by Huxham et al. [7]. Subsequently, the effects of various parameters on the performance of the studied device will be presented, as obtained by the systematic application of the systematic linearized model.
First, for the device operating in a uniform stream U = 0.5 m/s, and the foil (with dimensions and geometrical parameters as above in Section 3) oscillating at reduced frequency k = 0.1 with self-pitching amplitude θ 0 = 50 • , in Figure 11 the effective angle of attack and the angular response of the system are presented, as calculated from the present time-domain model using solid lines and from the linearized frequency domain model by using dashed lines, respectively. The corresponding comparison concerning the input and output power of the device as a function of time is presented in Figure 12. The results show very good agreement indicating that the simplified frequency domain model can provide very good predictions, at least for the cases of geometrical data and operating conditions, where the omitted non-linearities do not expect to have strong effects.  In the case examined the dimensionless heave amplitude is h/c = 0.39, and is found to be a little smaller than the time-domain prediction but still in good agreement with the  In the case examined the dimensionless heave amplitude is h/c = 0.39, and is found to be a little smaller than the time-domain prediction but still in good agreement with the experimental data provided in [7]. Also, the amplitude of the effective angle of attack is In the case examined the dimensionless heave amplitude is h/c = 0.39, and is found to be a little smaller than the time-domain prediction but still in good agreement with the experimental data provided in [7]. Also, the amplitude of the effective angle of attack is found to be α 0 = 36.6 • which is a little smaller than the value α 0 = 37.54 • obtained by the present time domain simulation. Similarly, the amplitude of the response of the system is found to be ϕ 0 = 7.51 • and 7.85 • by the frequency and the time domain analysis, respectively, and the above justify also the small differences observed in the power calculations presented and compared in Figure 12. In the case examined the results obtained by the time domain model are closer to the experimental data, as reported in the work by Huxham et al. [7].
In the examined case, the power extraction and the power input to drive the foil self-pitching mechanism, plotted in Figure 12 as obtained by the frequency domain model, are found to be 1.27 W and 0.572 W, respectively and, thus, the estimated performance for the particular operating condition is 21.38%. The latter small discrepancy is due to the difference in the predicted responses, and the estimation of the performance is smaller than the corresponding predicted value from the time-domain model which is 26.1%. On the other hand, the experimental value is reported to be 23%, which lies between the two calculated values; see references [7]. The predictive capability of the present linearized model in the frequency domain is further illustrated in Figure 13 where the calculated dimensionless heave amplitude (h/c) is presented against the self-pitching amplitude θ 0 , for several values of the reduced frequency k = 0.025, 0.05, 0.1, 0, 0.15, 0.2, and compared with the experimental data reported by Huxham et al. [6] and shown by using symbols. respectively, and the above justify also the small differences observed in the power calculations presented and compared in Figure 12. In the case examined the results obtained by the time domain model are closer to the experimental data, as reported in the work by Huxham et al. [7]. In the examined case, the power extraction and the power input to drive the foil self-pitching mechanism, plotted in Figure 12 as obtained by the frequency domain model, are found to be 1.27 W and 0.572 W, respectively and, thus, the estimated performance for the particular operating condition is 21.38%. The latter small discrepancy is due to the difference in the predicted responses, and the estimation of the performance is smaller than the corresponding predicted value from the time-domain model which is 26.1%. On the other hand, the experimental value is reported to be 23%, which lies between the two calculated values; see references [7]. The predictive capability of the present linearized model in the frequency domain is further illustrated in Figure 13 where the calculated dimensionless heave amplitude (h/c) is presented against the self-pitching amplitude 0 θ , for several values of the reduced frequency 0.025,0.05,0.1,0,0.15,0.2 k = , and compared with the experimental data reported by Huxham et al. [6] and shown by using symbols.
The above results indicate that the present simplified method could represent relatively well the dynamic response and performance of the considered system and, by taking advantage of its analytical structure, it could be very conveniently applied to study the effects of various parameters on the performance of the device. In particular, a contour map of efficiency as calculated by the present linearized method is presented in   The above results indicate that the present simplified method could represent relatively well the dynamic response and performance of the considered system and, by taking advantage of its analytical structure, it could be very conveniently applied to study the effects of various parameters on the performance of the device. In particular, a contour map of efficiency as calculated by the present linearized method is presented in Figure 14. Systematic results for U = 0.5 m and damping coefficient C = 29.5 modelling the PTO are presented as a function of pitching amplitude θ 0 and reduced frequency k, and compared with experimental data from Huxham et al. [7]. The results obtained by frequency-domain model are indicated with solid lines and the experimental data using circles.
(a) (b) It is observed in Figure 14 that the present model predicts very well the dependence of the energy device performance from the basic operating parameters for low and moderate values of the reduced frequency up to k = 0.16, and the optimal value k = 0.1 is consistent with numerical research of foils undergoing a prescribed pitch and heave; see also Zhu [15].
For higher values of reduced frequency, significant non-linear phenomena take place limiting the predictive capability of the present model, requiring applications of more sophisticated methods including modelling of dynamic stall effects. However, we observe in Figure 14 that in the subregion of parameters where efficient performance of the device occurs the present model still provides good predictions and could be exploited for the optimal selection of parameters concerning the design and operation of the system.
Following this remark, in the sequel a parametric study of the device is undertaken, and results are presented and discussed concerning the effect of the location of the pivot axis for the self-pitching motion of the foil and for the same value of the damping coefficient modelling the PTO. In Figure 15 , and the maximum performance is calculated to be 26.5%, as in the case presented in Figure 14. Moreover, the positioning of the pivot axis is preferably set at a distance of about 35% of the chord from the leading edge. The small increase of performance for this value of the pivot axis is attributed to the small reduction of input power required due to the corresponding reduced moment for p/c = 0.3-0.35, which is probably closer to the center of It is observed in Figure 14 that the present model predicts very well the dependence of the energy device performance from the basic operating parameters for low and moderate values of the reduced frequency up to k = 0.16, and the optimal value k = 0.1 is consistent with numerical research of foils undergoing a prescribed pitch and heave; see also Zhu [15].
For higher values of reduced frequency, significant non-linear phenomena take place limiting the predictive capability of the present model, requiring applications of more sophisticated methods including modelling of dynamic stall effects. However, we observe in Figure 14 that in the subregion of parameters where efficient performance of the device occurs the present model still provides good predictions and could be exploited for the optimal selection of parameters concerning the design and operation of the system.
Following this remark, in the sequel a parametric study of the device is undertaken, and results are presented and discussed concerning the effect of the location of the pivot axis for the self-pitching motion of the foil and for the same value of the damping coefficient modelling the PTO. In Figure 15, a contour map of efficiency is presented as function of the reduced frequency k and the pitching amplitude θ 0 for two other locations of pivot axis at p/c = 0.3 ( Figure 15a) and p/c = 0.35 (Figure 15b). It is seen that good efficiency is achieved for a range of reduced frequencies k = 0.06-0.11, and the maximum performance is calculated to be 26.5%, as in the case presented in Figure 14. Moreover, the positioning of the pivot axis is preferably set at a distance of about 35% of the chord from the leading edge. The small increase of performance for this value of the pivot axis is attributed to the small reduction of input power required due to the corresponding reduced moment for p/c = 0.3-0.35, which is probably closer to the center of hydrodynamic pressure of the foil in this case for the system in operating conditions around the optimal values. hydrodynamic pressure of the foil in this case for the system in operating conditions around the optimal values. We conclude that the present frequency domain analysis is a useful tool that can reproduce, at a first level of approximation, the behavior of the system, providing useful results for the optimal selection of parameters at the preliminary design stage.

Conclusions
In this work, the performance of a biomimetic stream energy harvester is investigated through numerical modelling. The biomimetic device includes a rotating, vertically mounted flapping wing, supported by an arm linked at a pivot point on its mid-chord section. Activated by an enforced self-pitching motion, the system performs angular oscillations around the vertical axis in parallel flow. For the unsteady analysis of the above flow energy harvester, a semi-3D model based on unsteady hydrofoil theory is developed and the results are verified by comparison to experimental data and a 3D boundary element method based on vortex rings, showing that the present method provides good results in the case of oscillating foils of simple rectangular planform shape. Next, a linearized version of the model was derived in the frequency domain and systematically applied to study the effect of various parameters concerning the performance optimization of the system. Future work includes the investigation of incident flow disturbances on the performance of the studied device [15,16]. Also, the application of 3D methods for the investigation of more general fish-like wing configurations, including interactions of multiple devices in array arrangement [17], and the effects of additional boundaries like the free-surface and the seabed, will be exploited for the study of the operations and performance optimization of this type of renewable energy device in the marine and coastal environment.  We conclude that the present frequency domain analysis is a useful tool that can reproduce, at a first level of approximation, the behavior of the system, providing useful results for the optimal selection of parameters at the preliminary design stage.

Conclusions
In this work, the performance of a biomimetic stream energy harvester is investigated through numerical modelling. The biomimetic device includes a rotating, vertically mounted flapping wing, supported by an arm linked at a pivot point on its mid-chord section. Activated by an enforced self-pitching motion, the system performs angular oscillations around the vertical axis in parallel flow. For the unsteady analysis of the above flow energy harvester, a semi-3D model based on unsteady hydrofoil theory is developed and the results are verified by comparison to experimental data and a 3D boundary element method based on vortex rings, showing that the present method provides good results in the case of oscillating foils of simple rectangular planform shape. Next, a linearized version of the model was derived in the frequency domain and systematically applied to study the effect of various parameters concerning the performance optimization of the system. Future work includes the investigation of incident flow disturbances on the performance of the studied device [15,16]. Also, the application of 3D methods for the investigation of more general fish-like wing configurations, including interactions of multiple devices in array arrangement [17], and the effects of additional boundaries like the free-surface and the seabed, will be exploited for the study of the operations and performance optimization of this type of renewable energy device in the marine and coastal environment.

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

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