A DMST-based tool to establish the best aspect ratio, solidity and rotational speed for tidal turbines in real sea conditions

A Double Multiple Stream Tube (DMST) routine to predict the performance of cross-flow hydrokinetic turbines in real environments is presented, along with a site assessment application to identify the most efficient turbine aspect ratio, solidity and configuration (single, or paired) for a selected area of the Northern Adriatic Sea. The peculiarity of this DMST tool is its 3D character, since it allows to reproduce the vertical distribution of the torque generated by the turbine. To this end, correlations for fluid dynamic phenomena, based on high-fidelity fully CFD simulations, were implemented. The marine circulation code SHYFEM is adopted to obtain velocity profiles for a half lunar cycle period. The sites with the highest mean kinetic power were identified. The DMST routine is equipped with an iterative process able to establish which rotational speed maximizes the power output. Indeed, a spatially non-uniform velocity profile requires to determine the flow velocity more suitable to obtain the rotational speed via Tip Speed Ratio (TSR) definition. To this end, the section of the blades working at optimal TSR varies from top to bottom, until the maximum power is reached. It works as a virtual Maximum Power Point Tracking system able to adapt the turbine operating conditions for the different turbine geometries, and for changes in flow conditions. The results show that for the case study, the performance curve shape influences the optimal TSR blade section: the latter is often located in the upper part of the turbine for the low solidity, whereas, for high solidity turbines, in the bottom half part.


