Modelling of a Flow-Induced Oscillation, Two-Cylinder, Hydrokinetic Energy Converter Based on Experimental Data

: The VIVACE Converter consists of cylindrical oscillators in tandem subjected to transverse ﬂow-induced oscillations (FIOs) that can be improved by varying the system parameters for a given in-ﬂow velocity: damping, stiffness, and in-ﬂow center-to-center spacing. Compared to a single isolated cylinder, tandem cylinders can harness more hydrokinetic energy due to synergy in FIO. Experimental and numerical methods have been utilized to analyze the FIO and energy harnessing of VIVACE. A surrogate-based model of two tandem cylinders is developed to predict the power harvesting and corresponding efﬁciency by introducing a backpropagation neural network. It is then utilized to reduce excessive experimental or computational testing. The effects of spacing, damping, and stiffness on harvested power and efﬁciency of the established prediction-model are analyzed. At each selected ﬂow velocity, optimization results of power harvesting using the prediction-model are calculated under different combinations of damping and stiffness. The main conclusions are: (1) The surrogate model, built on extensive experimental data for tandem cylinders, can predict the cylinder oscillatory response accurately. (2) Increasing the damping ratio range from 0–0.24 to 0–0.30 is beneﬁcial for improving power efﬁciency, but has no signiﬁcant effect on power harvesting. (3) In galloping, a spacing ratio of 1.57 has the highest optimal harnessed power and efﬁciency compared with other spacing values. (4) Two tandem cylinders can harness 2.01–4.67 times the optimal power of an isolated cylinder. In addition, the former can achieve 1.46–4.01 times the efﬁciency of the latter. (5) The surrogate model is an efﬁcient predictive tool deﬁning parameters of the Converter for improved energy acquisition.


Introduction
The present trend towards clean renewable energy is very strong [1][2][3]. It is gaining momentum worldwide and it is unlikely to slow down despite the current low cost of fossil fuels. Consequently, the exploration and development of renewable energy sources are urgent. Marine energy, as a widely distributed sustainable energy source, with the characteristics of strong stability and predictability, and high-energy density, has been studied by a considerable number of scholars [4][5][6]. One of the major forms of marine energy is horizontal hydrokinetic energy in ocean currents, tides, and rivers [7]. In last decade, research on the utilization of ocean current energy for power generation is widespread [8,9], with the main power generation device being the underwater turbine [10,11]. Underwater turbines, which generate electricity by rotation of their impeller induced by the incoming flow, require a high flow velocity (2.0-2.5 m/s) to start. At most sea locations, however, the current velocity is low, with an average velocity of about 0.4-1.5 m/s and a maximum velocity generally not exceeding 2.5 m/s. Typical river flows are 1.0-1.5 m/s. the relationship of influencing factors and power have not been directly built, but it is solved by the second-order linear differential equation and power equation [33], meaning that the probability of errors will increase. Therefore, this paper used the programming software to write neural networks instead of the "Isight 5.0", and directly searched for the relationship between input variables and power. The BPNN and RBFNN both can build the nonlinear relationship between input variables and power. However, the structure of BPNN is simpler than that of RBFNN when solving the problems with the same precision requirements. Moreover, the number of hidden layer neurons of RBF neural network is much higher than that of BPNN when there are more training samples, which makes the complexity and computation of RBFNN increase greatly. Based on above two reasons, we selected the BPNN to establish the surrogate model of VIVACE Converter with two tandem cylinders.
In this paper, four different surrogate-based models, three for two tandem cylinders with different spacing and one for an isolated single cylinder are established by implementing the back propagation neural network method. This modeling implementation is based on sample training and testing of the experimental data. Based on the established surrogate model of the VIVACE Converter, this paper explores and analyzes the influences of harnessing damping ratio, spring stiffness, and spacing ratio on the harnessed power and the corresponding efficiency of the VIVACE Converter. Extending the selected range of harnessing damping ratio and spring stiffness is conducive to enhancing the oscillating responses, especially for the power efficiency. In addition, a tandem cylinder with reasonable combinations of various input parameters can harness more power at higher efficiency than that of the isolated single-cylinder due to the positive synergy between cylinders in close proximity.

