Hydrodynamic Analysis of a Multibody Wave Energy Converter in Regular Waves

A performance assessment of wave power absorption characteristics of isolated and multiple wave energy converter (WEC) rotors was presented in this study for various wave-heading angles and wave frequencies. Numerical hydrodynamic analysis of the WEC was carried out using the three-dimensional linear boundary element method (BEM) and nonlinear computational fluid dynamics (CFD). Experimental results were used to validate the adopted numerical models. Influence with and without power take-off (PTO) was estimated on both isolated and multiple WEC rotors. Furthermore, to investigate the interaction effect among WECs, a q-factor was used. Incorporation of viscous and PTO damping into the linear BEM solution shows the maximum reduction focused around peak frequency but demonstrated an insignificant effect elsewhere. The q-factor showed both constructive and destructive interactions with the increase of the wave-heading angle and wave frequencies. Further investigation based on the prototype WEC rotor was carried, and calculated results of the linear BEM and the nonlinear CFD were compared. The pitch response and q-factor of the chosen wave frequencies demonstrated satisfactory consistency between the linear BEM and nonlinear CFD results, except for some wave frequencies. Estimated optimal time-averaged power using linear BEM show that the maximum extracted power close to the zero wave-heading angle around the resonance frequency decreases as the wave-heading angle increases. Overall, the linear BEM on the extracted power is overestimated compared with the nonlinear CFD results.


Introduction
Commercial development of many applications in the field of wave energy converter (WEC) is motivated to deploy in arrays for creating WEC farms. Understanding the hydrodynamic behavior in such arrays can determine the optimum configuration of multiple WECs and maximize wave power extraction. Wave approach angle and large-scale interference among WECs are important in multiple WEC spacing of hydrodynamic interaction of waves. A q-factor can be used to quantify the interaction among WEC, which is defined as the ratio of the total response/power extracted from a multiple of n WECs to n times the response/power extracted from an isolated WEC, where n is the number of WECs. Interaction effects may be constructive (q-factor > 1) or destructive (q-factor < 1) depending on the configuration of multiple WECs and sea conditions.

Mathematical Model and Experimental Setup
The configuration of multiple WEC rotors in the present study is shown in Figure 2. Three WEC rotors with a width (W) in the y-direction are considered. A fixed coordinate system is taken at the center of rotation (O) with a vertical upward z-axis, and the wave propagates from the negative x-direction to the positive x-direction. The rotor tip is placed in the "negative" x and approach wave angle (β) is set to 0°, 30°, and 60° to allow waves to hit the WEC rotor. Each rotor rotates around its own center of rotation. Geometry details of principle particulars; moment of inertia in x-, y-, and z-directions; and hydrostatic properties are listed in Table 1. With the presented background, the modeling of multiple WEC rotors is considered an extension of the framework of isolated WEC rotors carried out in previous studies [27][28][29][30][31]. The interaction of multiple WEC rotors is investigated in regular waves on the basis of linear boundary element method (BEM) (using WAMIT) and nonlinear CFD (using Star-CCM+). Each WEC rotor has one degree of freedom that includes pitch alone along the y-axis, as shown in Figure 1. The rest of this paper is organized as follows. The mathematical model and experimental details of a WEC rotor are outlined in Section 2. The CFD methodology and numerical settings use to simulate the wave and its interaction with multiple WEC rotors are discussed in Section 3. The validation, results, and discussion are presented in Section 4. Finally, conclusions of this study are drawn in Section 5.

Mathematical Model and Experimental Setup
The configuration of multiple WEC rotors in the present study is shown in Figure 2. Three WEC rotors with a width (W) in the y-direction are considered. A fixed coordinate system is taken at the center of rotation (O) with a vertical upward z-axis, and the wave propagates from the negative x-direction to the positive x-direction. The rotor tip is placed in the "negative" x and approach wave angle (β) is set to 0 • , 30 • , and 60 • to allow waves to hit the WEC rotor. Each rotor rotates around its own center of rotation. Geometry details of principle particulars; moment of inertia in x-, y-, and z-directions; and hydrostatic properties are listed in Table 1.