Introduction
Issues such as sustainable development and environmental impact of human activities have raised increasing concern in recent decades. Indeed, in 2015 almost two hundred countries adopted 14 ambitious goals (United Nation Report 2020) concerning both the protection and well-being of the planet and the quality of human life. From an energetic point of view, the need to produce clean energy is crucial, as reported in goal number seven. In this context, tidal energy can be an adequate energy source for the near future. Indeed, it shows many advantages compared to other renewable resources. Since dependent from circulation, it can take advantage of forecasting tools, therefore, being predictable. Moreover, offshore farms have no visual impact on the system and their effect on the neighbouring ecosystem can be monitored and mitigated. The exploitation of sea areas can be limited using Cross Flow Turbines (CFTs): these devices are able to reach higher values of power production per unit of sea area, compared to the more known Horizontal Axis Turbine (HAT). Indeed, CFTs benefit of favourable fluid dynamic conditions if placed close to each other (Zanforlin and Nishino 2016;Zanforlin 2018). Furthermore, wake recovery in CFTs is faster than HATTs, so a shorter distance between devices is possible (Kinzel et al. 2012). All these aspects contribute to increase the power densities per unit of horizontal area, which is one of the most penalizing aspects of renewable energies, thus limiting the exploited sea area. The know-how of HATs is wide due to the broad background coming from the wind energy field, but CFTs show further advantages. Some components (such as the gearbox and the generator) can be placed above the sea level if floating platforms are adopted, with the advantage of more easily accessible components for maintenance. The platforms motion due to the waves does not significantly affect turbine efficiency, still good in skewed flow (Ferreira et al. 2006;Orlandi et al. 2015). Moreover, tidal energy is characterized by the reversal of the flow direction. CFTs have the ability to work independently of flow direction, whereas HATs need yaw system to set the rotor plane perpendicular to the incoming flow in order to harness the energy. This makes the system more complex and consequently less reliable. On the other hand, CFTs exhibit lower starting-torque and efficiency. These inconveniences can be mitigated by introducing an acceptable complexity into the system by means of a blade-pitching mechanism (Coiro et al. 2005;Chougule and Nielsen 2014;Zhang et al. 2015).
CFTs are characterized by unsteady phenomena during the rotation of the blades. The angle of attack changes cyclically, and dynamic stall can occur. Furthermore, the blades pass through a disturbed flow along the downstream path. This behaviour requires high spatial and temporal resolution for 2D or 3D CFD to properly reproduce the flow, too computationally expensive for geometrical optimization studies. Moreover, experimental solutions seem not to be efficient, especially if considering constrains on turbine prototype sizes and maximum Reynolds numbers achievable. Results obtained with Reynolds number ranges far from those in real operating conditions cannot be extended for general conclusions. In Bachant and Wosnik (2016) an experimental study on CFTs is presented: their results show that, for the considered turbine geometry, performance become Re-independent for Reynolds numbers of about 10 6 (based on the diameter). In Mason-Jones et al. (2012) a CFD study on HATs shows that a-dimensional parameters become independent from Re, when the latter reaches values higher than 5 × 10 5 .
In this scenario, analytical models seem to be the most efficient solution for geometrical optimization studies. Here we present a simple methodology that makes use of a Double Multiple Stream Tube model (DMST) to predict CFTs performance in order to evaluate different geometrical solutions. The DMST model is the equivalent of the BEM theory, originally born for HAT, applied on CFT. The BEM theory is the application of the momentum theory to the blade element theory. The latter evaluates forces acting on blade element using 2D considerations, starting from the relative velocity on the airfoil/hydrofoil and the angle of attack. Hence, the BEM theory supposes that changes in axial momentum (for a certain annulus or thin disc, respectively, for CFT and HAT) are solely due to forces acting on blade elements. Various modifications and corrections can be taken into consideration to improve the model with 3D effects. The first application of the momentum theory to CFTs is attributed to Templin (1974), whereas the adoption of multiple stream tubes is first addressable to Strickland (1975). In literature we can find works, which make use of more or less elaborate DMST models. The complexity lies in the number and quality of the sub-models embedded to reproduce fluid dynamic phenomena. In Sanvito et al (2021) we find a recent example of DMST model application, equipped with various submodels taking into account: a variable velocity profile due to the atmospheric boundary layer, skewed flow and/ or rotor tilt, the dynamic stall, the Glauert correlation for high loaded conditions, the flow curvature correction, the tip losses, the supporting strut and pole effects.
In our work the DMST model was implemented as a MATLAB routine. The model has on the inside several submodels to reproduce unsteady phenomena, such as: dynamic stall, tips losses and deviation of streamlines approaching the turbine. We wanted to insert the turbine into real flow conditions: to this end we used, as input for the DMST routine, flow data obtained with the open source 3D circulation model SHYFEM developed by CNR-ISMAR in Umgiesser et al. (2014Umgiesser et al. ( , 2018, open source code freely available on the web pages http://www.ismar.cnr.it/shyfem and https:// github.com/SHYFEM-model). The purpose of this work is to provide a tool capable of establishing very quickly some geometrical and operational parameters of a CFT. This must be done before proceeding with more computationally expensive, and more accurate fluid dynamic simulations. This allows the designer to decide which turbine is most suitable for a particular site, setting some variables by considering the real environment in which the turbine will operate. The structure of this article is here summarized: in section two we describe the methodology adopted and the used tools; section three describes the setup of the tools and validation; section four shows the results, which are discussed in section five; section six summarizes some conclusions.

Turbine geometries
In this study, after a preliminary identification of suitable sites to place turbines, we evaluated different geometries of H-Darrieus turbine with 3 blades and based on NACA0018 hydrofoil. We considered two values for the frontal area 25 m 2 and 50 m 2 , to evaluate which is the best configuration among the following: a couple of small devices or a single device with the same overall frontal area. We assigned four values for the aspect ratio (AR, which is the ratio between the height H and diameter D of the turbine) for each area. Two values for the solidity are evaluated, a low value, σ L , equal to 0.0637 and a high value, σ H , equal to 0.159. Th solidity σ is defined as where N B is the number of blades and c is the chord. All the geometric configurations considered, for each solidity, are summarized in Table 1.

DMST model
For this study we used the DMST model developed at the University of Pisa , further information are available in Appendix 1. The DMST model is based on a simplified approach, which does not take into account some peculiar phenomena that characterize the CFTs behaviour. Hence, special sub-models are necessary to reproduce the effect of these phenomena. The dynamic stall is due to the periodical variation of the angle of attack during blades' rotation that characterized CFTs. This can cause an increase in the angle of attack and the born of the Leading Edge Vortex (LEV), which enhances suction on the hydrofoil. Therefore, the maximum lift coefficient achieved with dynamic stall conditions, is higher than in static conditions. A further increase in the angle of attack α makes the LEV move toward the trailing edge, causing a quick drop in performance. Decreasing the angle of attack a hysteresis is observable in the C L − α curve. This behaviour is reproduced in the MATLAB routine through the dynamic stall sub-model developed at the University of Pisa and described in Rocchio et al. (2020).
An important fluid dynamic loss is related to the helical vortex shed at the blade tip. This phenomenon, which consists in a flow circulation over the tip, occurs at the ends of any lifting body (airplane wings, turbine and propeller blades), since, due to a pressure difference between the two sides of the body, the flow tends to pass, where the fluid dynamic resistance is lower. This causes efficiency to drop at tips, since part of the fluid is not available to produce lift force. To reproduce the effects of tip losses in the spanwise direction we used the k(μ) factor, which represents the ratio between the local power coefficient along the blade, C P (μ), and the power coefficient of the mid plane of the turbine (where μ 0). Defining C P (μ) as where P is the generated power, ρ the fluid density and z the vertical thickness of the considered turbine layer. The last phenomenon, which must be considered is the streamlines deviation. The flow tends to pass sideways to the turbine, because the latter represents an obstacle (Goward Brown et al. 2017). To implement this behaviour we used the f factor, which is the arctangent of flow deviation angle with respect to the undisturbed direction. We modified velocity approaching the turbine considering W d , where U d represents the x component, whereas V d the y component. Relationships are the following: where U ∞ is the undisturbed flow velocity, a the axial induction factor defined as a W d /U ∞ , and θ is the azimuthal position.
In Ansys-Fluent CFD simulations, used for the validation and calibration step, we adopted a hydrofoil cambered on the circumferential trajectory of the blade, as shown in Fig. 1. This hydrofoil shape combined with flow curvature effects (Migliore et al. 1980) behaves as a symmetrical hydrofoil in straight flow. Therefore, we run the DMST model with curvature sub-model inactive, and used NACA0018 airfoil polars reported in Sheldahl and Klimas (1981).

Sites assessment
We used flow data, coming from the SHYFEM code, as input for the DMST routine. The flow data refer to the Northern Adriatic Sea (latitude from 44.5 to 45.5 and longitude from 12 to 13), and cover a period from 7 to 21th of February 2014 with 337 hourly outputs. The grid was characterized by a spatial discretization of about 600 m for both horizontal directions. In the vertical direction we had 17 levels with variable thickness: from 1 m at the top, to 2 m moving towards the bottom.
To identify adequate sites to place turbine, we first needed to discriminate locations based on depth. We considered to leave 2 m from turbines tips and the free surface and also at least 2 m from the turbine bottom and the seabed. An appropriate distance from the free surface is necessary to have a proper wake development, and to benefit from the bypass flow (Birjandi et al. 2013;Kolekar and Banerjee 2015). Indeed, the presence of the turbine causes the flow to deviate, and consequently a by-pass flow is observable at the top and the bottom of the device. Hence, at the same time a proper distance from the seabed is indispensable to mitigate shear stress. This phenomenon must be limited, because erosion is one of the most evident environmental impact of such energy exploitation mechanism (Ramírez-Mendoza et al. 2020; Gillibrand et al. 2016).
An energetic criterion was then adopted: we searched the maximum available power over the entire area, and accepted sites characterized by at least the 50% of the maximum power found.
The DMST routine represents a preliminary study able to identify suitable sites, geometries and, as it will be explained in the following paragraph, the optimal rotational speed of the device in single or paired configurations. However, it does not allow to predict the fluid dynamics interactions between nearby turbines nor the effects of the turbine on the environment (momentum reduction, turbulence). For this reason, this approach can be used as the first phase of a more complex site assessment study, aimed at the identification of the site and the planning of the layout of a tidal farm. The subsequent phases should be carried out through a marine circulation code within which a turbine model routine is embedded and used in "two-way" mode, as already done in Pucci et al. (2020), Thiébot et al. (2020). Therefore, in this first phase it is indispensable to minimize computation time. To this end, each of the ten chosen sites is characterized with a unique velocity profile. We calculated the average power in time P av for each site and we extracted a velocity profile U av representative of the average power as shown below: We considered 337 h, which represents half of a lunar cycle. Furthermore, this method, which requires only one calculation to be carried out with the U av velocity profile, allows to predict the energy production with limited errors compared to the hour-by-hour calculation (this aspect will be discuss in section four and five).

Optimal rotation speed
To evaluate which rotational speed assign to the turbine, we implemented an iterative process inside the DMST routine. It consists of two phases, the first of which begins by dividing the turbine into 51 vertical layers, of variable thickness depending on turbine height (the thickness z goes from a minimum of 0.08 m to a maximum of 0.19 m for the different turbine geometries). Then, the code assigns at each layer the optimal Tip Speed Ratio (TSR) value (on the basis of the C P −α curve of Fig. 2), and evaluates ω from the following equation: The value of the chosen ω, among the 51 plausible values, is the one which allows to maximize the power of the turbine. This concludes the first phase. Then, to assure we have identified the best omega value, we check that it does not coincide with the value assigned to the bottom or top layer of the turbine. If, on the other hand, it coincides with an end layer, we proceed with the second phase, which consists in iteratively increasing or decreasing omega (depending on whether the optimal omega coming from the first phase has been found for the top or for the bottom layer), using a small delta-omega equal to 5% of the optimal TSR until the turbine power continues to increase.
The advantage of this procedure is that the first phase is extremely rapid, and that in the vast majority of cases it is sufficient (in the case study described in the following, it was enough to resort to the second phase in less than 10% of the simulations).
The choice of setting 51 layers is in analogy with the number of layers considered in the DMST routine. Indeed, in Deluca (2018) a sensitivity analysis for the vertical discretization was done. 51 layers was a good compromise between required calculation time and results accuracy. Indeed, considering the power coefficient as output, an error of 2.5% was found with respect to results obtained with 2001 layers. The height was variable also for Deluca (2018): a turbine with radius equal to 3.16228 m and three different AR (0.9, 1.35 and 1.8) were considered. However, it must be considered that different turbine geometries means different diameters and consequently different Reynolds numbers. Variations in Reynolds numbers could significantly affect the value of the actual best value of TSR for each layer. This is particularly true for low solidity turbines, as shown in Fig. 2 (obtained without applying the correction for tip losses, in order to emphasize only the effect of Re).
From the results of Fig. 2 we can deduce that, in case of low solidity, the optimal TSR is higher than 2.65, and that for the turbines with smaller diameters (lower Re values) the optimal TSR could even be higher than the last value taken into consideration (that was 2.85).
What we have done in the search of ω is a kind of Maximum Power Point Tracking (MPPT) strategy. MPPT usually refers to an electrical system capable of adapting turbine operating conditions to varying wind or flow velocity in order to maximize the power production. Some of this electrical equipment work i with TSRstrategies (Carriveau 2011). Wind speed can be measured with anemometer or can be estimated with particular algorithm. This method is useful for wind turbine, where the rotor region is far from the boundary layer, and the velocity profile can be considered flat, equal to the undisturbed flow speed at each vertical position. However, it can be inadequate in other applications, such as wind turbines intended for city contests, or tidal turbine, which are characterized by turbulent flow (Carpman 2010). For tidal turbine, which often exploit the energy resource in shallow water, we cannot overlook the shape of the velocity profile. Therefore, it becomes important for ω evaluation, the vertical height, where we assign the TSR value and consequently the height of the reference undisturbed flow velocity.

Sub-models setup and validation
In this section, we will show the validation process used for each sub-model embedded in the DMST model. It should be noted that the use, where possible, of 2D CFD simulations for the validation, is driven by: the need to isolate the effect of a particular phenomenon, and using 2D simulations threedimensional effects are undoubtedly eliminated, and the need to save time, since the aim of our research is to provide a quick tool to establish turbine parameters.

Dynamic stall sub-model
A tuning process for each solidity was necessary to adapt the model to the chosen hydrofoil. In Fig. 3, we can compare performance obtained from our CFD 2D Ansys Fluent simulations (red in the graphs) and DMST tuning (blue in the graphs), both carried out considering U ∞ 1.75 m/s and R 3 m. Good agreements emerge between the efficiency of a  Pucci et al. (2021) single blade during a rotation ( Fig. 3a, b, respectively, for σ L and σ H ). Similar trends are also noted for the total C P − TSR curve (Fig. 3c). The grids adopted for the 2D CFD simulations of Fig. 3 are unstructured and only made of quad-elements; they were generated with the software Ansys-Icem by means of the technique "patch-dependent". The overall validation of the 2D CFD model has been described in a previous article (Pucci et al. 2020) in comparison to the experimental data by Bravo et al. (2007).

Fluid dynamic tip losses sub-model
The values of the correction factor k, which reproduces tip losses, were extracted from 3D CFD simulations carried out with the commercial software Ansys-Fluent v.19. Before analysing the behaviour of k for the turbines of interest, an essential description of the numerical model and its setup is here reported. Ansys-ICEM was used to generate multiblock structured 3D grids, with the addition of O-grids to thicken the distribution of cells in the areas of greatest interest and at the same time to improve their quality. Two grid levels are used to simulate the blade rotation via the sliding mesh method: a fixed sub-grid with the outer dimensions of the flow domain and a rotating sub-grid including the turbine blades. All around the blades the grid is very fine to make sure that y + at the walls remains below 0.4, following the work of Maître et al. (2013), who analysed the effect of y + , realizing that averaged y + > 1 causes a pressure drag overestimation in turbines exposed to significant flow separation, as happens for high solidity hydrokinetic turbines. Some grid details are shown in Fig. 4 (which refers to a turbine with σ H and AR 2).
The moving and the fixed domains are joined together by means of an interface boundary condition (BC). The other BCs are: velocity inlet for the inlet face; pressure outlet for the exit; symmetry for the top and the bottom boundaries. The whole turbine was simulated in case of low AR, and to save computation time, in case of high AR only half turbine was simulated, since blade are very long. In this way, the fixed domain has~2.6 × 10 6 hexahedral cells, while the rotating domain has~12 × 10 6 and 9.6 × 10 6 cells in case of low and high AR, respectively. To model the turbulence, the k-ω SST (Shear Stress Transport) was adopted Fig. 3 Comparison between CFD and DMST "Stall Tuning" efficiency of a single blade during a rotation at the optimal TSR for σ L and σ H , respectively (a, b). C P − TSR CFD and DMST curve for both solidities (c) (Menter 1994;Wilcox 2008); this model is widely used in the simulation of wind and tidal cross-flow turbines, since it is considered appropriate in case of flow characterized by strong adverse pressure as happens in cross-flow turbines, especially when operating at low TSR. The algorithm for the velocity-pressure coupling is SIMPLEC (Jang et al. 1986). About the spatial discretization scheme, the Least Squares Cell-Based (LSCB) is set for gradient; pressure interpolation, turbulent kinetic energy and specific dissipation rate formulations are based on second-order schemes. Temporal discretization is also based on a second-order implicit method. The convergence criteria for each time-step was 1 × 10 −4 for the residuals of continuity, velocity components, turbulence kinetic energy and specific dissipation rate. When using a sliding mesh, to obtain a satisfactory numeric convergence, the time-step should not be greater than the time required by the mobile interface to advance one cell. Thus, in order to consider the smallest cell at the interface, a time-step corresponding to 0.5°of revolution is adopted; this value is also in accordance with Balduzzi et al. (2016). The reader can find details about the validation of our 3D CFD model in previous works: in Zanforlin and Deluca (2018) some C P sensitivity analyses to grid refinement and number of revolutions are discussed, together with a comparison of the C P − TSR curve with the experimental data by Maitre et al. (2013), in Pucci et al. (2020) the model is validated versus experimental profiles of the velocity measured in the turbine wake by Vergaerde et al. (2020).
Four simulations were performed assuming two solidities (σ L and σ H ) and two aspect ratios (AR 2/3 and 2), while the turbine diameter was fixed at 6 m. Given these ARs and coherently with Zanforlin (2018), to achieve steady and then reliable results, 15 and 10 revolutions have been done, respectively, for the high and the low AR cases. The free stream (or undisturbed) velocity is 1.75 m/s. Since the lift force and hence the torque generation strongly depend on the pressure difference between the pressure and suction sides of the blade, it is useful to analyse the spanwise distribution of the static pressure especially when the conditions of attack angle and incoming flow velocity are more favourable, as it happens in upwind when θ is about 100°. Figure 5 shows the pressure distribution on the blade surfaces for the low σ and low AR turbine, when blade_1, 2 and 3 have angular position θ of 100°, 220°and 340°.
Focusing the attention on blade_1, it can be observed that the pressure on the blade pressure-side is uniform for about the 60% of the length, but it progressively drops approaching the tips making impossible the power extraction. This happens, because the flow approaching the turbine deviates in the spanwise (i.e., the vertical) direction, due to the turbine blockage, the finiteness of the blade length and the tendency of the flow to move towards lower pressure regions.
In Fig. 6, the z-components of the flow velocity are depicted on a vertical plane parallel to the free stream and cutting blade_1, for the low σ turbine. The positive and negative values (red and blue spots) calculated upstream and downstream the blade tip suggest the helical motion of the flow approaching the final part of the blade span. The flow "leakage" around the tip decreases the pressure difference between the suction and pressure sides, reducing lift. Moreover, tip vortices imply a localized huge pressure drag increase. As a result, performance drastically drops at the blade tip. Moreover, it should be noted that the phenomenon is not confined at the tip but propagates along the span in relation to the vortex strength, which in turn depends on σ and TSR. Since the turbines of Fig. 6 have the same σ and TSR, the blade absolute length affected by tip losses is expected to be quite similar, as proved by the vertical extension, in the blade proximity, of the most significant z-velocities (for instance, greater of 0.2 m/s, corresponding to yellow-red). However, by comparing these lengths to the actual blade lengths, it can be seen that they correspond to 25% of span in case of high AR and 40% in case of low AR; therefore, the relative incidence of tip losses is higher the lower the AR.
The small difference (in absolute terms) between the tip-leakage flow-rates of two turbines is also confirmed by Fig. 7a), showing the z-velocities on the horizontal plane passing through the tips. Indeed, the extensions of the red and blue stripes in correspondence to the passage of blade_1 Fig. 6 Flow z-velocity on a vertical plane parallel to the free stream and cutting the blade_1 (θ 100°), in case of low σ turbine (high AR on the left, low AR on the right). Only half turbine is shown in upwind can be considered as indicative of the tip vortex strength, and it can be observed that this extension is not very dissimilar for the two cases depicted.
To exemplify the effect of the turbine solidity on tip losses, we can compare Fig. 7a, b, the latter showing the z-velocity field in case of the high-σ and high-AR turbine. The extension of the red-blue stripes, and therefore, the tip-leakage flow-rate, is significantly higher for high σ turbine of Fig. 7b.
The spanwise distribution of the dimensionless C P , k(μ), predicted by CFD for the turbines working at their optimal TSR (1.4 and 2.75 for the high and low σ , respectively), are reported in Fig. 8. Consistently with what has been observed up to this moment, it can be seen that the performance deterioration due to tip losses increases (occurring closer to the blade mid-plane) as σ increases and AR decreases.
The k(μ) values of Fig. 8 were used to multiply the C P (μ 0) achieved by means of the DMST, regardless of the particular TSR and size of the turbine. For AR different from 2/3 and 2, the k(μ) coefficients have been assigned by linear interpolation.

Streamlines deviation sub-model
To reproduce the effect of streamlines deviation we adopted the already described f factor, obtained from 2D CFD simulations with Ansys Fluent. As can be seen from Fig. 9, deviation is more evident with higher TSR (maintaining constant the solidity), because it means higher rotational speed, and consequently a turbine which is less permeable to the flow. At the same time, higher solidity at almost the same TSR generates greater deviation (compare Fig. 9 right column and left column). We extrapolated the f factor for discretized θ and TSR values from the CFD simulations, and via interpolation, we obtained the deviation factor for each azimuthal position and each TSR in the DMST model.

Results
The site assessment procedure consisted in deleting those areas (black coloured regions in the bathymetric map of Fig. 10a which do not allow to respect bathymetric constraints (not even in the case of the turbine with the shortest blades).
From the energetical evaluations, we can see that the most of the energy content is located in front of the delta Po River, and at the three inlets connecting the Venice lagoon to the Adriatic Sea (Fig. 10b shows power available on the surface layer). We finally identified 10 suitable sites (Fig. 10c).
Then, we applied the DMST tool to the 10 sites. Below we summarize some results coming from the iterative optimization process. In Figs. 11 and 12, turbines with frontal area 25 m 2 and 50 m 2 are, respectively, considered. For both figure, we compare site 5, 6, 8 and 9 (respectively, in a, b, c and d), which are the sites that allow to test the majority of turbine geometries. The graphs show the vertical extension of the turbine (dashed lines) by varying the AR, and the reference height for the available power along blade length H PR (line with triangles in Figs. 11 and 12). The latter variable is defined as the height weighted on available power. The marine area considered in the study, is characterized by velocity profiles such as to make H P R fall into the upper part of the turbine, and H PR changes with an almost linear trend as the AR varies. For σ L the section of blade working at optimal TSR (green line with squares) lays always above H PR . Whereas, for σ H the optimal section (red line with rhombus) lays often under the H PR line. It is worth noting that we refer to optimal section, as the one showing the maximum local C P along blade.
In Fig. 13, we can compare the power coefficient along blade length (for site 7) both with and without considering tip losses. The figure shows the lower solidity on the left side, Fig. 7 a Flow z-velocity on the horizontal plane passing through the blade tips, in cases of low-σ turbine (high AR on the left, low AR on the right); b flow z-velocity on the horizontal plane passing through the blade tips, for the case of high-σ and high-AR turbine  whereas the higher solidity on the right side. Obviously, tip losses have major effects on lower ARs. However, from the diagram of Fig. 13 "with tip losses" it is not easy to assess at a glance which AR has the highest average C P along the blade. However, it is worth noting that the C P is not sufficient to determine which geometry achieves the maximum power output per unit area. Indeed, the power coefficient indicates the performance with respect to the available energy resource. The extracted power also depends on the available power, which is greater in surface flow layers everywhere in the sea area. Thus, geometries with short blades, i.e., low AR, are favoured. Probably, this fact justifies higher specific powers found for low AR, as shown in Table 2.
As previously anticipated, we characterized each site with a velocity profile representative of the average power content. This choice also makes it possible to evaluate the energy production during the entire period of work. Indeed, comparing the energy production obtained with the unique calculation and the energy production with hourly calculation, differences are very limited. In the hourly calculation, we assumed that the turbine does not work when the 2/3 of the blade are reached by a flow slower than the cut in velocity (we considered a value of 0.5 m/s). In Table 3, we summarize energy productions obtained both with the unique calculation (AV subscript) and with hourly calculation for random sites and turbine configurations.
The averaged calculation always overestimates the hourly calculation, but differences are limited to 4-5%. The adoption of a time averaged velocity profile would lead to sensible differences in results. In Fig. 14, we can compare quantities based on time averaged velocity profile (black curves) and quantities based on time averaged available power (green curves). The dashed lines represent generated power along a blade, whereas the dot lines represent the velocity profile. Figure 14 is obtained considering site 7, Area 50 m 2 and AR 1.11. We can see the huge difference in behaviour given by the two criterions. Since the great accordance in energy productions results from Table 3, time averaged velocity profile would surely lead to errors in characterizing sites.

Discussion
From Figs. 11 and 12, we have highlighted trends for the section working at optimal conditions (i.e., where C P (μ) is the highest). The behaviour is mostly due to the tuning curve C P − TSR. Considering Fig. 3, we can see that the tuning curve for σ H is almost symmetrical around the optimal TSR value 1.4. Indeed, for variation limited to ± 15% (i.e., TSR 1.2-1.6) the curve is almost specular around the maximum value, but always with higher efficiency for lower TSR. For further increases, ± 30% (i.e., TSR 0.98-1.82) the drop in performance is more relevant in the left side of the curve. Otherwise, the σ L tuning curve is not so specular around the optimal TSR value: there is a great drop in efficiency as the TSR decreases. For the same above mentioned percentage variations from the reference value of 2.75, we have, respectively, TSR range of 2.34-3.16 and 1.92-3.57. In both cases we have higher efficiency for high TSR. Indeed, looking at Figs. 11 and 12, the optimal section for σ L tends to rise as the AR increases. Greater ARs lead to increase in blade height, and to avoid stall condition in the upper part of the blade, the optimal section moves towards tip, getting farther away from H PR . This is due to limits in the tuning process, indeed from literature knowledge the expected behaviour is a quicker change in performance for high solidity turbines instead of low solidity (Rezaeiha et al. 2018). Hence, considering this limitation in the C P − TSR curve, we can explain results for Fig. 12 Comparison between the height of the section working at optimal TSR for σ L (green line with squares) and σ H (red line with rhombus) with respect to the reference height for available power (black line with triangles). The considered turbines have a frontal area of 50 m 2 and variable AR. The graphs show results for site 5, 6, 8 and 9, respectively, in (a-d) the optimal section as follows: velocity profiles considered in this work are representative of real environmental conditions, so they are not uniform profile, but there is a vertical gradient. For this reason changes in speed along blade length causes changes in TSR. The latter increases with depth. Hence, for σ L the optimal section lays always in the half upper part of the turbine, to privilege a TSR range shifted towards high TSR, where higher performances are observed (Fig. 15).
Whereas, for σ H for limited TSR percentage range variation, the DMST routine privileges lower TSR values by assigning the optimal TSR value often in the lower half part of the turbine. Therefore, the trend for the optimal TSR section is primarily guided by the C P − TSR (Figs. 11 and 12).
Considering now the turbine geometry, we have higher power production with lower solidities (as shown in Table 2). It is well known that, both for HATs and for CFTs, the turbine solidity influences the performance and the TSR value at which the peak C P is obtained. Indeed, the higher the solidity, the lower the TSR of optimal operation (Jamieson 2011). However, the effect of solidity on the peak value of C P is not monotonous but, for given values of Re and the number of blades, there is a solidity that maximizes the C P . For solidities lower than this value the performance deteriorates because of the power losses due to viscous friction, which are roughly proportional to the cube of TSR. On the other hand, for higher solidities the performance is penalized by some phenomena that typically affect CFTs: the onset of stall, since a lower optimal TSR and, therefore, a lower rotation speed implies an increase in the attack angle; the interference between the blade and the wake released by the preceding blade; the fact that the more the chord is extended the more the blade generate torque almost exclusively during the upwind path. Our results, show C P slightly lower in the case of higher solidity. This is in accordance with literature, since experimental and numerical studies proved that for three bladed CFTs the Fig. 13 Single blade efficiency along the vertical direction for site 7 turbines. On the left side we consider σ L , whereas on the right σ H . At the top of the figure we have blade efficiency without considering tip losses, whereas on the bottom we have included tip losses effect power coefficient peak is reached for solidity equal to 7.1% (Hyun et al. 2010), 7.2% (Bedon et al. 2012), 9.5% (Rezaeiha et al. 2018) and 9.8% (Sagharichi et al. 2018) that are relatively low values. The only exception is Du et al. (2019) work, where the optimal sigma was 12.9%. Therefore, in agreement with literature, it is reasonable to conclude that for our work too the peak of power coefficient arises at σ values relatively low, and closer to σ L (6.37%), instead of σ H (15.9%). This can be the explanation of better performance obtained with σ L .
From Table 2, we notice that the higher power productions occur for the lowest AR value and frontal area 25 m 2 . This means that tips losses were not decisive in this case study, major advantages with low AR are due to: greater available kinetic power, lower TSR variation along blade length and higher Re (for a fixed value of frontal area Re increases with decreasing AR).
In Brusca et al. (2014) a similar approach is presented, carried out with a Multiple Stream Tube model: they use the same airfoil NACA0018. An iterative process is at the basis of their study too. The authors conclude that a lower AR leads to increase in diameter and consequently in Reynolds number. The enhanced Reynolds number causes the performance to rise. We can share this conclusions for our work too. Furthermore, we can highlight also differences caused by the frontal area value. Indeed, Brusca et al. (2014), have considered a uniform wind velocity profile, which is a reasonable assumption in the wind field. Therefore, their conclusions are mostly oriented in the influence of Reynolds number without considering the effect of the frontal area. Differently, in our study having considered realistic flow velocity profile, characterized by vertical gradient, the influence of the frontal area is evident. Indeed, even if the turbine "Area 50 m 2 AR 0.67" reaches the highest Reynolds number, this can not overcome the penalty of exploiting a slower flow. Comparing the same AR value for the two different frontal area, the lower frontal area reaches always higher power production, because these devices exploit a flow with higher kinetic potential. For this reason, with the same AR, it is preferable to install a couple of counter-rotating small devices, instead of a single device with the same overall frontal area. Furthermore, a pair of devices can benefits of other fluid dynamics mechanisms able to enhance efficiency (Zanforlin and Nishino 2016;Zanforlin 2018;Dabiri 2011). Moreover, two devices are capable to naturally mitigate back-torque on supports, especially if a floating platform is adopted (Vergaerde et al. 2020). In our work we have taken into consideration tip losses effects, whereas this effect is not mentioned in Brusca et al. (2014).
In Hunt et al. (2020) we find a criticism to one aspect of the Brusca et al. (2014) research. They underline that in the study, variation in AR are reached by varying diameter and also solidity, changing consequently the chord to radius ratio. They assert that as chord to radius and solidity decrease, the optimal TSR increases causing the Reynolds number to locally increase, and so also performance. Moreover, changes in chord to radius values lead to flow curvature phenomenon (Migliore et al. 1980), which requires empirical corrections. In our study we have considered two "turbines families", and for each of them the solidity is kept fixed with variable AR, and therefore, it is also fixed the chord to radius ratio. Furthermore, for the tuning process we have used 2D and 3D CFD simulations carried out with a hydrofoil profile cambered on the circumferential trajectory of the blade in order to delete the curvature effect (as already shown in Fig. 1). For this reason we have neglected the virtual camber effect in the DMST model. Thus, we can reasonably conclude that for each solidity "family", performance are guided by the Reynolds number (beside the C P − TSR curve shape in the TSR range available on blade length), whereas for the power production the influence of the quality of the energy resource available is predominant. Hunt et al. (2020), have carried out experiments maintaining the blockage, the Reynolds numbers and the Froude numbers at a fixed value. They also use the same NACA0018 hydrofoil. They notice that no significant changes in efficiency occur by varying AR. However, they evaluate two different Reynolds number values, i.e., 4.27 × 10 4 and 2.03 × 10 4 , both of which are very low compared to those under real operating conditions. The above mentioned Reynolds number was calculated using as characteristic velocity, the undisturbed flow velocity, instead of using the tangential blade velocity. However, even for TSR greater than one, using the tangential blade velocity would not lead to significant enhancement in Reynolds number. Therefore, it is not possible to extend results, obtained for low Reynolds number, to realistic operating conditions.
In Li et al. (2013), we find a study carried out with the panel method: this method is based on incompressible and potential flow assumptions together with a free vortex wake, and it is able to reproduce forces acting on surfaces. They evaluate performance of VAT for different AR and solidity. They keep the diameter constant, and vary the height and chord. What emerges is that, keeping the solidity set to a fixed value, increases in AR lead to higher optimal TSRs. This is in accordance with our preliminary analysis shown in Fig. 2a. Furthermore, they found that greater AR lead to higher peak of C P . This aspect is apparently in contrast with our conclusions, but it can be explained considering that their main stream velocity is uniform, they do not consider a variable velocity profile along the vertical direction. In this way, devices with higher AR benefit from a lower impact of tip losses on performance.

Conclusions
In this study we use a DMST model implemented in a MAT-LAB routine, and equipped with some sub-models in order to reproduce dynamic and 3D phenomena. The peculiarity of our work is in the sub-models validation: for the first time, it makes use of high fidelity 2D and 3D CFD simulations instead of complex and expensive experimental tests. The input flow data for the DMST routine, are obtained with the SHYFEM marine code, thus to reproduce realistic environmental flow conditions. Moreover, by adopting realistic velocity profile with vertical gradients, the routine is adequate to evaluate the rotational speed for turbine in turbulent flow. This is the case of tidal turbines, but also wind turbines in city contexts. The optimal rotational speed is found by means of an iterative process, which assigns the optimal TSR to various blade sections, progressively from top to bottom. The chosen speed is the one which maximizes power production. The results show that for the area of study, the optimal section of the low solidity turbines is located in the upper part of the turbine, thus favouring high TSRs, where greater efficiency occur. Whereas, for high solidity turbines, the optimal section lies often in the bottom half part of the turbine. This is mostly due to the TSR variation along blade, and hence is linked with the C P −TSR curve shape. The combined use of DMST routine and SHYFEM code, allows quick sites assessment, turbine geometric evaluation, turbine working conditions, and forecasts in energy production. Hence, it is suitable for preliminary studies, in order to examine some parameters before proceeding with a broader research, which implies computationally more expensive simulations.
Funding Open access funding provided by Università di Pisa within the CRUI-CARE Agreement.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecomm ons.org/licenses/by/4.0/.

DMST fundamentals
The DMST model consist in an iterative process to assign the a value, which is the axial induction factor. The process starts by assuming an a value. By using it, we calculate the relative velocity to the blade, the angle of attack, drag, lift and extrapolate the analytical expression for the thrust coefficient C X . At the same time C X is evaluated using the theoretical formulation coming from the actuator disk approach (Burton et al. 2011). If these two results match, the a value hypothesized is correct, otherwise we repeat the calculation (Fig. 16). The turbine is virtually divided in several horizontal planes of thickness z, and each plane is divided into multiple streamtube which are double, an upstream tube and a downstream tube as schematically shown in Fig. 17a. Double stream-tubes allow to identify an induction factor for each half tube: we will refer with subscript 1 for the upstream tube and 2 for the downstream tube.
Assuming the a values we can write the balance equations: Combining the above equations with the a values, defined as a 1 U 1 U ∞ and a 2 U 2 U e , we obtain the following velocities relationships: Defining thrust coefficients, starting from the total thrust force F X , we obtain: . 17 a Generical horizontal plane of the turbine for DMST reference; b schematic representation of velocity involved in DMST calculations which implies a theoretical formulation for C X , i (with i 1, 2) equal to However, this formulation shows few correspondence with empirical data for a i < 0.6. For this reason, the following relationship provided by Spera (2009) is used: To obtain the analytical formulation for the thrust coefficient we consider a generic horizontal plane of the turbine located at z quote (Fig. 17b). The following relationship between velocities can be noticed: U n (θ ) U 1 sin(θ ) ( 2 0 ) U t (θ ) ω R + U 1 cos(θ ) (21) From the local angle of attack and the Reynolds number, respectively, α tan −1 (U n (θ )/U t (θ )) and R e (θ ) cW (θ) ν , we can evaluated the lift coefficient C L and the drag coefficient C D . By further geometrical consideration it follows that C t C L sin(α) − C D cos(α) C n C L cos(α) + C D sin(α) C X C t cos(θ ) − C n sin(θ ) C Y C t sin(θ ) + C n cos(θ ) We compare C X , i obtained from Eqs. 18 and 25 inside the iterative process.