Experimental Facility: LTFSW Channel
All tests were conducted in the Low Turbulence Free Surface Water (LTFSW) Channel in the MRELab of the University of Michigan. In the LTFSW Channel, 10,000 gallons (37,854 lt) of freshwater are recirculated at speed up to 1.5 m/s realized by using an impeller powered by a 20 hp induction motor. The test section is 2.44 m long and 1.0 m wide and is made of transparent plexiglass, which allows for visualization of both cylinders in FIO using argon lasers of 5 W and adding aluminum oxide particles of 100 µm. Figure 1 shows the VIVACE Converter with a two cylinder-oscillator system in the Channel. The water depth of the test section is 1.17 m and the maximum flow velocity is kept lower than 1.35 m/s to prevent damage due to the galloping instability.

Cylinders with Distributed Roughness
Passive turbulence control (PTC) [23] consists of two roughness strips on the surface of the smooth cylinder. It was introduced in the VIVACE Converter and extensively studied in the MRELab [19,20,24,27] in order to harvest more marine hydrokinetic energy from the steady flow using FIO of cylinders. As shown in Figure 2, the total height of the PTC consists of double-sided tape at the bottom, sandpaper in the middle, and coarse sand on top. In addition, the total height of PTC is set to be about equal to the boundary layer thickness of the fluid around the cylinder. Attaching the PTC on the surface of the cylinder results in geometric asymmetry of the cross-section of the oscillator inducing the galloping, which is characterized by unstable oscillation and high energy acquisition. A secondary benefit is that PTC fixes the separation point in VIV, which for a smooth cylinder oscillates by about 5 • to 10 • . This seemingly minute detail is very important and difficult to predict numerically. A negative effect of PTC is that the surface roughness reduces the vortex strength in VIV and the response amplitude. This is corrected by introducing smooth PTC in the form of step surface elements like simple tape.

Vck System
Physical springs and dampers were applied to the early VIVACE Converter models in the MRELab. For systematic experiments, however, to test different combinations of spring stiffness and damping ratio, replacing the physical springs and dampers requires extensive design, fabrication, calibration, and alignment time. Besides that, the nonlinear viscous damping of a Converter is not easily modeled by a physical damper. To compensate for the above two challenges, the virtual spring-damping (Vck) system was developed and improved by the MRELab and was used to emulate physical springs and dampers. In 2011, Lee et al. [26] developed the first generation of Vck. In 2013, Sun et al. [39] developed the embedded second-generation Vck system, which reduced the controller response time by three orders of magnitude by eliminating LabView. It used digital instead of analogue signal and the damping model simplified significantly. The Vck system is a controller with a feedback loop that uses a servo motor encoder to track the position and speed of the cylindrical oscillator and to provide the necessary torque to simulate spring stiffness, as well as the required linear viscous damping. Further, it can emulate any nonlinear spring and damper function and its parametric values.

Mathematical Model of Harnessed Power and Harnessing Efficiency
To reasonably and systematically analyze the performance of the VIVACE Converter, the harnessed power P harness and harnessed efficiency η harness are proposed and introduced as two important evaluation metrics. As expressed in Equation (1), the harnessed power is positively correlated with the oscillating system mass, harnessing damping ratio ζ harness, and natural frequency of oscillator in water f n,water , respectively. A detailed description of ζ harness and f n,water is provided in [23].
where A and f osc are, respectively, the oscillating amplitude and frequency of the cylinder. The harnessed power can be calculated by averaging the instantaneous harnessed power. For a linear system, since Equation (1) depends explicitly on the oscillating mass of the system, m sys , power can be calculated via the recorded amplitude and the frequency of oscillation. For either linear or nonlinear damping, the instantaneous power can be written as follows: P harness (t) = c harness . y 2 (t) (2) where . y(t) is the oscillation speed of the cylinder as measured by the Vck system. Equation (2) does not depend on the specific way the added mass is calculated or even the value of the added mass explicitly; it cancels out, thus, eliminating the need to deal with modeling of the added mass [13].
The power of fluid-flow can be calculated as: where A max represents the maximum amplitude, ρ is the fluid density, U is the incoming flow velocity, and D and L are the diameter and length of the cylinder, respectively. The Betz Limit is the theoretical maximum efficiency of power extracted from the flowing fluid. Therefore, the converter efficiency can be defined as:

Modeling of the Converter
In this section, four different surrogate models of the VIVACE Converter are established based on the BP neural network method by using the experimental data measured for flow velocity, spring stiffness, and harnessing damping ratio values uniformly distributed within the specified range [23,35]. These four models correspond to the isolated single-cylinder device and three double-cylinder devices with tandem spacing ratios of L/D = 1.57, 2.00, and 2.57. The surrogate models built in this section can construct the complex nonlinear mathematical relationship between input variables (U, spring stiffness K and ζ harness ) and objective variable (P harness or η harness ). After the establishment of the prediction model, the power harvesting of the VIVACE Converter can be predicted by inputting appropriate and reasonable values of the input variables, which are not limited to the tested experimental data. The principle of BP neural network is described in Section 3.1. The construction and validation of surrogate models are presented in Section 3.2.

Principle of BP Neural Network
The BP neural network is a multi-layer feedforward network trained by an error inverse propagation algorithm. It is one of the most widely used neural network models [36,37]. As illustrated in Figure 3, a BP neural network is mainly composed of the input layer, middle layer, and output layer. The middle structure can be designed as single layer or multiple hidden layers. The node numbers of the neural network of the above three layers is set as p (1 ≤ i ≤ p), s (1 ≤ j ≤ s), r (1 ≤ k ≤ r), respectively. The neurons between the adjacent layers are connected in the form of weights and biases. In addition, the weights or biases of the input layer to the hidden layer and the hidden layer to the output layer are, respectively, w ij and w jk or a j and b k . The BP neural network method consists of the forward-propagation of input signal and back-propagation of error. For the forward-propagation of the input signal, the neurons of the input layer receive the input signal x i and transmit it to the neurons of the hidden layer. The hidden layer is the internal information processing layer and is responsible for information transformation. The output signal G j of the hidden layer can be expressed as in Equation (5): where f is the nonlinear transfer function, also called the activation function. Then, the hidden layer transmits the signal G j to the output layer. The output layer completes the signal processing and outputs it. The output signal Q k is calculated by the Equation (6) as follows: The error E that exists between the expected output D k (i.e., experimental data) and predicted output Q k can be expressed as follows: In the process of the back-propagation, the error will continuously correct the network connection weights and biases of two adjacent layers by backward propagation along with the output layer, hidden layer, and input layer, respectively. The process will be completed when the output error E satisfies the step requirement or reaches the maximum training time of the neural network.

Construction and Validation of Surrogate Model
The correctness of the construction of the surrogate model is largely determined by the selection of input parameters. For the VIVACE Converter, the most easily changeable variable parameters are the incoming flow velocity, harnessing damping ratio, and spring stiffness, which would directly affect the acquisition of harnessed power and the formation of power efficiency. Due to the reason that there are only three groups of optional spacing ratios, the spacing ratio (L/D = 1.57, 2.00, 2.75) is not included in the input variables. In addition, the output objective is harnessed power measured experimentally in reference [23]. The modeling of power efficiency can be realized by utilizing the same approach as used for the power harvesting model. Nevertheless, in the following discussion, only the surrogate model of harnessed power is considered. The steps for the construction and validation of the surrogate model in order of occurrence are: determination of neural network structure, sample training, and sample testing.