Mathematical Model and Experimental Setup
The configuration of multiple WEC rotors in the present study is shown in Figure 2. Three WEC rotors with a width (W) in the y-direction are considered. A fixed coordinate system is taken at the center of rotation (O) with a vertical upward z-axis, and the wave propagates from the negative x-direction to the positive x-direction. The rotor tip is placed in the "negative" x and approach wave angle (β) is set to 0°, 30°, and 60° to allow waves to hit the WEC rotor. Each rotor rotates around its own center of rotation. Geometry details of principle particulars; moment of inertia in x-, y-, and z-directions; and hydrostatic properties are listed in Table 1. This study primarily deals with the hydrodynamic response of oscillating multiple WEC rotors within the framework of the linear BEM model under the assumption of irrotational flow, inviscid and incompressible fluid, linearized kinematic and dynamic boundary conditions on the free surface, and smaller body motions compared with their characteristic length.  This study primarily deals with the hydrodynamic response of oscillating multiple WEC rotors within the framework of the linear BEM model under the assumption of irrotational flow, inviscid and incompressible fluid, linearized kinematic and dynamic boundary conditions on the free surface, and smaller body motions compared with their characteristic length.
All variables can be expressed as complex amplitude multiplied by e −iωt under these assumptions. Angular oscillation amplitude ξ j as a harmonic function of time is expressed as follows: Hydrodynamic response of the M th WEC excited by waves with amplitude "a" angular frequency ω along the pitch direction (denoted subscript 5) can be expressed as where • n is the total number of WECs (n = 3 in the present study), • J M5 is the mass moment of inertia of the M th WEC in the pitch direction, The energy that can be extracted from the WEC rotor is captured by the PTO. The constant PTO (B PTO M5 ) damping obtained from the isolated WEC rotor was used for multiple WEC rotors. Free decay tests were conducted on the isolated WEC rotor to obtain the viscous damping moment from laboratory tests in Jeju National University [31]. Logarithmic decrement method is used to obtain the damping ratio. The viscous damping can be estimated by subtracting radiation damping from actual damping. The actual damping is defined as the product of the damping ratio and critical damping where ω N is the undamped natural frequency of the isolated WEC rotor. Viscous damping values from experiments, hydrodynamic coefficients from the linear BEM analysis, and PTO damping from Equation (5) were used to solve the motion Equation (2).
The time-averaged M th WEC power based on linear theory can be expressed as follows: where B PTO M5 is the optimal power take-off given by M th WEC rotor satisfying ∂P opt M5 (ω)/∂ B PTO M5 = 0, ξ M5 a (ω) is the pitch response amplitude operator (RAO) (from Equation (2)) and are expressed as ω N is the undamped natural frequency of the pitch direction, is given by The interaction factor known as q-factor was used to evaluate the interaction effect of multiple WEC rotors. The energy yield can be quantified because the extracted power is directly proportional to the square of the response amplitude. Interactions are constructive with respect to energy yield when q-factor > 1, and the energy diminishes if q-factor < 1 because destructive effects appear. The q-factor response amplitude from "M" WECs ( ξ M5 a (ω)) can be expressed as follows: where ξ 5 ,iso a (ω) is the isolated WEC rotor pitch RAO. A full-scale prototype WEC rotor (see Table 1) was proposed for real-time testing at the west coast of Jeju Island, South Korea. Initial validation of numerical models on the model-scale WEC rotor is compared with the experimental results prior to the analysis of the prototype WEC rotor. Properties of the 1:11 scale-down model of the WEC rotor experimentally investigated in Jeju National University are presented in Table 1. The experimental setup of the model made of acrylic material is illustrated in Figure 3. The experimental setup consists of the following:

•
The wave tank is equipped with a hydraulic piston-type wave maker at one end to generate waves and a wave absorber at the other end.

•
The WEC rotates due to wave excitation around a fixed-axis shaft (y-direction) located below the waterline, and both ends of the shaft are connected to vertical rods, which are supported from the fixed frame along the wave flume. • A ball bearing mechanism is provided for smooth rotation in the pitch motion.

•
An image processing technique is used to measure the pitch motion of the WEC rotor.

Computational Fluid Dynamics
The problem domain is solved using three-dimensional implicit, unsteady, and incompressible Navier-Stokes equations. The computational domain is discretized into a small number of finite control volumes with a grid, which uses integral form of conservation equation. Quantities in governing equations are decomposed into their mean (Ui, P) and fluctuating ( ' , ' i u p ) components, where Reynolds-averaged Navier-Stokes (RANS) equations along with the continuity are expressed as follows:

Computational Fluid Dynamics
The problem domain is solved using three-dimensional implicit, unsteady, and incompressible Navier-Stokes equations. The computational domain is discretized into a small number of finite control volumes with a grid, which uses integral form of conser-vation equation. Quantities in governing equations are decomposed into their mean (U i , P) and fluctuating (u i , p ) components, where Reynolds-averaged Navier-Stokes (RANS) equations along with the continuity are expressed as follows: where the velocity field is U i (x, y, z, t), P is pressure, ρ is the constant density of the fluid, δ ij is Kronecker delta function and υ is the kinematic viscosity of the fluid. The WEC is allowed to move by only a single degree-of-freedom in the pitch around the y-axis where the rigid body motion was carried out using dynamic fluid body interaction. The motion equation of the M th WEC rotor is expressed as follows: where I wec M is the moment of inertia of the M th WEC rotor, ω M is the angular velocity of the M th WEC rotor, n wec M is the resultant moment acting on the M th WEC rotor, and E d M is the external damping moment of the M th WEC rotor.
The computational domain is set in the Cartesian coordinate system Oxyz with the origin located at the COR of the center WEC rotor 2, where the x-direction is the propagation of the wave from the negative x-axis to the positive x-axis and positive z-axis pointing upwards. For prototype WEC rotor see Table 1, three-dimensional numerical wave tank (NWT) is used to generate desired waves, and the flow is controlled by boundary conditions assigned to different regions in the domain, as shown in Figure 4a. The velocity inlet is located at the upstream side of the model, while the pressure outlet boundary is found in the downstream side; the top is a pressure outlet and all other boundaries are wall boundaries except the y = 0 symmetry. Wave forcing was applied at inlet and outlet boundaries with 0.5 × λ, where λ is the wave length. The forcing is performed using fifth-order Stokes wave theory. The middle portion in the vicinity of the WEC rotor of the computational domain was chosen as 2 × λ. Fifth-order Stokes waves were generated to represent realistic regular waves. Water depth is fixed at 40 m.
A gap equal to the width of the WEC rotor between the WEC rotor 3 and the wall in the y-direction was considered for all test cases. An overset meshing strategy was selected to address the WEC response. An overset mesh typically contains a WEC rotor superimposed on a background mesh containing the surrounding geometry and data are interpolated between the two meshes [32]. Multiple overset regions are created for each WEC rotor in the present study. Only half of the computational domain was chosen due to the symmetrical configuration of the xz-plane, and this cannot be implemented when the wave-heading is different from zero. Component grids are allowed to overlap and flow field information is transmitted from one grid to another through the overlap region. A mesh size that matches the overset boundary must be maintained to use this method effectively because conservation issues will occur when large discrepancies exist in the mesh size the grid technique. The cell size along x-(∆x/∆z) and y-(∆y/∆z) directions was set to 4 to minimize conservation errors (Figure 4a,b). The cell size in overset regions and the background was matched using overlapping reason. Linear interpolation technique is used to transfer the information between the two grids. The cell size around the WEC rotor is reduced to 15% of the basic cell size of the overset region (base size = 1.0 m), and three cells within the prism layer were generated. All − y+ wall along with the two-equation low-Re k − ε model was adopted for simulations. Automatic convective Courant-Friedrichs-Levy (CFL) was utilized to stabilize the numerical solution and choose a suitable time step. The time integration is based on a second-order implicit scheme with 10 inner iterations per time step is adopted. The volume of fluid (VOF) is used to handle the free surface interface between the two immiscible fluids (water and air) and to improve the interface tracking high-resolution interface capturing, along with interface momentum dissipation, is used.

Validation
Adopted numerical models of the linear BEM and nonlinear CFD were validated with the experimental results, and for the investigation a single WEC rotor without PTO damping was used. The advantage of any numerical model can be nullified if unsuitable settings are selected. The first and second authors performed an in-depth sensitivity analysis on the single WEC rotor and compared the linear BEM and nonlinear CFD results with the findings in the literature and experimental results in Poguluri and Bae and Poguluri et al., respectively [27,28]. The details of the NWT used for validation can be found in [28].

