Development of a Computational Model for Investigation of and Oscillating Water Column Device with a Savonius Turbine

: The present work aims to develop a computational model investigating turbulent ﬂows in a problem that simulates an oscillating water column device (OWC) considering a Savonius turbine in the air duct region. Incompressible, two-dimensional, unsteady, and turbulent ﬂows were considered for three different conﬁgurations: (1) free turbine inserted in a long and large channel for veriﬁcation/validation of the model, (2) an enclosure domain that mimics an OWC device with a constant velocity at its inlet, and (3) the same domain as that in Case 2 with sinusoidal velocity imposed at the inlet. A dynamic rotational mesh in the turbine region was imposed. Time-averaged equations of the conservation of mass and balance of momentum with the k – ω Shear Stress Transport (SST) model for turbulence closure were solved with the ﬁnite volume method. The developed model led to promising results, predicting similar time–spatial-averaged power coefﬁcients ( C P ) as those obtained in the literature for different magnitudes of the tip speed ratio (0.75 ≤ λ ≤ 2.00). The simulation of the enclosure domain increased C P for all studied values of λ in comparison with a free turbine ( Case 1 ). The imposition of sinusoidal velocity ( Case 3 ) led to a similar performance as that obtained for constant velocity ( Case 2 ).


Introduction
The energy demand will increase by more than 1.0% per year up to 2040, increasing gases emission rates [1]. Moreover, the costs of commodities for energy generation from fossil fuels have increased significantly, leading to economic difficulties, risks associated with energy security, and geopolitical conflicts around the world [2]. Considering this scenario, there is a growing search for a better comprehension of the development of technologies and the use of devices and economic impacts of different renewable sources of energy such as wind, solar, geothermal, and ocean energy [2][3][4][5][6][7][8][9][10][11].
One of the important sources of renewable energy with high potential, but not frequently explored worldwide, is the conversion of ocean energy into electricity [12][13][14][15], i.e., wave energy conversion. Despite several signs of progress in technological development, there is no dominant main operational principle. Several devices have been proposed and investigated based on various ways to convert wave energy such as point absorbers, attenuators, oscillating surge converters, overtopping, submerged plates, and oscillating water columns (OWCs) [12][13][14]. The efficiency and survivability of the wave energy converters (WECs) are important issues to make them more competitive and viable [16]. In this context, the OWC device has advantages as its simplicity and maintenance, e.g., the moving parts are located outside of the water, increasing the lifetime material of the power take-off (PTO) system, and the structures of buildings are robust [12][13][14]. Therefore, several studies and prototypes using the OWC as the main operational principle have been developed around the world: Sakata-Japan (60 kW), Mutriku-Spain (296 kW), Pico-Portugal (400 kW), Tofteshallen-Norway (500 kW), Islay island-Scotland (500 kW), and Lewis island-Scotland (4.0 MW) [17][18][19][20][21][22].
Important experimental works have sought to improve the comprehension of the fluid dynamic behavior of water/air flow in the OWC device and investigate the influence of some parameters over its performance. For instance, experiments in the laboratory and large-scale domains analyzed the influence of the inclination of the frontal wall, entrance areas of the OWC chamber, and water depth on the device efficiency, reflection, and loading of an OWC for different wave conditions [23][24][25]. Recently, the experimental progress extended to understand the hydrodynamic of the fluid flow into dual chambers OWC [26,27].
The numerical simulation of OWC devices has also been worth investigating. Several works have been performed since the development of computational models to represent the main operating principle of the device and the investigation of several parameters regarding the performance. For the former purpose, the representation of fluid flow in a laboratory and large-scale devices has been done without considering the effect of a turbine, as in Maciel et al. [28]. Other studies have considered it using orifice plates, plate-baffle, obstacles, or actuator disks to simulate the head loss caused by the turbine over the airflow in the hydropneumatic chamber and air duct of the device [29][30][31][32][33]. Recently, an interesting approach employed the numerical simulation of water and air in the chamber and considered the effects of Wells and impulse turbines by means of analytical thermodynamic models [34,35]. All models mentioned above have been used to obtain recommendations about parameters such as the depth and inclination of the frontal wall, height and length of the chamber, the diameter of the turbine, ramp placed in the seabed below the OWC chamber, and, recently, the design of multiple coupled chambers regarding device performance for different scales and wave conditions [24,32,[35][36][37][38][39][40]. It is also worth mentioning the efforts made to represent the sea state in a more trustworthy form. In this field, some studies modeled the irregular waves using a wave spectrum such as JONSWAP, and others obtained the sea state from geophysical models such as TOMAWAC and used this as an input for the modeling of a channel with the device to be investigated [41][42][43]. Some significant advancements regarding the PTO of OWC have also been reported in the literature. For instance, Britto-Melo et al. [44] numerically investigated the influence of the aerodynamic parameters of the Wells turbine and the influence of guide vanes and the bypass valve on the pressure drop, torque, and the overall performance. In this work, the conversion of pneumatic energy into electrical energy was estimated with a computational model based on the results extrapolated from aerodynamic tests on a scale model and empirical approximations for the generator losses. Recently, Rodríguez et al. [45] proposed a computational study to understand the behavior of OutFlow Radial (OFR) turbines in both direct and reverse modes, simulating an axisymmetric domain with the flow between a blade-to-blade arrangement of the turbine. The authors observed that the outer blade angle had poor performance in reverse mode despite the improvement of the global performance due to the rotor efficiency gain in direct mode. However, there are few studies related to the numerical simulation of the OWC device turbine considering the rotational domain with the intruded turbine rotor. Prasad et al. [46] performed similar work in this direction and developed a numerical model with a Savonius turbine immersed in a channel under regular wave flow, simulating a hydrokinetic turbine. The authors also investigated the influence of some parameters such as the submergence and the rotational speed for different blade entry angles regarding the rotor power.
Despite the several above-mentioned contributions, to the authors' knowledge, the development of computational models for the simulation of OWC devices considering the rotational turbine in the air duct is an approach still little explored in the literature. Recently, however, Liu et al. [47] presented the validation of an integrated three-dimensional numerical model considering an axial-flow impulse turbine coupled with an OWC inserted in a numerical wave tank (NWT). The present work aims to perform the first step in this direction. Initially, turbulent air flow over a free turbine inserted in a long and large channel (commonly used to represent the numerical modeling of wind turbines) was simulated to verify/validate the present model (Case 1). The effect of the tip speed ratio (λ) on the time and spatial averaged power coefficient (C P ) of a Savonius turbine was compared with numerical and experimental results available in the literature [48,49], investigating the tip speed ratios in the range 0.75 ≤ λ ≤ 2.00. An enclosure domain that mimics an OWC device was simulated considering a constant velocity imposed at the inlet of the domain (Case 2) and sinusoidal velocity that simulated the alternate flow in an OWC device (Case 3). For all cases, the influence of λ on the turbine performance and aerodynamic coefficients was investigated.
Moreover, incompressible, two-dimensional, unsteady, and turbulent flows with Re D = 867,000 were considered. URANS (Unsteady Reynolds-Averaged Navier Stokes) modeling was applied to all cases. Time-averaged equations of the conservation of mass, the balance of momentum, and transport equations of the k-ω SST model (used in the closure of turbulence) were solved with the finite volume method (FVM) [50][51][52][53][54][55], using the commercial code Ansys FLUENT 14.5 [56].

Mathematical Modeling
Incompressible, two-dimensional, unsteady, turbulent flows with constant thermophysical properties were considered. It is worth mentioning that one of the characteristics of turbulence is the three-dimensional structure of the flow [57]. Despite this fact, the present simulations represent most of the characteristics of turbulent flows, and the modeling properly predicts the transient behavior of parameters such as drag, lift, moment, and power in the Savonius turbine. Figure 1 illustrates the computational domain of Case 1 used for verification/validation of the present model. The domain consists of a long and large channel with an inserted free Savonius turbine. The dimensions are similar to those investigated in the work of Akwa et al. [49]. In this case, the fluid flow of air is caused by the imposition of a constant velocity (V ∞ ) at the inlet (left side surface). On this surface, there is also an imposed turbulent intensity of IT = u 2 /u = σ u /u = 1.0%, where u' is the fluctuation of the velocity field, u is the time-averaged velocity field, and σ u represents the variance of u'. At the exit of the channel (right side surface), there is an imposed null gauge pressure (p g = 0 atm). At the upper and lower surfaces, there is a free slip and impermeability boundary condition (also called symmetry). In the turbine region, there is an imposed constant angular velocity (n) in the gray region named the "Rotational Domain" simulating the effect of wind action over the turbine. In the turbine walls, there is an imposed no-slip and impermeability boundary condition (u = v = 0 m/s) related to the rotational domain. Figure 1 illustrates the details of the turbine region with the geometric variables used to design it. Table 1 presents the parameters of the fluid flow, thermo-physical properties, dimensions of the computational domain, and turbine variables used here for the four different tip speed ratios of the rotor (λ = nD/2V ∞ ) of λ = 0.75, 1.00, 1.25, and 2.00; n is the angular velocity of the turbine and D is the turbine diameter. For the unsteady analysis, there was a time interval of t f = 3.5 s, with the last 1.75 s being analyzed for computation of the drag, lift, momentum, and power coefficients.   The magnitude of the imposed velocity was defined in such a way to have a value of 7.0 m/s at the inlet of the air duct (in the contraction from the chamber to the air duct), leading to a Reynolds number in the turbine similar to that reached for Case 1. For Case 3, the sole difference in comparison with Case 2 was the imposition of a sinusoidal function that mimicked the oscillating behavior in the OWC chamber:

Description of the Studied Cases
where Hw = 0.4 m, Tw = 0.875 s, allowing the reproduction of the piston-type movement generated by the regular waves incidence over the OWC. The range of magnitudes was limited by (−1.4 m/s ≤ V(t) ≤ 1.4 m/s). The use of real configurations of a sea wave, for example, leads to long periods of simulation, requiring a high computational effort. As  The magnitude of the imposed velocity was defined in such a way to have a value of 7.0 m/s at the inlet of the air duct (in the contraction from the chamber to the air duct), leading to a Reynolds number in the turbine similar to that reached for Case 1. For Case 3, the sole difference in comparison with Case 2 was the imposition of a sinusoidal function that mimicked the oscillating behavior in the OWC chamber: where H w = 0.4 m, T w = 0.875 s, allowing the reproduction of the piston-type movement generated by the regular waves incidence over the OWC. The range of magnitudes was limited by (−1.4 m/s ≤ V(t) ≤ 1.4 m/s). The use of real configurations of a sea wave, for example, leads to long periods of simulation, requiring a high computational effort. As the purpose of this case is to compare the imposition of the sinusoidal velocity profile with the case with constant velocity to investigate the effect of oscillating flow over the device parameters and performance, the idealized imposed velocity variation was adequate for the desired investigation. For Case 2 and Case 3, the same four magnitudes of the tip speed ratio studied in the verification/validation case (λ = nD/2V 1 ) of λ = 0.75, 1.00, 1.25, and 2.00 were investigated. It is worth mentioning that a mean velocity was measured in the air duct and before the turbine (V 1 ) for the calculation of λ.  In the present work, results the influence of λ over the drag, lift, moment, and power coefficients (C d , C l , C T , and C P ) are based on [57]: (3) where F d is the drag force (N), F l is the lift force (N), T is the rotor moment (N·m), A r is the projected area of the Savonius rotor (A r = D·W), W is the depth of the domain in the z-direction (m), V is the upstream turbine velocity (V = V ∞ for the Case 1 and V = V 1 for Case 2 and Case 3) (m/s), r is the radius of the rotor (m), P is the turbine power (W), and P available is the available power of the wind upstream of the rotor (W). It is worth mentioning that the use of Equation (5) is valid only for the prediction of power coefficients in turbines subjected to open flow conditions (seen in tidal or wind power devices). In OWC devices, the air flow is driven by the pneumatic power, i.e., the pressure drop in the device must also be taken into account. Therefore, for power coefficients predicted for the enclosure domain of Case 2 and Case 3, the CP is calculated It is worth mentioning that the use of Equation (5) is valid only for the prediction of power coefficients in turbines subjected to open flow conditions (seen in tidal or wind power devices). In OWC devices, the air flow is driven by the pneumatic power, i.e., the pressure drop in the device must also be taken into account. Therefore, for power coefficients predicted for the enclosure domain of Case 2 and Case 3, the C P is calculated as: where ∆p is the pressure drop between the OWC chamber and the exit of the chimney (Pa), V is the air velocity at the OWC turbine-duct (m/s) (V = V 1 ), and Q is the volumetric flow rate of the air (m 3 /s). The time-averaged magnitudes of the drag, lift, moment, and power coefficients are obtained as follows:

Governing Equations of Turbulent Flows
For all simulations, the modeling of incompressible, two-dimensional, unsteady, and turbulent flows is given by the time-averaged conservation equation of mass and balance of momentum in the x and y directions and can be written as [50,51]: where x and y are the spatial coordinates (m), u and v are the velocity components in the x and y directions, respectively (m/s), p is the pressure (N/m 2 ), µ is the dynamic viscosity (kg/m·s), µ t is the turbulent viscosity (kg/m·s), and the overbar represents the time-averaged operator. For the k-ω SST closure model, the turbulent viscosity (µ t ) is [52,53]: The transport equations of turbulent kinetic energy (k) and its specific dissipation rate (ω) are as follows: ∂ω ∂t where P k is a function that prevents the turbulence generation in stagnation regions, i represents the direction of fluid flow (i = 1 represents the x direction and i = 2 represents the y direction), β = 0.09, α 1 = 5/9, β 1 = 3/40, σ k = 0.85, σ w = 0.5, σ 2 = 0.44, β 2 = 0.0828, σ k2 = 1, σ w2 = 0.856 are ad hoc constants used in [52]. F 1 and F 2 are blending functions defined as follows: In Equation (13), the term CD kω is calculated as follows:

Numerical Modeling
The solution of the governing equations was performed with the FVM using the commercial package Ansys Fluent TM [54][55][56]. To tackle the advective terms of the balance of momentum and transport of k and ω, the second-order upwind interpolation function was employed. The algorithm SIMPLE (Semi-Implicit Method for Pressure Linked Equations) was used for pressure-velocity coupling. The simulations were considered to have converged when the residuals for continuity, the balance of momentum, k and ω transport equations were less than 10 −5 . Moreover, the maximum number of iterations per time step was 200. Concerning the time advancement, there was an implicit time advancement scheme and a fixed time step of ∆t = 1.75 × 10 −3 s. All simulations were performed using desktops with six core Intel ® Core ™ i7 5820K @ 3.30 GHz processors and 16 Gb of RAM memory. The processing time for the simulation of t = 3.5 s of physical time was nearly 20 × 10 3 s.
Concerning spatial discretization, hybrid triangular and rectangular finite volumes were used with the domain in simulations of Case 1, Case 2, and Case 3. Figure 3a-c illustrates the mesh generated with software GMSH [58] for the free Savonius turbine (Case 1), the configuration similar to the OWC device (Case 2 and Case 3), and a detail of the mesh in the blades of Savonius turbine, respectively. In detail, it is possible to observe a region around the blades with the refined rectangular mesh. The dimensions of the rectangular volumes were defined as a function of a grid independence study and the parameter for representation of the boundary layer (y + ), which must be y + ≤ 1.0 in the walls; y + was defined as follows [51][52][53]: where y is the normal distance to the wall (m), τ w is the surface tension in the wall (N/m 2 ), ν is the kinematic viscosity (m 2 /s), and u τ is the friction velocity (m/s). For the grid independence study, four different meshes were simulated, and the results for the time-averaged power coefficient for the free Savonius rotor (Case 1) with Re D = 867,000, λ = 1.25, and bucket overlap ratio of R S = s/c = 0.15 are presented in Table 2. The mesh was considered independent when the relative difference between the results of C P for two successive grids met the criterion given by: where j represents the result obtained with the coarser mesh, and j + 1 represents the result obtained with the next successive refined mesh. Based on the results of Table 2, the mesh with 369,653 volumes was used in the independent grid. The same parameters used in this mesh were applied to the spatial discretization of the domain of Case 2 and Case 3.
parameter for representation of the boundary layer (y ), which must be y ≤ 1.0 in the walls; y + was defined as follows [51][52][53]: where y is the normal distance to the wall (m), τw is the surface tension in the wall (N/m²), ν is the kinematic viscosity (m²/s), and uτ is the friction velocity (m/s).

Results and Discussion
This section is divided into two parts, the verification/validation of the developed computational model (Case 1) and the investigation of the influence of the enclosure model and imposition of the sinusoidal velocity inlet over the aerodynamic and performance coefficients (Case 2 and Case 3). Figure 4 illustrates the instantaneous drag, lift, and moment coefficients as a function of time for Case 1 with Re D = 867,000, λ = 2.0, and R S = 0.15. For the first instants of time, mainly for t ≤ 0.5 s, the coefficients had a strong variation due to the incidence of the fluid flow and the imposition of angular velocity in the rotational domain region. Therefore, the present model must be used to predict power in the turbine when the flow is stabilized, and the first cycles of rotation of the turbine were disregarded for the analysis of coefficients and power. Here, for the computation of time-averaged parameters, only the results in the range of time 1.75 s ≤ t ≤ 3.5 s were used. Despite the complexity of the fluid flow, the results also demonstrated a regular oscillation in the magnitudes of C d , C l , and C T , which had similar behavior to that previously obtained in Akwa et al. [49]. Moreover, the crest and cave magnitudes were also similar, showing the generation of regular wakes of vortices behind the turbine.    Figures 6a and 7a, the velocity and pressure fields noticed were generated by the imposition of angular velocity in the rotational domain by the computational model. This behavior does not represent the real condition of fluid flow over the turbine since the turbine should not be in motion before the incidence of the fluid flow. Therefore, the present computational method must be used when the flow is stabilized, which happens a few cycles later at the beginning of the fluid flow (as shown in Figure 4 for monitoring the coefficients). As the time advances, mainly for t > 1.0 s, Figures 6c-e and 7c-e show an increase in the pressure field magnitude on the concave side of the advancement blade and a pressure drop as a consequence, being the main reason for the pressure drag force that drives the turbine. The results also show the fluid flowing in the region between the two blades. The results also show the generation of wakes behind the rotor and vortices generated in the tip region of the blade. It is also  Figure 5 shows the comparison of power coefficients (C P ) as a function of the tip speed ratio (λ) obtained with the present computational model and the numerical predictions of Akwa et al. [49] obtained with other commercial code (Star-CCM+) also based on the FVM and the experimental results of Blackwell et al. [48]. The results predicted with the present method are in close agreement with those previously obtained in the literature, verifying and validating the method used here. Even for λ = 1.00 and 2.00, where the highest differences between the present results and the experimental ones were obtained, the deviations were lower than the uncertainty of the experiment. The results indicated that the highest magnitudes of C P were reached at the lowest values of λ. With the increase in λ, the magnitude of C P had a slight decrease in the range 1.00 ≤ λ ≤ 1.25 and a step decrease for λ ≥ 1.25. For λ = 2.00, negative magnitudes of C P were obtained, indicating that the device supplies energy for the fluid flow and not the contrary.    Figures 6a and 7a, the velocity and pressure fields noticed were generated by the imposition of angular velocity in the rotational   Figures 6a and 7a, the velocity and pressure fields noticed were generated by the imposition of angular velocity in the rotational domain by the computational model. This behavior does not represent the real condition of fluid flow over the turbine since the turbine should not be in motion before the incidence of the fluid flow. Therefore, the present computational method must be used when the flow is stabilized, which happens a few cycles later at the beginning of the fluid flow (as shown in Figure 4 for monitoring the coefficients). As the time advances, mainly for t > 1.0 s, Figures 6c-e and 7c-e show an increase in the pressure field magnitude on the concave side of the advancement blade and a pressure drop as a consequence, being the main reason for the pressure drag force that drives the turbine. The results also show the fluid flowing in the region between the two blades. The results also show the generation of wakes behind the rotor and vortices generated in the tip region of the blade. It is also worth mentioning that the behavior found here is similar to that described in previous literature, e.g., in Prasad et al. [45] and Blackwell et al. [49]. worth mentioning that the behavior found here is similar to that described in previous literature, e.g., in Prasad et al. [45] and Blackwell et al. [49].

The Results of the Savonius Turbine Inserted in an OWC Domain (Case 2 and Case 3)
For the turbulent air flow in the enclosure domain, the instantaneous aerodynamic coefficients (C d , C l , and C T ) were obtained as a function of time to better understand the behavior of the turbine in an OWC domain. Figure 8 depicts the coefficients for the flow with Re D = 867,000, λ = 2.0, and R S = 0.15 for a constant imposed velocity at the inlet (Case 2). In general, the results show an average increase in C d (black line) for all instants of time compared to the free turbine simulations (Case 1) due to the insertion of the turbine in the enclosure domain. The increase in the drag coefficient (C d ) also led to an augmentation of the moment coefficient (C T ) represented in red color in a similar form reached for the C d , while the lift coefficient (C l ) did not suffer important modifications in its mean magnitude. The results also indicate that the transient behavior of the coefficients was strongly modified compared to Case 1, showing C d and C l with sharp peaks and smooth troughs, i.e., each cycle did not behave in a sinusoidal form as noticed in Figure 4 for the free turbine configuration. It is worth mentioning that, for other magnitudes of λ, similar behavior for the instantaneous coefficients was obtained. Therefore, for the sake of brevity, the instantaneous coefficients for other magnitudes of λ are not presented.
The results also indicate that the transient behavior of the coefficients was strongly modified compared to Case 1, showing Cd and Cl with sharp peaks and smooth troughs i.e., each cycle did not behave in a sinusoidal form as noticed in Figure 4 for the free turbine configuration. It is worth mentioning that, for other magnitudes of λ, similar behavior for the instantaneous coefficients was obtained. Therefore, for the sake of brevity, the instantaneous coefficients for other magnitudes of λ are not presented. As previously mentioned, the air flow in the OWC chamber is subjected to the pistontype oscillatory motion of the water column (hydropneumatic chamber). In order to simulate this effect, the results of instantaneous power coefficient obtained for Case 3 and Case 2, considering the same conditions (ReD = 867,000, λ = 2.0, and RS = 0.15), are presented in Figure 9. It is important to reinforce that, for prediction of CP for Case 2 and Case 3 Equation (6) was used instead of Equation (5), and the time-averaged magnitudes were predicted with Equation (7). The results reveal a strong similarity between the instantaneous CP for both cases, indicating that the imposition of sinusoidal velocity at the inlet of the domain did not have a significant influence over the pattern of the CP investigated here, which is not intuitively expected since the mean imposed momentum decreased for Case 3 in comparison with Case 2. The results also demonstrate a slight increase in the power coefficients when the sinusoidal velocity was imposed, with a As previously mentioned, the air flow in the OWC chamber is subjected to the pistontype oscillatory motion of the water column (hydropneumatic chamber). In order to simulate this effect, the results of instantaneous power coefficient obtained for Case 3 and Case 2, considering the same conditions (Re D = 867,000, λ = 2.0, and R S = 0.15), are presented in Figure 9. It is important to reinforce that, for prediction of C P for Case 2 and Case 3, Equation (6) was used instead of Equation (5), and the time-averaged magnitudes were predicted with Equation (7). The results reveal a strong similarity between the instantaneous C P for both cases, indicating that the imposition of sinusoidal velocity at the inlet of the domain did not have a significant influence over the pattern of the C P investigated here, which is not intuitively expected since the mean imposed momentum decreased for Case 3 in comparison with Case 2. The results also demonstrate a slight increase in the power coefficients when the sinusoidal velocity was imposed, with a difference of nearly 9.0% on average. Possible explanations for the behavior found here include the synchronization of velocity augmentation with the rotation of the turbine and the imposition of the rotational domain, as well as the acceleration of the fluid caused by the variation of the imposed velocity at the inlet for Case 3, which did not happen for Case 2. The magnitude of the imposed mean velocity at the inlet of Case 2 and Case 3 and acceleration at the inlet for Case 3 are illustrated in Figure 10 to make this visualization easy. It is worth mentioning that, despite different imposed inlet velocities, the transient fields of velocity and pressure had only slight differences. Future investigations should be performed with other magnitudes of amplitude and periods of imposed velocity to corroborate this hypothesis and the development of other models where the inertia moment of the turbine is taken into account.
2. The magnitude of the imposed mean velocity at the inlet of Case 2 and Case 3 and acceleration at the inlet for Case 3 are illustrated in Figure 10 to make this visualization easy. It is worth mentioning that, despite different imposed inlet velocities, the transient fields of velocity and pressure had only slight differences. Future investigations should be performed with other magnitudes of amplitude and periods of imposed velocity to corroborate this hypothesis and the development of other models where the inertia moment of the turbine is taken into account. and Case 2 were almost the same, only the fields for one of the cases are illustrated. As previously mentioned, for the free turbine case (see Figures 6 and 7), the initial time steps should be disregarded due to the artificial imposition of the angular velocity in the turbine region. In spite of that, when the flow was stabilized, the model properly represented the physical problem. In turn, one can note in Figures 11 and 12 an increase in the pressure difference between the concave and convex sides of the advancement blade for Case 2 and Case 3 compared to Case 1, which explains the augmentation of the pressure drag in comparison with Case 1. Moreover, the pressure difference between the upstream and downstream regions of the turbine was also augmented. This influence over the water oscillation was not investigated, but the results pointed out that this aspect is worth investigation when the turbine is taken into account in the problem in order to avoid the restriction of the water column, mainly when coupling with a wave channel is performed. The results also indicated the variation of the pressure magnitude in the chamber for different instants of time, explaining the non-symmetric differences between the peaks and troughs with respect to the mean magnitudes of the coefficients. The momentum of the fluid flow intensified due to the insertion of the turbine in the enclosure domain. Moreover, the wakes generated in the turbine and secondary vortices in the tip of the blades could not be spread in the spanwise direction of the main flow due to the limitation imposed by the air duct walls. Therefore, as expected, the insertion of the turbine in the enclosure domain affected the fluid dynamic behavior of the flow considerably.
To summarize, Figure 13 shows the effect of λ on the time-averaged magnitudes of , ,   As the behaviors of Case 3 and Case 2 were almost the same, only the fields for one of the cases are illustrated. As previously mentioned, for the free turbine case (see Figures 6 and 7), the initial time steps should be disregarded due to the artificial imposition of the angular velocity in the turbine region. In spite of that, when the flow was stabilized, the model properly represented the physical problem. In turn, one can note in Figures 11 and 12 an increase in the pressure difference between the concave and convex sides of the advancement blade for Case 2 and Case 3 compared to Case 1, which explains the augmentation of the pressure drag in comparison with Case 1.
Moreover, the pressure difference between the upstream and downstream regions of the turbine was also augmented. This influence over the water oscillation was not investigated, but the results pointed out that this aspect is worth investigation when the turbine is taken into account in the problem in order to avoid the restriction of the water column, mainly when coupling with a wave channel is performed. The results also indicated the variation of the pressure magnitude in the chamber for different instants of time, explaining the non-symmetric differences between the peaks and troughs with respect to the mean magnitudes of the coefficients. The momentum of the fluid flow intensified due to the insertion of the turbine in the enclosure domain. Moreover, the wakes generated in the turbine and secondary vortices in the tip of the blades could not be spread in the spanwise direction of the main flow due to the limitation imposed by the air duct walls. Therefore, as expected, the insertion of the turbine in the enclosure domain affected the fluid dynamic behavior of the flow considerably.     To summarize, Figure 13 shows the effect of λ on the time-averaged magnitudes of C d , C l , C T and C P for Case 2 and Case 3, and Figure 14 illustrates the same effect considering only C P for Case 1, Case 2, and Case 3. The results indicate that the differences between the values of the aerodynamic and power coefficients obtained for Case 2 and Case 3 were not significant, with the highest difference being lower than 10.0%. The effect of λ on the coefficients was also similar for both cases. Concerning the magnitude of C d , it decreased in the range 0.75 ≤ λ ≤ 1.00, and after this point, the magnitude was almost constant. The C l and C T magnitudes decreased with the increase in λ. However, for C T the decrease was more significant only in the range 1.00 ≤ λ ≤ 1.25, being almost constant in other intervals of λ. For the C P , there was an increase with the augmentation of λ. For Case 3, for example, the magnitude increased from C P = 0.2922 when λ = 0.75 to C P = 0.5054 when λ = 2.0. Figure 14 indicates that in the region 0.75 ≤ λ ≤ 1.25, a similar trend of C P was noticed when Case 1 and the enclosure cases (Case 2 and Case 3) were compared. However, for λ = 2.0, contrary to the behavior noticed for Case 1, Case 2 and Case 3 did not have a reduction of C P for the highest magnitude of λ investigated. Therefore, the results indicated that changes in the domain where the turbine is placed are important to define its application range and the effect of λ over the problem performance. Further studies should be performed to define the range of application of λ for the enclosure domain.

Conclusions
The present work developed a computational model to investigate turbulent flows in enclosure domains with an inserted Savonius turbine simulating the air flow in an OWC/WEC device. Initially, a free turbine configuration (Case 1) was investigated in order to perform the verification/validation of the present method. Then, two different study cases were investigated (Case 2 and Case 3) with constant and sinusoidal velocity imposed at the domain inlet. For all cases, ReD = 867,000, RS = 0.15, and four different magnitudes

Conclusions
The present work developed a computational model to investigate turbulent flows in enclosure domains with an inserted Savonius turbine simulating the air flow in an OWC/WEC device. Initially, a free turbine configuration (Case 1) was investigated in order to perform the verification/validation of the present method. Then, two different study cases were investigated (Case 2 and Case 3) with constant and sinusoidal velocity imposed at the domain inlet. For all cases, ReD = 867,000, RS = 0.15, and four different magnitudes

Conclusions
The present work developed a computational model to investigate turbulent flows in enclosure domains with an inserted Savonius turbine simulating the air flow in an OWC/WEC device. Initially, a free turbine configuration (Case 1) was investigated in order to perform the verification/validation of the present method. Then, two different study cases were investigated (Case 2 and Case 3) with constant and sinusoidal velocity imposed at the domain inlet. For all cases, Re D = 867,000, R S = 0.15, and four different magnitudes of λ were studied (λ = 0.75, 1.00, 1.25, and 2.00). For the simulation of incompressible, two-dimensional, transient, and turbulent flows, the time-averaged equations of mass conservation, the balance of momentum, and transport equations of the k-ω SST model were solved with the FVM. The proposed computational model was verified and validated using a comparison of C P with the numerical and experimental results of the literature [48,49], reproducing the effect of λ on C P for the turbulent flow over a free Savonius turbine. Moreover, the velocity and pressure fields demonstrated that the driven force of the Savonius turbine was dominated by the pressure difference between the concave and convex sides of the advancement blade, which the present model adequately predicted.
After the verification/validation of the computational model, new recommendations were reached for the turbulent flows over the Savonius turbine inserted in the enclosure domain representing an OWC. The results demonstrated that the new configuration had a strong influence over the behavior and magnitudes of the instantaneous coefficients such as C d and C T , also affecting the instantaneous values of C P . For the investigated values of λ, the results indicated that the insertion of a turbine in the enclosure domain led to an overall augmentation of C P for all values of λ. In addition, the effect of λ on C d , C l , C T , and C P was strongly affected by the domain change, indicating that the geometric configuration can be important to define the range of applicability of the turbine and the design of the PTO in the OWC device. The comparison between the imposition of sinusoidal velocity (Case 3) and constant velocity (Case 2) at the domain inlet led to similar performance and aerodynamic coefficients. It is important to highlight that due to the two-dimensional approach, some characteristics of the air flow phenomenology through the Savonious turbine cannot be completely revealed from the obtained results, being a limitation of the proposed computational model. Future investigations should be performed considering a model where the inertia moment is taken into account to verify if this behavior continues to be observed. It is recommended to investigate other parameters, such as the effect of the overlap and spacing between the blades (s and a) and other tip speed ratios (λ) on the aerodynamic and performance of the OWC and the imposition of irregular velocity variation to represent the real sea state movement. It is also recommended to investigate the coupling between the present model and the wave channel for adequate simulation of the interaction between the oncoming waves, the structure of the device, and air flow over the turbine inserted in the OWC air duct.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author. The data are not publicly available due to privacy reasons.