Range of Input Variables
The selection range of input experimental design parameters is constrained by the geometry of the model converter ( Figure 1), engineering applications, and mechanical properties of cylinder oscillators.
Incoming flow velocity (U): The incoming flow velocity is selected in the range from 0.39 to 1.31 m/s with an interval of 0.04 m/s, which covers to the actual ocean conditions of water current [40] and the experimental condition of Sun et al. [23,35].
Harnessing damping ratio (ζ harness ): In order to trigger FIO easier in VIV branches and harvest higher power in galloping, the values of ζ harness chosen vary from 0.00 to 0.24 in intervals of 0.04.
Spring stiffness (K): The five spring stiffness values are evenly chosen from 400 N/m to 1200 N/m in increments of 200 N/m. Low stiffness with low damping induce high oscillating amplitude, resulting in reaching the amplitude limit of the water channel. Further, the highest value of K is chosen depending on the strength limit of the timing belts of the Vck controller system. The input experiment-design parameters and its constraint range for four different surrogate models are uniform and summarized in Table 1.

Neural Network Structure
The number of input layer nodes and output layer nodes for all the four different models are selected as 3 and 1, respectively, according to the number of input and output variables. Accordingly, in order to avoid the phenomenon of excessive training time caused by multiple hidden layers, a single hidden layer is set in this paper. Additionally, for the process from input layer to hidden layer and from hidden layer to output layer, the selected activation function was set as 'purline' and 'tansig', respectively.
The determination of the number of hidden layer nodes is of pivotal importance for establishing a correct prediction model by using the BP neural network. Too few nodes result in insufficient training of the BP neural network, meaning that the prediction error of the neural network is unreasonably large. On the other hand, too many nodes will increase the learning time of neural network, leading to the phenomenon of "overfitting". That is to say, even if the prediction accuracy of sample training is high, the prediction error of other samples may be large. Selecting the experimental data for spacing ratio L/D = 2.00 as the input or output parameter of BP neural network and introducing the root-mean-square error (RMSE, between predicted data and actual experimental data) as an evaluative criterion, the relationship between the RMSE of BP neural network and the number of hidden layer node is established. The corresponding curve of RMSE vs. the number of hidden layer nodes for the isolated single cylinder was added in Figure 4. Figure 4 shows that the minimum RMSE values of two tandem cylinders and isolated cylinder appear when l = 30 and l = 24, respectively. Accordingly, the number of the hidden layer nodes is selected in this paper. To analyze the influence of the spacing ratio on FIO of two tandem cylinders, the node number of l = 30 is selected for all three spacing-ratio cases.

Sample Training
For sample training, the total sample set of experimental data (i.e., different combinations of U, K, and ζ harness ) for three double-cylinder devices and for isolated single cylinder are 816 and 809, respectively. The relative parameters including the number of hidden layer nodes and activation function were explored and determined in Section 3.2.2. The correlation coefficient R 2 [41], which can be expressed as Equation (8), was introduced to evaluate the modeling quality and fitting result of the surrogate model. The closer R 2 is to 1, the higher the modeling quality is: As shown in Equation (8), SSE represents the sum of the squares of errors between the predicted output Q k and the actual experimental data, and SST is the sum of the squares of the difference between the actual experimental data and its mean. As shown in Figure 5, it can be observed that the values of R 2 for the four different models are all higher than 0.97 and close to 1. Moreover, the RMSE values are far smaller than the corresponding mean value of the experimental data M. Therefore, the modeling quality and fitting accuracy of effectiveness are both high.