Validation
Adopted numerical models of the linear BEM and nonlinear CFD were validated with the experimental results, and for the investigation a single WEC rotor without PTO damping was used. The advantage of any numerical model can be nullified if unsuitable settings are selected. The first and second authors performed an in-depth sensitivity analysis on the single WEC rotor and compared the linear BEM and nonlinear CFD results with the findings in the literature and experimental results in Poguluri and Bae and Poguluri et al., respectively [27,28]. The details of the NWT used for validation can be found in [28].
The WEC rotor in the linear BEM is discretized with 1200 panels using Multi-surf, and only half of the model was used for analysis due to the symmetry of the xz-plane. The damping ratio of 0.0849 and the corresponding viscous damping value of 4.113 N·m·s/rad were obtained from the experimental free decay test conducted by Kim et al. [31]. The influence of domain, mesh size, mesh moving strategy around the WEC rotor, and associated domain boundary conditions in nonlinear CFD simulations were referred to in [28] and omitted in the present paper. Tests are conducted with the following fixed wave parameters: wave frequency range: 3.25-6.0 rad/s, wave amplitude (a): 0.005 m, and water depth (h): 0.6 m. Non-dimensional pitch RAO *, ( where L is the characteristic length of the model equal to 1) in the experiments was plotted along with the linear BEM and nonlinear CFD results for different wave frequencies, as shown in Figure 5. The chosen value of the damping ratio was based on the acceptable match between the maximum pitch motion response of the WEC rotor and the experimental results; see [31]. Consistency between the linear BEM results and the experimental and nonlinear CFD findings is only possible when incident, diffracted, and radiated waves are linear and the same with the present chosen wave frequencies. The means of three to four cycles of the steady solution from the nonlinear CFD and experimental results were used to obtain the single value of the pitch response. The present linear BEM and nonlinear CFD results are consistent with the experimental results except at wave frequencies close to the resonance region. The nonlinear CFD and experimental results slightly overestimate the non-dimensional RAO compared with the linear BEM but the overall consistency was relatively acceptable. Furthermore, the time history of the pitch response of the WEC rotor is compared with the experimental results for two wave frequencies of 3.70 and 3.98 rad/s (see Figure 6). The pitching response of the WEC rotor was acceptable for four cycles compared with the experimental results. The experimental pitch response was then affected by wave reflections from both ends of the wave flume and eventually led to differences in amplitude and phase. In the foregoing sections, additional investigations on the prototype isolated and multiple WEC rotors will be carried out by retaining most of the numerical settings of the model unless otherwise mentioned. The WEC rotor in the linear BEM is discretized with 1200 panels using Multi-surf, and only half of the model was used for analysis due to the symmetry of the xz-plane. The damping ratio of 0.0849 and the corresponding viscous damping value of 4.113 N⋅m⋅s/rad were obtained from the experimental free decay test conducted by Kim et al. [31]. The influence of domain, mesh size, mesh moving strategy around the WEC rotor, and associated domain boundary conditions in nonlinear CFD simulations were referred to in [28] and omitted in the present paper. Tests are conducted with the following fixed wave parameters: wave frequency range: 3.25-6.0 rad/s, wave amplitude (a): 0.005 m, and water depth (h): 0.6 m. Non-dimensional pitch RAO *, where L is the characteristic length of the model equal to 1) in the experiments was plotted along with the linear BEM and nonlinear CFD results for different wave frequencies, as shown in Figure  5. The chosen value of the damping ratio was based on the acceptable match between the maximum pitch motion response of the WEC rotor and the experimental results; see [31]. Consistency between the linear BEM results and the experimental and nonlinear CFD findings is only possible when incident, diffracted, and radiated waves are linear and the same with the present chosen wave frequencies. The means of three to four cycles of the steady solution from the nonlinear CFD and experimental results were used to obtain the single value of the pitch response. The present linear BEM and nonlinear CFD results are consistent with the experimental results except at wave frequencies close to the resonance region. The nonlinear CFD and experimental results slightly overestimate the non-dimensional RAO compared with the linear BEM but the overall consistency was relatively acceptable. Furthermore, the time history of the pitch response of the WEC rotor is compared with the experimental results for two wave frequencies of 3.70 and 3.98 rad/s (see Figure 6). The pitching response of the WEC rotor was acceptable for four cycles compared with the experimental results. The experimental pitch response was then affected by wave reflections from both ends of the wave flume and eventually led to differences in amplitude and phase. In the foregoing sections, additional investigations on the prototype isolated and multiple WEC rotors will be carried out by retaining most of the numerical settings of the model unless otherwise mentioned.    Figure 7 shows the pitch RAO of multiple WEC rotors by changing wave frequencies from 0.5 rad/s to 2.0 rad/s with an increment of 0.02 rad/s for different wave-heading angles β = 0°, 30°, and 60°. Line denotes the isolated WEC rotor and symbols represent multiple WEC rotors. The viscous damping is 25.0 kN⋅m⋅s/rad (see [31]). The spacing between rotors was tested for 10 m, which is the center-to-center distance. Notably, the pitch RAO of isolated and multiple WEC rotors

