*2.3. Spatial and Temporal Sensitivity Study*

In order to obtain an adequate level of spatial and temporal discretization, a thorough sensitivity study was made in this work, based on the papers of Balduzzi et al. [13–15]. This was done to ensure the best compromise between accuracy and computation time. Specifically, it was suggested in [13–15] that, as the tip speed ratio was reduced, the grid refinement must be increased in order to adequately capture the higher vorticity gradients related to the increasing unstable operating conditions [13–15]. Furthermore, the temporal discretization must be reduced accordingly to the grid refinement so as the Courant number was below 10, thus providing the optimal error damping properties in the implicit scheme solver used in the present work. Therefore, due to the highly unstable condition related to the very low operating tip speed ratios of the micro rotor analyzed (between 0.025 and 0.29 in Table 1), a specific spatial and temporal sensitivity study was made. The study was carried out for each of the three turbulence models. In this way, the sensitivity of the turbulence models to the grid refinements was evaluated as well. Nine sensitivity tests for each rotor operating condition were made. Three meshing levels and three time steps were tested for each of the three turbulence models. The sensitivity was tested for the maximum (λ = 0.29) and the minimum (λ = 0.025) tip speed ratio in order to find the best compromise for all the intermediate simulations. A global amount of 54 simulations were carried out in this sensitivity study. In this way, an optimal spatial and temporal discretization level was found to be valid for all the simulated operating conditions.

The three grid refinement levels used are reported in Table 2, where M1 is the coarsest mesh and M3 is the finest one. To refine the mesh, the number of nodes on the airfoil was increased and the cell sizing of the rotating domain was reduced accordingly. The same cell sizing of the rotating ring was imposed to the inner circle to have a very fine discretization leading to an accurate transport of the vortices detached from the upwind blades. This constant fine mesh reduced the numerical diffusion on the upwind blade wakes, thus improving the accuracy of the downwind blade-vortex interactions. Exactly the same refinement was imposed to both the sliding interfaces so as to reduce the interpolation errors at the non-conformal sliding mesh. The boundary layer of the blades was discretized using quadrangular layers of elements. The quadrangular elements allowed for a more accurate resolution of the boundary layer compared to the triangular elements used for the rest of the domain.

**Table 2.** Grid independence study meshing characteristics.


The three meshes were developed with a different number of quadrangular layers and growth rates as reported in Table 2. The first layer height was fixed in such a way to guarantee a y<sup>+</sup> < 1 for all the simulations, as required by the turbulence models [29,31–33]. A sizing function, which limited the maximum dimension of the elements, was used in the stationary outer domain as well. The same growth rate as that used for the inflation layers was applied to the whole domain, thus allowing for a gradual mesh coarsening from the rotating interfaces to the boundaries.

In Figure 4, some details of the Mesh M2 are shown. It is worth remarking that Mesh M2 and Mesh M3 did satisfy the Grid Reduced Vorticity (GRV) criterion proposed by Balduzzi et al. [14]. GRV is a quantitative parameter for a qualitative a priori evaluation of the spatial discretization accuracy. GRV is defined as a vorticity normalized with respect to the local grid size, and thus gives an estimate of the velocity variation within a single element. Therefore, it represents the capability of the mesh itself of correctly computing the gradient flow features. The mesh sensitivity analysis and the evaluation of GRV are therefore strictly related: grid independent results are obtained when the discretization error becomes irrelevant, i.e., when GRV is sufficiently small [14]. In the present paper, it was verified that Mesh M2 and Mesh M3 had GRV < 1%, as recommended by [14], while Mesh M1 presented slightly larger values. Furthermore, all three grids satisfied the LES filter constraint imposed by the use of the DDES model, which required to have cell dimension lower than 1/30 the cord length [35].

**Figure 4.** Details of the rotor region (**a**), airfoil region (**b**) leading edge (**c**) and trailing edge (**d**) for mesh M2.

The temporal discretization must be defined so that all the relevant time scales of the flow were resolved. For this purpose, the time step dimension was chosen in such a way to try a wide set of Courant numbers. The Courant number was defined as:

$$\mathbf{C}\mathbf{o} = \mathbf{V} \,\,\Delta\mathbf{t} / \Delta\mathbf{x} \tag{1}$$

As reported in [14], for VAWTs, V is the peripheral velocity of the airfoil, Δx is the average distance between two cell centroids on the airfoil wall and Δt is the time step. The Courant Number expresses the ratio between the temporal time step (Δt) and the time required by a fluid particle, moving with V velocity, to be convected throughout a cell of dimension Δx.

Since V was the tangential velocity, in order to obtain similar Courant numbers for the different grids and rotational speeds, the time step must be varied together with the tip speed ratio. This obviously involved the angular step being kept constant by reducing the time step with increasing tip speed ratios. In Table 3, the Courant numbers, obtained for different combination of grids and time steps, are shown. Only the finest grids with the largest time step presented slightly high Courant numbers while, in all the other combinations, the Courant numbers were well below 10. This was in order to obtain the best error damping properties, as recommended in [14]. The adaption of the time step was made specifically for all the simulated operating conditions, thus ensuring the same angular step and approximately the same Courant number. In this way, the spatial and temporal discretization was suitable for the specific spatial and temporal scales in each simulation.

The results of the sensitivity study are shown in Figure 5. The charts present the trend of the time-averaged torque coefficient for the various combinations of grids, angular steps and turbulence models. The time-averaged torque coefficient is plotted as a function of the angular steps corresponding to the time steps in Table 3. The number of grid elements is reported on the horizontal axis. Figure 5a,b refers to the DDES model results at maximum (a) and minimum (b) tip speed ratio, respectively. Figure 5c,d refers to the RANS SST transition model while Figure 5e,f shows the RANS SST k-ω model results. The sensitivity analysis demonstrates that the grids M2 and M3 with angular steps Δθ = 0.035◦ and Δθ = 0.01◦ lead to results which are in very close proximity to each other, for both λmax and λmin. For the DDES turbulence model, the difference is approximately 1%. Therefore, mesh M2 with an angular step equal to Δθ = 0.035◦ was the best compromise for all the simulated range of tip speed ratio. Specifically, for all the other simulated operating conditions, the time step was adapted to obtain Δθ = 0.035◦.

**Figure 5.** *Cont.*

**Figure 5.** Average torque coefficient sensitivity analysis at λmax for DDES (**a**), SST Trans (**c**), SST kω (**e**) and at λmin for DDES (**b**), SST Trans (**d**), SST kω (**f**).


**Table 3.** Courant numbers for different grids and time steps at maximum and minimum λ.

It is already evident that both RANS turbulence models highly underestimate the average torque coefficient with respect to the DDES model. Even negative time-averaged torque coefficients are predicted at λmax. This would suggest that in highly unstable boundary layer conditions, the RANS turbulence models lead to wrong physical predictions of the flow-field. On the contrary, the DDES model results agreed very well with the experiments as demonstrated in the following sections.

#### *2.4. CFD SOLVER Settings*

The CFD solver setup is reported in Table 4. The ANSYS Fluent transient solver was used with a coupled algorithm for pressure-velocity coupling. The three aforementioned turbulence models were tested. Optimized local correlation parameters were used for the SST transition formulation both in URANS and in DDES. This optimization was carried out according to a previous work by the authors [36]. These local correlation parameters triggered and controlled the onset of the laminar to turbulent boundary layer transition within the transport equation for intermittency and momentum thickness Reynolds number. The number of iterations within each time step was set in order to ensure that all the residuals, within each sub-iteration, were well below 10–4. The turbulence boundary conditions were set according to wind tunnel data and literature suggestions [20]. The convergence criterion was to have a time-averaged torque coefficient variation lower than 0.1% between two consecutive revolutions [14]. This was achieved in about five to ten complete rotor revolutions. Once the convergence was reached, the data were sampled for two consecutive revolutions for the torque time-averaging. The simulated operating conditions are shown in Table 4. The simulations were carried out on a HP Z820 workstation with 24 available cores for parallel calculation and 128 GB of RAM memory. The calculation time per revolution was approximately 58 h with the SST k-ω model, 65 h with the Transition SST model and 71 h with the DDES model.


**Table 4.** Main CFD solver settings.