Sample Testing
For sample testing, a random combination of K and ζ harness (K = 1200 N/m, ζ harness = 0.16) was selected to calculate the variation trend of power acquisition with flow velocity so as to verify the predictive ability and modeling accuracy of the four different surrogate models developed in this paper. It should be noted that experimental sample data for this pair of K and ζ harness have not been used for the sample training. Although the experimental data used for testing the accuracy of modeling did not reach 20% of the total experimental sample size, the testing samples selected is extensive and includes all flow velocity values. The flow velocity varies from 0.39 to 1.31 m/s in 0.04 m/s increments for a total of 21-23 velocity values. The purpose of doing so is that this paper mainly studies the power variation trend in the range of fluid-induced vibration at each flow velocity. It is essential to verify the energy at each flow velocity. That is, there's a rational in selecting the "random samples". Based on experience and reference [33] we concluded that the testing samples selected in this paper are reasonable and their number is adequate. A second abscissa axis, showing the reduced velocity U*, is provided in Figure 6 because it makes the results easier to interpret for a reader with background in VIV. For constant spring-stiffness K, the reduced velocity U* defined by Equation (9) is proportional to the flow velocity: where f n,water is the natural frequency of the oscillator in quiescent water. Figure 6 shows the comparison between the harnessed power as measured experimentally and as predicted by the surrogate model. It can be observed that, in these four cases, the corresponding power curves of prediction and experiment fit well. In addition, except for the VIV to galloping transition region, the relative error between prediction and experiment is low. Combined with the results in Figure 6, it can be concluded that the BP neural network method utilized in this paper is suitable for modeling of two tandem cylinders and isolated single cylinder in flow-induced oscillation. Based on the above analysis of sample training and sample testing, Figure 6 demonstrates that the present BP surrogate model can predict the corresponding harnessed power accurately within the range of the input parameters.    Figure 7, for tandem cylinders with L/D = 2.00, when the flow velocity is 0.584 m/s, the power is obviously not a monotonic function and has several peak values under the different combinations of K and ζ harness . In addition, with increasing flow velocity, the overall value of harnessed power increases but with some fluctuations. These two obvious phenomena are consistent with the experimental results [23], verifying the accuracy of the modeling. Besides that, the majority of P harness of tandem cylinders is higher than 2.0 times the corresponding P harness of the single-cylinder, indicating that in most instances, the interaction between two tandem oscillators has a positive influence on the improvement of harnessed power.

Prediction, and Discussion of Results
In Section 3, four predictive models of the VIVACE Converter with different combinations of oscillators were established and validated: One for an isolated single cylinder, and three for two tandem cylinders with spacing ratio L/D = 1.57, 2.00, 2.57. In this section, these four surrogate-based models are respectively utilized to predict and compare the harnessed power and power efficiency under different combinations of harnessing damping ratio, spring stiffness, and incoming flow velocity. It should be pointed out that the results shown in Section 4 are all the predictions of the established predictive models. Figures 9-11 show the distribution of predicted P harness or η harness with incoming flow velocity as the independent variable (0.39 m/s < U < 1.31 m/s), for two tandem cylinder with spacing ratio L/D = 1.57, 2.00, 2.57, respectively. The other parameters are ζ harness and K. In addition, the corresponding information for an isolated single cylinder can be found in Figure 12. Considering the structural feasibility and actual engineering situation, ζ harness selected has a wider range of 0.00-0.30 rather than 0.00-0.24. According to the experimental results in [23], when the spacing ratio is L/D = 1.57, 2.00, and 2.57, the corresponding optimal P harness appears in the cases of K = 1200 N/m, 800 N/m, and 400 N/m, respectively. Likewise, the optimal spring stiffness for an isolated single cylinder is K = 1200 N/m. Thus, in this section, these four stiffness values are selected respectively as the fixed values for the four different surrogate models established above. This reduces the variable and controls parameters to better analyze the impact of ζ harness on P harness . The following results can be observed: Nevertheless, the same phenomenon that occurs with a single-cylinder is that for L/D = 2.57, the dropping tendency of P harness becomes obvious in the transition region.