Linear BEM Results
predicts a similar hydrodynamic behavior with a range of wave frequencies. The isolated WEC rotor at a small wave-heading angle predicts increased RAO within the range of frequencies considered, while peak amplitudes are slightly affected in the case of multiple WEC rotors. Peak pitch RAO variations according to the wave frequency and wave-heading angle were observed, and the main resonance frequency was within the range of 1.2-1.23 rad/s (see Table 2). The influence of PTO damping on the WEC RAO is illustrated in Figure 8. Equation (3) shows that the estimated minimum PTO damping of the WEC rotor is 34.6 kN⋅m⋅s/rad. Figure 8 clearly shows the maximum reduction focused around the peak frequency although the influence is in-significance as moving away from the peak. The PTO damping has a strong influence on peak frequencies at all wave-heading angles (see Table 2). Figure 9 shows the variation of the q-factor (given by Equation (6)) as a function of the wave-heading angle and frequency for with and without PTO damping. Interaction effects are constructive (see the dark red contour in Figure 9) within the range of (0.5 < ω < 1.08 and β ≤ 60°) and (1.34 < ω < 1.82 and β < 30°) for without PTO and (0.52 < ω < 1.1 and β ≤ 60°) and (1.33 < ω < 1.77 and β < 30°) for with PTO but destructive otherwise.  Figure 7 shows the pitch RAO of multiple WEC rotors by changing wave frequencies from 0.5 rad/s to 2.0 rad/s with an increment of 0.02 rad/s for different wave-heading angles β = 0 • , 30 • , and 60 • . Line denotes the isolated WEC rotor and symbols represent multiple WEC rotors. The viscous damping is 25.0 kN·m·s/rad (see [31]). The spacing between rotors was tested for 10 m, which is the center-to-center distance. Notably, the pitch RAO of isolated and multiple WEC rotors ξ 5,iso a , ξ M5 a predicts a similar hydrodynamic behavior with a range of wave frequencies. The isolated WEC rotor at a small waveheading angle predicts increased RAO within the range of frequencies considered, while peak amplitudes are slightly affected in the case of multiple WEC rotors. Peak pitch RAO variations according to the wave frequency and wave-heading angle were observed, and the main resonance frequency was within the range of 1.2-1.23 rad/s (see Table 2). The influence of PTO damping on the WEC RAO is illustrated in Figure 8. Equation (3) shows that the estimated minimum PTO damping of the WEC rotor is 34.6 kN·m·s/rad. Figure 8 clearly shows the maximum reduction focused around the peak frequency although the influence is in-significance as moving away from the peak. The PTO damping has a strong influence on peak frequencies at all wave-heading angles (see Table 2). Figure 9 shows the variation of the q-factor (given by Equation (6)) as a function of the wave-heading angle and frequency for with and without PTO damping. Interaction effects are constructive (see the dark red contour in Figure 9) within the range of (0.5 < ω < 1.08 and β ≤ 60 • ) and (1.34 < ω < 1.82 and β < 30 • ) for without PTO and (0.52 < ω < 1.1 and β ≤ 60 • ) and (1.33 < ω < 1.77 and β < 30 • ) for with PTO but destructive otherwise.  • Each cell represents (peak frequency (rad/s), peak pitch RAO(rad/m)).

Nonlinear CFD Results
NWT simulations were performed for a fixed spacing (center-to-center spacing = 10 m), zero wave-heading angle, wave height = 1.5 m (target wave height), and three wave frequencies of 1.09, 1.20, and 1.32 rad/s. The time history of angular velocity and pitch response of isolated and multiple WEC rotors with and without PTO are presented in Figures 10 and 11 and Table 3. The time history shows that the WEC responds nonlinearly to incoming waves. The angular velocity and pitch response signal demonstrated that crests and troughs are affected by the applied PTO while the phase slightly increases over time. These changes were significant in the case of ω = 1.32 rad/s and a similar trend was observed with multiple WEC rotors but with reduced magnitudes.

Nonlinear CFD Results
NWT simulations were performed for a fixed spacing (center-to-center spacing = 10 m), zero wave-heading angle, wave height = 1.5 m (target wave height), and three wave frequencies of 1.09, 1.20, and 1.32 rad/s. The time history of angular velocity and pitch response of isolated and multiple WEC rotors with and without PTO are presented in Figures 10 and 11 and Table 3. The time history shows that the WEC responds nonlinearly to incoming waves. The angular velocity and pitch response signal demonstrated that crests and troughs are affected by the applied PTO while the phase slightly increases over time. These changes were significant in the case of ω = 1.32 rad/s and a similar trend was observed with multiple WEC rotors but with reduced magnitudes.  Table 3 presents the pitch RAO and q-factor comparison between the linear BEM nonlinear CFD results. The pitch RAO of wave frequencies (1.02, 1.20, and 1.32 r without PTO shows the satisfactory overall consistency between the linear BEM nonlinear CFD solutions, except for the significant variation in the case of 1.20 rad/s linear BEM results demonstrate higher destruction in the q-factor compared with nonlinear CFD findings. The overall consistency of the pitch RAO with PTO was sonable, whereas variations were higher than 42% when the linear BEM and nonl CFD results were compared at 1.32 rad/s. The variation of the q-factor was high a wave frequency of 1.20 rad/s. The free surface wave profile colored by the vertical w particle displacement in the NWT is shown in Figure 12. Large differences in the surface wave elevation at the rear and front sides of rotors, phase lag among rotors the interaction of radiated waves allow for a rigorous mixing of the flow at the cent rotation (see Figure 12). The WEC rotor 2 demonstrated stronger and more distinct ing than WEC rotors 1 and 3. The violent liquid motion resulted in large slamming f and significant changes without PTO at ω = 1.20 rad/s and with PTO at 1.32 rad/s pared with other cases. However, compared with the nonlinear CFD, the linear results typically led to the overestimation of findings in most of the tested cases bec the violent interaction is neglected.  Table 3 presents the pitch RAO and q-factor comparison between the linear BEM and nonlinear CFD results. The pitch RAO of wave frequencies (1.02, 1.20, and 1.32 rad/s) without PTO shows the satisfactory overall consistency between the linear BEM and nonlinear CFD solutions, except for the significant variation in the case of 1.20 rad/s. The linear BEM results demonstrate higher destruction in the q-factor compared with the nonlinear CFD findings. The overall consistency of the pitch RAO with PTO was reasonable, whereas variations were higher than 42% when the linear BEM and nonlinear CFD results were compared at 1.32 rad/s. The variation of the q-factor was high at the wave frequency of 1.20 rad/s. The free surface wave profile colored by the vertical water particle displacement in the NWT is shown in Figure 12. Large differences in the free surface wave elevation at the rear and front sides of rotors, phase lag among rotors, and the interaction of radiated waves allow for a rigorous mixing of the flow at the center of rotation (see Figure 12). The WEC rotor 2 demonstrated stronger and more distinct mixing than WEC rotors 1 and 3. The violent liquid motion resulted in large slamming forces and significant changes without PTO at ω = 1.20 rad/s and with PTO at 1.32 rad/s compared with other cases. However, compared with the nonlinear CFD, the linear BEM results typically led to the overestimation of findings in most of the tested cases because the violent interaction is neglected.    Figure 13 illustrates the extracted power characteristics of isolated and multiple WEC rotors as a function of the wave frequency and wave-heading angle. The dotted line denotes the isolated WEC rotor and the contour represents multiple WEC rotors using linear BEM. The absorbed power presented in kW and based on Equation (3) corresponds to the optimal time-averaged extracted power of each wave frequency. Peak values are summarized in Table 4. The evident influence of the natural period on the extracted power in both isolated and multiple WEC rotors is illustrated in Figure 13. The maximum time-averaged optimum extracted power is distributed around the zero-heading angle close to the resonance frequency and continuously decreases as the wave-heading angle increases. Table 5 lists the extracted power at specific frequencies where the nonlinear CFD simulation results are compared with the linear BEM results. According to the linear BEM results, the maximum extracted power in both isolated and multiple WEC rotors were obtained at ω = 1.20 rad/s, followed by ω = 1.32 rad/s and ω = 1.09 rad/s. Nonlinear CFD simulations showed that the maximum extracted power was observed at a similar wave frequency of linear BEM results, followed by ω = 1.09 rad/s and ω = 1.32 rad/s.  Figure 13 illustrates the extracted power characteristics of isolated and multiple WEC rotors as a function of the wave frequency and wave-heading angle. The dotted line denotes the isolated WEC rotor and the contour represents multiple WEC rotors using linear BEM. The absorbed power presented in kW and based on Equation (3) corresponds to the optimal time-averaged extracted power of each wave frequency. Peak values are summarized in Table 4. The evident influence of the natural period on the extracted power in both isolated and multiple WEC rotors is illustrated in Figure 13. The maximum timeaveraged optimum extracted power is distributed around the zero-heading angle close to the resonance frequency and continuously decreases as the wave-heading angle increases. Table 5 lists the extracted power at specific frequencies where the nonlinear CFD simulation results are compared with the linear BEM results. According to the linear BEM results, the maximum extracted power in both isolated and multiple WEC rotors were obtained at ω = 1.20 rad/s, followed by ω = 1.32 rad/s and ω = 1.09 rad/s. Nonlinear CFD simulations showed that the maximum extracted power was observed at a similar wave frequency of linear BEM results, followed by ω = 1.09 rad/s and ω = 1.32 rad/s.  Figure 13 illustrates the extracted power characteristics of isolated and multiple WEC rotors as a function of the wave frequency and wave-heading angle. The dotted line denotes the isolated WEC rotor and the contour represents multiple WEC rotors using linear BEM. The absorbed power presented in kW and based on Equation (3) corresponds to the optimal time-averaged extracted power of each wave frequency. Peak values are summarized in Table 4. The evident influence of the natural period on the extracted power in both isolated and multiple WEC rotors is illustrated in Figure 13. The maximum time-averaged optimum extracted power is distributed around the zero-heading angle close to the resonance frequency and continuously decreases as the wave-heading angle increases. Table 5 lists the extracted power at specific frequencies where the nonlinear CFD simulation results are compared with the linear BEM results. According to the linear BEM results, the maximum extracted power in both isolated and multiple WEC rotors were obtained at ω = 1.20 rad/s, followed by ω = 1.32 rad/s and ω = 1.09 rad/s. Nonlinear CFD simulations showed that the maximum extracted power was observed at a similar wave frequency of linear BEM results, followed by ω = 1.09 rad/s and ω = 1.32 rad/s.    Note:

Conclusions
A preliminary assessment of multiple WEC rotors was conducted, and the effects of the wave-heading angle and wave frequency with and without PTO damping were quantified using linear BEM along with the isolated WEC rotor. Laboratory experiments were carried out to validate the adopted numerical models on a scale-down isolated WEC rotor. Nonlinear CFD simulations were carried out on prototype WEC for three different wave frequencies at a fixed wave-heading angle to verify and compare the results obtained from the linear BEM solution. The following conclusions can be drawn from the present study:

•
The range of tested wave frequencies showed that the small wave-heading angle of the isolated WEC rotor increases the RAO while peak amplitudes are slightly affected in multiple WEC rotors with the linear BEM.

•
The maximum reduction of the pitch RAO is focused around the peak frequency but insignificant elsewhere with the application of the PTO damping system.

•
The q-factor demonstrated a constructive interaction in the range of (0.5 < ω < 1.08 and β ≤ 60 • ) and (1.34 < ω < 1.82 and β < 30 • ) for without PTO and (0.52 < ω < 1.1 and β ≤ 60 • ) and (1.33 < ω < 1.77 and β < 30 • ) for with PTO but destructive otherwise with the linear BEM. • Investigation based on the prototype isolated and multiple WEC rotors, the pitch RAO, and the q-factor with without PTO at chosen wave frequencies demonstrated satisfactory overall consistency between the linear BEM and nonlinear CFD results, except for significant differences at 1.20 rad/s without PTO and 1.32 rad/s with PTO. These differences were qualified by nonlinear CFD simulations because the violent liquid motion at the center of the rotation leaves large slamming forces on the WEC rotor and leads to reduced and distinct responses from other wave frequencies.

•
The estimated maximum time-averaged optimal power is distributed close to zeroheading angles around the resonance frequency and continuously decreases as the wave-heading angle increases.

•
The linear BEM results showed overestimated extracted power compared with the nonlinear CFD due to the absence of nonlinear effects. Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.
Data Availability Statement: Data available in a publicly accessible repository.