Effect of Harnessing Damping Ratio on Power and Efficiency
In galloping, for all spacing ratios, P harness continually increases with the increase of velocity. The above results are highly consistent with the experimental results in [23], indicating that these surrogate models have high accuracy. It should be pointed out that the area covered in gray in Figures 9-12 represents the transition region. (b) For L/D = 1.57, Figure 9a shows that, with the increase of damping within the range of 0.20-0.30, the variation of P harness is negligible. Figure 9b illustrates that, in the VIV range, increasing the damping from 0.20 to 0.30 will translate the curve of η harness to higher velocity and is followed by a significant decrease in η harness . The above phenomenon is owing to the negative interaction between two tandem cylinders in the sense that the higher damping limits the oscillating amplitude of the downstream cylinder, which is then shielded by the upstream cylinder. As shown in Figure 9b, the optimal η harness occurs at the upper branch of VIV, and as the damping ratio increases from 0.2 to 0.3, the optimal efficiency increases as well. For instance, the optimal efficiency changes from 40.82% ( Figure 10b shows that, except for the initial VIV branch, in all selected FIO regions, as ζ harness increases within the range 0.04 < ζ harness < 0.30, the corresponding η harness increases. As the harnessing damping ratio increases from 0.20 to 0.30, the increment rate of harnessing efficiency will reach up to 59.52% at U = 0.71 m/s, and the maximum harnessing efficiency is 68.78% at U = 0.63 m/s.  Figure 11a. Here, the optimal harnessed power changes from 31.37 W (K = 400 N/m, ζ harness = 0.20) to 27.54 W (K = 400 N/m, ζ harness = 0.30) at U = 1.31 m/s. In contrast, it is found that increasing the harnessing damping ratio from 0.20 to 0.30 is followed by an increase of power efficiency in the whole FIO region except at the end region of galloping. For instance, the efficiency changes from 52.09% (K = 400 N/m, ζ harness = 0.20) to 64.92% (K = 400 N/m, ζ harness = 0.30) at U = 0.43 m/s. Due to the low spring stiffness (K = 400 N/m), the initial velocity of power harvesting is lower than 0.39m/s. Besides that, the amplitude at low velocity has a relative low value, resulting in high efficiency. As mentioned above, for two tandem cylinders with spacing ratio L/D = 2.57 and K = 400 N/m, it can be observed that increasing the damping ratio from 0.20 to 0.30 will decrease the harnessed power; however, it enhances the power efficiency of the VIVACE Converter. (e) For an isolated cylinder in general, with the increase of spacing ratio, the interaction between the upstream cylinder and downstream cylinder becomes weaker. Therefore, the changes of curves of harnessed power between the spacing ratio of L/D = 2.57 and the isolated single-cylinder cases are similar. However, in this section, the selected fixed stiffness for isolated cylinder and two cylinders with tandem spacing ratio L/D = 2.57 are different. Thus, a large difference in oscillating response exists between these two cases. As we can see in Figure 12a, except for the initial galloping region, the value of P harness decreases with the increase of harnessing damping ratio (0.2 < ζ harness < 0.30). Moreover, as ζ harness increases, the reduction between two adjacent damping ratios increases. Because of the large decrease in P harness , the corresponding η harness decreases as well when the damping ratio increases. (f) Comparing the power efficiency between the four cases (spacing ratio of L/D = 1.57, 2.00, 2.57 and isolated single cylinder), it is observed that the efficiency of tandem cylinders is much larger than that of the isolated cylinder. For L/D = 2.00, the optimal efficiency for K = 800 N/m and ζ harness = 0.04-0.30, reaches 1.91 times the hydrokinetic power efficiency of one isolated cylinder. As the spacing ratio increases from 1.57 to 2.00 to 2.57, the optimal efficiency increases from 50.51% to 68.78% to 64.92%, indicating that the smaller spacing of 1.57 induces a negative effect on the power efficiency due to the strong interaction between the two cylinders.     Figure 13, shows the two-dimensional diagrams of prediction-results of P harness with incoming flow velocity (0.39 m/s < U < 1.31 m/s), for the three spacing values of L/D = 1.57, 2.00, 2.57 and the isolated single-cylinder case, respectively. K is the parameter in Figure 13. In order to prevent excessive load on the timing belt of the Vck system, the spring stiffness was allowed to expand from 400-1200 N/m to 400-1600 N/m only. As demonstrated by Sun et al. [23], the harnessing damping ratio of 0.20 is conducive to obtaining the highest values of P harness for the cases of these three spacing ratios. Therefore, ζ harness = 0.20 is chosen as a fixed value in order to reduce the variable and control parameter so as to better analyze the impact of K on P harness . The following results can be observed:

Effect of Spring Stiffness on Power and Efficiency
(a) As shown in Figure 13a-d, in the VIV region, with the increase of flow velocity, the value of P harness increases. An important finding is that, for the spacing ratio L/D = 1.57 and 2.00, the transition region is not apparent due to the close spacing between the cylinders. Contrary, the spacing ratio L/D = 2.57 and the isolated singlecylinder have a relatively obvious transition. Therefore, the correctness of surrogate models has verified again.  Figure 13b, considering the spring stiffness, all the curves of P harness show a similar tendency. Moreover, in galloping, two local peaks of P harness appear at the K = 800 N/m and 1000 N/m, respectively. This phenomenon of K = 800 N/m is in good agreement with the experimental data [23]. Almost in all cases, in the VIV upper branch and the galloping region, the harnessed power shows the downward trend with increase of K in the range 1200 N/m < K < 1600 N/m. Nevertheless, in galloping, the decrement of P harness between adjacent stiffness values is small. In consequence, for a given harnessing damping ratio of 0.20, increasing K from 1200 N/m to 1600 N/m is adverse to the value of P harness .  Figure 13d, an isolated single-cylinder with high spring stiffness cannot induce VIV at low flow velocity due to the high natural frequency of the oscillator in quiescent water, f n,water . This is due to the strong restoring force provided by the hard spring. For the same reason, as K increases from 600 N/m to 1600 N/m, the initial flow velocity of power harvesting increases as well. This increase is expected to follow the square root of K as in f n,water. Of course, the curve of P harness of the isolated single cylinder will shift to higher velocity. Therefore, in the VIV region, expanding the selected range of K from 400-1200 N/m to 400-1600 N/m has no effect on the increase of P harness . On the contrary, in the transition region and initial galloping, expanding the selected range of K is beneficial for enhancing the power harvesting. The range of synchronization in VIV shifts to the right proportionally to the square root of K while galloping initiation is independent of K. When U = 0.83-1.11 m/s, as K increases from 1200 N/m to 1600 N/m, P harness increases. In the transition region, the optimal harnessed power can reach up to 10

Optimal Power and Efficiency Envelope for Two Cylinders in Tandem
To further analyze the performance of this Current Energy Converter (CEC) consisting of two cylinders in tandem, the above established four surrogate models are used to find the optimal harnessed power and optimal power efficiency. This is done for each available flow velocity: i.e., U = 0.39-1.31 m/s, in increments of 0.04 m/s. It should be pointed out that the optimal values of P harness and η harness can be found by introducing the traditional method of grid searching based on the modeling of power and efficiency in a series of spring stiffness and harnessing damping ratio values. The spring stiffness and harnessing damping ratio are, respectively, selected in the ranges of K = 400-1200 N/m in 10N/m intervals and ζ harness = 0.00-0.30 in 0.005 increments. Selecting too large an interval would reduce the accuracy of the optimal value. On the other hand, selecting too small an interval would result in unnecessary calculations as the power and efficiency of the original experimental data were not measured in such small increments. The optimal power and efficiency for each velocity are shown in Figure 14a,b, respectively. The corresponding optimal parameters are presented in Appendix A, Tables A1 and A2.  Table A1, the combination of lower stiffness with lower harnessing damping ratio stimulates the onset of VIV for lower velocity U < 0.39 m/s, which is consistent with the conclusion in Ref. [23]. Additionally, for spacing L/D = 2.00 and 2.57, increasing the flow velocity is accompanied by an increase in power efficiency. The efficiency exhibits a local maximum after which it decreases as shown in Figure 14b. (e) To better analyze the synergy effect of the two tandem cylinders, the optimal harnessed power and efficiency of the isolated single cylinder are also modeled, and are added in Figure 14a,b, respectively. The corresponding optimal parameters for an isolated cylinder are presented in Appendix A, Table A3. Specifically, the green line in Figure 14 represents the maximum optimal value at each velocity for three cases of the tandem cylinder. As shown in Figure 14a, it can be observed that at almost every selected flow velocity, the optimal harnessed power of two tandem cylinders has a higher value than that of the isolated cylinder. More concretely, the two cylinders in tandem can harness 2.01-4.67 times the power of the isolated cylinder at the selected range of flow velocity except for the upper branch of VIV. Further, the two tandem cylinders can achieve 1.46-4.01 times the efficiency of the isolated one.
Based on the analysis in Section 4.1, in order to further verify that the wider range of damping ratio (0.00 < ζ harness < 0.30) has better performance than that of the original range (0.00 < ζ harness < 0.24) [23], a similar surrogate model as the one developed in Section 4.3 is established by replacing the wider range of ζ harness by the narrower one. Figure 15a shows that extending the range from 0.00-0.24 to 0.00-0.30 has little impact on the optimal harnessed power for each flow velocity, except for the spacing L/D = 2.00. For L/D = 2.00, the maximum increment will reach up to 5.30 W at U = 0.95 m/s. However, this operation for ζ harness obviously increases the optimal power efficiency, especially in the upper branch of the VIV region and a transition region between VIV and galloping, as shown in Figure 15b. Therefore, extending the harnessing damping ratio from 0.00-0.24 to 0.00-0.30 has a positive influence on optimal power harvesting and optimal power efficiency.

Conclusions
The harnessed power and converting efficiency of the VIVACE hydrokinetic energy Converter with two tandem cylinders were modeled by a surrogate-based model established by the BP neural network method. The design parameters are flow velocity, harnessing damping ratio, spring stiffness, and spacing ratio. Modeling is based on the experimental data from tests conducted in the LTFSW channel of the MRELab of the University of Michigan. The major conclusions can be summarized as follows: (1) The structure of BPNN is simpler than that of RBFNN when solving the problems with the same precision requirements, and the number of hidden layer neurons of RBF neural network is much higher than that of BPNN when there are more training samples, which makes the complexity and computation of RBFNN increase greatly. Therefore, this paper selected the BPNN to establish the surrogate model of the VIVACE Converter.
Due to synergy between two tandem cylinders, the VIVACE Converter with two tandem cylinders has a better performance of power harvesting compared to that of the isolated single cylinder. Analysis of the optimal harnessed power and efficiency shows that the two tandem cylinders can harness 2.01-4.67 times the power of the isolated cylinder. In addition, two tandem cylinders can achieve 1.46-4.01 times the efficiency of the isolated one.
The results provide valuable insight into the parametric dependence of power and efficiency of a two-oscillator converter on spacing, spring-stiffness, damping ratio, and flow velocity. What may be even more valuable though, is the potential implementation of the developed model in engineering applications. That is to say, at each flow velocity, selecting appropriate optimal spacing ratio, spring stiffness and harnessing damping ratio can induce the optimal harvested power or efficiency in the actual current. Specifically: (a) During field-tests or in commercial use of a two-cylinder VIVACE Converter, flow speed changes. To maintain near-optimal performance, parameters should be adjusted. The derived model can provide real-time guidance to that effect. (b) The model can be tested against field-data and be modified or expanded accordingly. (c) This Converter has been studied experimentally and with field-tests for over a decade using many different parameters and various size oscillators. One of the most powerful ways of increasing its harnessed power or its harnessing efficiency is an adaptive change of parameters. For example, in [18] the spring stiffness was adjusted based on response amplitude. In [25], the damping was adjusted based on the speed of the oscillating cylinder. (d) One of the greatest advantages of the VIVACE Converter is that it is fish-friendly because it operates on alternating lift like fish rather than steady lift like wings. Nevertheless, until recently [6], there was no physics-based mathematical model of the alternating lift force in VIV or galloping, which is required for proper control. The model developed in this paper serves as a guide for a numerical control model. Institutional Review Board Statement: Not applicable.

Informed Consent Statement:
Informed consent was obtained from all subjects involved in the study.

Conflicts of Interest:
The authors declare no conflict of interest.
Appendix A Table A1. Parameters for optimal harnessed power for two tandem cylinders with spacing: L/D = 1.57, 2.00 and 2.57.   Table A2. Parameters for optimal converted efficiency for two tandem cylinders with spacing: L/D = 1.57, 2.00 and 2.57.   Table A3. Parameters of optimal harnessed power and efficiency for single cylinder.

Parameters for Optimal P harness
Parameters for Optimal η harness