Capture Power Prediction of the Frustum of a Cone Shaped Floating Body Based on BP Neural Network

How to improve the power generation of wave energy converters (WEC) has become one of the main research objectives in wave energy field. This paper illustrates a framework on the use of back propagation (BP) neural network in predicting capture power of the frustum of a cone shaped floating body. Mathematical model of single floating body is derived, and radius, semi-vertical angle, mass, submergence depth, power take-off (PTO) damping coefficient, and stiffness coefficient are identified as key variables. Commercial software ANSYS-AQWA is used for numerical simulations to obtain hydrodynamic parameters, and then capture power is calculated by these parameters. A database containing 100 samples is established by Latin hypercube sampling (LHS) method, and a simple feature study is conducted. A BP neural network model with high accuracy is designed and trained for predictions based on built database. The results show that forecasting results and desired outputs are in great agreement with error percentage not greater than 4%, correlation coefficient (CC) greater than 0.9, P value close to 1, and root mean square error (RMSE) less than 139 W. The proposed method provides a guideline for designers to identify basic parameters of the floating body and system damping coefficient.


Introduction
Wave energy converters (WEC), a new type of energy extractor with little pollution, are expected to be a reliable alternative to the current generation method. There are two stages for an oscillating body WEC transforming wave energy into other forms of energy like electricity. A floating body is firstly required to capture the wave energy induced by a wave's motion. Then the moving floating body drives a generator to generate power. An intact oscillating body WEC system is generally composed of a moving floating body, a power take-off (PTO) system, and an anchor chain, etc. At present, the conversion efficiency of WECs is relatively low, so the main research objective is to improve the power generation of a specific device.
One method is to design a different floating body's shape, and the shape is usually irregular curved surface. McCabe [1] researched the optimization of the shape of a wave energy collector to improve energy extraction by genetic algorithms, and a benchmark collector shape was identified. Colby [2] used evolutionary algorithms to optimize the ballast geometry and achieved 84% improvement in power output. Fang [3] designed a mass-adjustable float, and a new optimization calculation method was proposed. Multifreedom buoys have been also proposed in [4][5][6]. They can translate or rotate in more than one freedom, so more wave energy can be captured. Another means is to design an innovative PTO system. Reabroy [7] proposed a novel floating device integrated with a Under three assumptions, the following forces act on the floating body namic forces (excitation and radiation force); hydrostatic buoyancy; PTO dam rigid restoring force. According to the theory of fluid mechanics and Newto governing equation of the floating body can be expressed as follows: where M0 represents mass; z(t) represents the heave displacement; fE represent force; fR represents radiation force; fS represents hydrostatic buoyancy; fPTO repr damping force; fK represents rigid restoring force.
The excitation force imparts on the floating body by the incoming wave. It mation of the Froude-Krylov force fFK and the diffraction force fD, so it can be follows: The radiation force is induced by the floating body's motion and can be d into an added mass term and a radiation damping term [25], so it can be ex follows: Under three assumptions, the following forces act on the floating body: hydrodynamic forces (excitation and radiation force); hydrostatic buoyancy; PTO damping force; rigid restoring force. According to the theory of fluid mechanics and Newton's law, the governing equation of the floating body can be expressed as follows: where M 0 represents mass; z(t) represents the heave displacement; f E represents excitation force; f R represents radiation force; f S represents hydrostatic buoyancy; f PTO represents PTO damping force; f K represents rigid restoring force. The excitation force imparts on the floating body by the incoming wave. It is the summation of the Froude-Krylov force f FK and the diffraction force f D , so it can be written as follows: The radiation force is induced by the floating body's motion and can be decomposed into an added mass term and a radiation damping term [25], so it can be expressed as follows: where A M and B C are the added mass and radiation damping in the vertical direction, respectively.
The hydrostatic buoyancy, induced by seawater static pressure, is the resultant force of gravity and buoyancy. It is a force that restores the structure to hydrostatic equilibrium and is linear with the heave displacement of the floating body. It can be written as: where ρ is seawater density; g is acceleration of gravity; A W is water cross area of the floating body. In this paper, the value of PTO damping is relatively large, and the heave displacement of the floating body is small. As a result, it is assumed that the water cross area of the floating body does not change. It is the section where the water line is located when the floating body is in still water. Therefore, the hydrostatic buoyancy can be expressed as: where D is the diameter of the floating body. The energy conversion system can be simplified to a linear spring damping system, so the PTO damping force is where c is the damping coefficient of the PTO system. The rigid restoring force is proportional to the heave displacement, and it can be written as follows: where k is stiffness coefficient. Reformulate Equation (1) through Equations (2), (3), (5)- (7): ..
Apply Fourier transform to Equation (8) and obtain another governing equation in the frequency domain. It is where ω is the wave frequency; j is imaginary unit; Z(ω) and F E (ω) are functions of displacement and excitation force in the frequency domain, respectively. In the frequency domain, the excitation force can be expressed by the product of the unit excitation force and the incident wave amplitude [35]. It is Equation (9) can be rewritten as follows: Formula (11) is a typical damped and forced vibration equation, so the natural frequency and damping factor can be expressed as below: From Equations (12) and (13), the natural frequency and damping factor of a given WEC change over added mass and added damping.
According to Equation (11), the heave response in the frequency domain is The average power in one wave period, captured by the floating body with heave motion, can be written as the product of damping force and vertical velocity. The work done by damping force is the energy absorbed by the floating body, so the mean capture power is The mean capture power reaches the maximum when the following conditions are met.
This stiffness and damping are called the best stiffness and the best damping, respectively. When B c > 0, the natural frequency, damping factor, and displacement are The max capture power is

Latin Hypercube Sampling
The sampling method is of great importance in experimental design. A good sampling method can result in more reasonable sample distribution, leading to a better model with higher accuracy. In this paper, a Latin hypercube sampling (LHS) method is utilized to generate sample points. Different from random sampling, LHS has a high efficiency of space filling by maximizing the stratification of each edge distribution, which improves the uniformity.
According to Equation (15), the factors that determine the capture power under given wave conditions are PTO damping coefficient c, system stiffness coefficient k, wave exciting force F E , float mass m, added mass A M , and added damping B C . Added mass, added damping, and wave exciting force are related to the geometry and submergence depth of the floating body. The geometric features of the floating body depend on radius R, semivertical angle α, and mass m. As a result, four main geometric parameters, including radius R, semi-vertical angle α, mass m, and submergence depth d, plus two system parameters, PTO damping coefficient c and stiffness coefficient k, are selected as key variables that affect the capture power. The sample space of six key variables are defined as follows: A database covering 100 sample points is established (see in Appendix A) and scatter diagrams of these samples are shown in Figure 2.
ing force FE, float mass m, added mass AM, and added damping BC. Added mass, added damping, and wave exciting force are related to the geometry and submergence depth of the floating body. The geometric features of the floating body depend on radius R, semivertical angle α, and mass m. As a result, four main geometric parameters, including radius R, semi-vertical angle α, mass m, and submergence depth d, plus two system parameters, PTO damping coefficient c and stiffness coefficient k, are selected as key variables that affect the capture power.
The sample space of six key variables are defined as follows: [2,3] [1. 5,3] [7000,8000] [5,25] [10000, 30000] A database covering 100 sample points is established (see in Appendix A) and scatter diagrams of these samples are shown in Figure 2. In Figure 2, each variable fills the whole sample space and the standard deviation of the value is small. It can reflect the relationship between the factor and the response in the six spaces. In Figure 2, each variable fills the whole sample space and the standard deviation of the value is small. It can reflect the relationship between the factor and the response in the six spaces.

Feature Study
Suitable feature study on the data set can give an insight to the correlation between the inputs and output. Pearson correlation analysis is conducted in this section to identify the correlation between six key variables and the capture power. Figure 3 shows the correlation coefficients in different wave situations. In this heatmap, a negative value means a negative correlation, and a positive value means a positive correlation. A large absolute value means a strong correlation.
(e) (f) In Figure 2, each variable fills the whole sample space and the standard deviation of the value is small. It can reflect the relationship between the factor and the response in the six spaces.

Feature Study
Suitable feature study on the data set can give an insight to the correlation between the inputs and output. Pearson correlation analysis is conducted in this section to identify the correlation between six key variables and the capture power. Figure 3 shows the correlation coefficients in different wave situations. In this heatmap, a negative value means a negative correlation, and a positive value means a positive correlation. A large absolute value means a strong correlation.  In general, radius, submergence depth, and damping show a strong correlation, while semi-vertical angle, mass, and stiffness behave a weak correlation. Besides, the correlation is different at different wave frequency. When the wave frequency is 0.53 and 0.81 rad/s, semi-vertical angle shows a negative correlation, while a positive correlation comes up at other frequencies. The similar situation happens on stiffness. Mass and damping behave a positive correlation, while radius and submergence depth show a negative correlation all the time.

Simulation Scheme
The structural schematic of the cone shaped floating body investigated in this study is shown in Figure 4.
The height of the cylinder part above the waterline is a constant, 0.5 m. In ANSYS Design Modeler, the 3D geometry with given parameters is created.
In this paper, ANSYS-AQWA, a commercial computation software based on potential flow theory, is utilized to calculate hydrodynamic parameters. The simulation process, including numerical modeling, parameters setting, mesh generation, and data post-processing, can be conducted in the graphical interface directly. The basic simulation steps for each sample are as follows: 1.
The moment of inertia and center of mass of the floating body are calculated in Static Structural module;

2.
Set solution environment in hydrodynamic diffusion module. The water line is at z = 0, the seawater depth is 200 m, and the surface area are 100 m × 100 m. Details of the point mass, additional damping, and additional hydrostatic stiffness are set according to the results obtained in Static Structural module and parameters in Table A1. In this study, the considered wave range is from 0 to 0.4 Hz, meaning that the wave circular frequency within 2.5 rad/s needs to be simulated. Therefore, the defeaturing tolerance and maximum element size are 0.5 m and 1 m, respectively. The maximum allowed wave frequency is 0.61 Hz in this scheme; 3.
Solve the model in the frequency domain and obtain Diffraction and Froude-Krylov force F e , added mass A M , and radiation damping B C .
ar. Sci. Eng. 2021, 9,656 behave a positive correlation, while radius and submergence depth show relation all the time.

Simulation Scheme
The structural schematic of the cone shaped floating body investigat is shown in Figure 4. The height of the cylinder part above the waterline is a constant, 0 Design Modeler, the 3D geometry with given parameters is created.
In this paper, ANSYS-AQWA, a commercial computation software tial flow theory, is utilized to calculate hydrodynamic parameters. The sim including numerical modeling, parameters setting, mesh generation, and cessing, can be conducted in the graphical interface directly. The basic s for each sample are as follows:  For each simulation, the given frequency range is divided into 52 frequency points. The mean power for each sample at each frequency is calculated. The results of sample 1 and sample 2 are shown in Figure 5.
With the increase in wave frequency, the capture power rises firstly and then drops steadily. For each sample, there is a unique optimal frequency in which the capture power can reach the maximum. The 100 samples' capture power are calculated so that they can be used as training set and test set for BP neutral network. Only two samples' results are presented in this figure. With the increase in wave frequency, the capture power rises firstly and then drops steadily. For each sample, there is a unique optimal frequency in which the capture power can reach the maximum. The 100 samples' capture power are calculated so that they can be used as training set and test set for BP neutral network. Only two samples' results are presented in this figure.

Theoretical Verification of Simulations
Falnes [36] illustrated that the maximum power that a heaving axisymmetric body can absorb is

Theoretical Verification of Simulations
Falnes [36] illustrated that the maximum power that a heaving axisymmetric body can absorb is where J is the wave energy flux; λ is the wavelength. For deep-water waves, λ = g/2π. J is where T and H are wave period and height, respectively. Budal's upper bound [36] gave another upper limit power that a submerged body with given volume V can absorb. It is where V is the volume of the submerged part. The point of intersection of two theoretical curves can be defined as (Tc, Pc). P c is In this study, Equations (23) and (25) are used to verify the validity of simulations. To make comparisons, the results are normalized by dividing P c . The three curves are shown in Figure 6.
It can be found that the capture power curves of two samples are in the area enclosed by curve P max , curve P u , and coordinate axes, which means the simulation scheme is accurate and reliable. All the samples are verified successfully and only two of them are demonstrated in this section. It can be found that the capture power curves of two samples are in the area enclosed by curve Pmax, curve Pu, and coordinate axes, which means the simulation scheme is accurate and reliable. All the samples are verified successfully and only two of them are demonstrated in this section.

BP Neural Network
The back propagation (BP) neural network is a kind of feedforward neural network trained by error back propagation algorithm. It is a most widely used form, and is composed of many nonlinear transformation units. This algorithm has a strong non-linear mapping ability and can simulate any nonlinear continuous functions with much higher accuracy theoretically. After the network is trained, the reflection between the inputs and outputs can be obtained and memorized. They are shown on the weights of each layer. BP neural network's structure is flexible, which means the number of layers and neurons can be changed according to research objectives. A BP neural network generally includes an input layer, one or two hidden layers, and an output layer. Full connections are applied between layers. More details about BP neural network can be seen in [37].

Neural Network Design
The first step to design a good neural network is to identify the number of hidden layers. A three layers neural network, which contains only one hidden layer, can simulate any reflection from n-dimensional inputs to m-dimensional outputs. Hence, a three-layer neural network with one hidden layer is selected in this study. Next, the nodes of each layer need to be identified. In this study, six key variables are selected, so the number of nodes in input layer is six. Only one parameter, capture power, needs to be predicted, so the number of nodes in output layer is one. Finally, the number of nodes in hidden layer needs to be identified. There is an empirical formula [38] that can be referred to identify the number of hidden nodes. l n m a = + + (27) where l, n, and m are the number of nodes in hidden layer, input layer, and output layer, respectively; a is an adjustment constant ranging from 1 to 10. In this paper, the number of hidden nodes is tested from 3 to 12 to identify the most suitable value. MSE is used to elevate the performance, and the results are shown in Table  1.

BP Neural Network
The back propagation (BP) neural network is a kind of feedforward neural network trained by error back propagation algorithm. It is a most widely used form, and is composed of many nonlinear transformation units. This algorithm has a strong non-linear mapping ability and can simulate any nonlinear continuous functions with much higher accuracy theoretically. After the network is trained, the reflection between the inputs and outputs can be obtained and memorized. They are shown on the weights of each layer. BP neural network's structure is flexible, which means the number of layers and neurons can be changed according to research objectives. A BP neural network generally includes an input layer, one or two hidden layers, and an output layer. Full connections are applied between layers. More details about BP neural network can be seen in [37].

Neural Network Design
The first step to design a good neural network is to identify the number of hidden layers. A three layers neural network, which contains only one hidden layer, can simulate any reflection from n-dimensional inputs to m-dimensional outputs. Hence, a three-layer neural network with one hidden layer is selected in this study. Next, the nodes of each layer need to be identified. In this study, six key variables are selected, so the number of nodes in input layer is six. Only one parameter, capture power, needs to be predicted, so the number of nodes in output layer is one. Finally, the number of nodes in hidden layer needs to be identified. There is an empirical formula [38] that can be referred to identify the number of hidden nodes.
where l, n, and m are the number of nodes in hidden layer, input layer, and output layer, respectively; a is an adjustment constant ranging from 1 to 10.
In this paper, the number of hidden nodes is tested from 3 to 12 to identify the most suitable value. MSE is used to elevate the performance, and the results are shown in Table 1. MSE reaches a minimum when the number of hidden nodes is 4, which is the optimal value for this case. The final BP neural network structure designed in this paper is shown in Figure 7. MSE reaches a minimum when the number of hidden nodes is 4, which is the opti value for this case. The final BP neural network structure designed in this paper is sho in Figure 7. According to the structure, the output bj of input layer can be expressed as follow where wj is the weight from the hidden layer to the output layer; θ' is the threshold v of the output layer.

Data Standardization and Neural Network Training
Before training, data standardization for individual features needs to be condu to improve training speed. The standardization formula used in this paper is where x is the standardized result; xmax and xmin are the maximum and minimum value the dataset, respectively. The standardized data have a distribution range between 0 1. The network training process is to adjust the weights and thresholds so that the v of loss function reduces to a minimum. The training parameters for this model are sho in Table 2.  According to the structure, the output b j of input layer can be expressed as follows: where w ij is the weight from the input layer to the hidden layer; x i is the input variable; θ j is the threshold value of the hidden layer. The output y of the BP neural network is where w j is the weight from the hidden layer to the output layer; θ' is the threshold value of the output layer.

Data Standardization and Neural Network Training
Before training, data standardization for individual features needs to be conducted to improve training speed. The standardization formula used in this paper is where x is the standardized result; x max and x min are the maximum and minimum values in the dataset, respectively. The standardized data have a distribution range between 0 and 1. The network training process is to adjust the weights and thresholds so that the value of loss function reduces to a minimum. The training parameters for this model are shown in Table 2. Tangent sigmoid function (tansig) is adopted for the hidden layer, and the linear function (purelin) is adopted for the output layer. In the training process, mean squared error is used as loss function. It is defined as where m is the number of samples;ŷ is the observed value; y is the real value. In this paper, the top 80 samples are defined as training set. This model is trained in MATLAB R2019a, and the trendline of MSE for training set is shown as Figure 8.
Tangent sigmoid function (tansig) is adopted for the hidd tion (purelin) is adopted for the output layer. In the training is used as loss function. It is defined as The training process is terminated at 234 epochs becau minimum (10 −7 ). The rest 20 samples are used to test, and the f de-standardized are shown in the next section.

Results and Discussion
In this section, six groups' forecasting data (ω = 0.53 ra rad/s, ω = 1.42 rad/s, ω = 1.76 rad/s, and ω = 2.09 rad/s) is giv common wave frequency. The desired outputs and forecas Figure 9 under different frequency. For each sample, the ou points can be predicted. 3000 5000 desired outputs de The training process is terminated at 234 epochs because the gradient reaches the minimum (10 −7 ). The rest 20 samples are used to test, and the forecasting results after being de-standardized are shown in the next section.

Results and Discussion
In this section, six groups' forecasting data (ω = 0.53 rad/s, ω = 0.81 rad/s, ω = 1.14 rad/s, ω = 1.42 rad/s, ω = 1.76 rad/s, and ω = 2.09 rad/s) is given because they are the most common wave frequency. The desired outputs and forecasting results are presented in Figure 9 under different frequency. For each sample, the output power at 52 frequency points can be predicted.
In Figure 9a, the deviation of five samples (85, 90, 94, 99, and 100), which are at the lowest position of the graph, are relatively large, with mean error about 60 W. In Figure 9c, the error of sample 81 is the largest, with approximately 500 W. The forecasting results of sample 95 and 96 are rather larger than desired outputs in Figure 9d,f, and the error of sample 89 is around 200 W in Figure 9f. The highest accuracy is at ω = 0.81 rad/s and almost all the forecasting points fit the desired points. In contrast, the worst result is at ω = 2.09 rad/s and there are five forecasting results deviating the desired outputs.

Results and Discussion
In this section, six groups' forecasting data (ω = 0.53 rad/s, ω = 0.81 rad/s, ω = 1.14 rad/s, ω = 1.42 rad/s, ω = 1.76 rad/s, and ω = 2.09 rad/s) is given because they are the most common wave frequency. The desired outputs and forecasting results are presented in Figure 9 under different frequency. For each sample, the output power at 52 frequency points can be predicted.  In Figure 9a, the deviation of five samples (85, 90, 94, 99, and 100), which are at the lowest position of the graph, are relatively large, with mean error about 60 W. In Figure  9c, the error of sample 81 is the largest, with approximately 500 W. The forecasting results of sample 95 and 96 are rather larger than desired outputs in Figure 9d,f, and the error of sample 89 is around 200 W in Figure 9f. The highest accuracy is at ω = 0.81 rad/s and almost all the forecasting points fit the desired points. In contrast, the worst result is at ω = 2.09 rad/s and there are five forecasting results deviating the desired outputs. To further verify the accuracy of the BP model, correlation coefficient (CC), root mean square error (RMSE), and error percentage are introduced in this section. They are defined as follows [29] where m is the number of forecasting results; t i is the desired value; y i is the output of the network; t and y are average values of desired and forecasting results, respectively. The significance analysis of ANOVA is also conducted in MATLAB R2019a, and the statistical parameters (after de-standardization) are listed in Table 3. The values of CC are greater than 0.9, meaning that the correlations with each group are well fitted. The values of RMSE do not exceed 140 W, and the error percentage is no more than 4%, indicating that desired outputs and forecasting results are reasonably fitted. All P values are close to 1, which means there is no significant difference between desired and forecasting outputs. These validation factors indicate that this model has a good prediction accuracy and meets the engineering requirement.

Conclusions
In this paper, capture power predictions of a specific shape floating body are attempted based on mathematical model, ANSYS-AQWA simulations, and BP neural network. The key variables are identified and the simulation scheme is proposed. A sample database is built by LHS and the corresponding power of each sample is calculated. In the end, a BP neural network, of which training set is from simulation results, is designed to predict the capture power at different wave frequency. Its performance and accuracy are also evaluated through statistical parameters.
According to the results, the conclusions can be given as follows: 1.
A mathematical model is constructed to identify the most important factors that affect the capture power. Four geometric parameters (radius, semi-vertical angle, mass, and submergence depth) and two system parameters (PTO damping coefficient and stiffness coefficient) are identified as key variables; 2.
A BP neural network with high accuracy is designed and it is used to predict the capture power. The error percentage of top five groups is less than 2.5%, and that of the last group is 4%. The values of CC are greater than 0.9 and that of RMSE are less than 80 W except for the third group, of which the value of RMSE is 138.7 W. The P values are close to 1. However, due to the error of simulations caused by commercial software, this method needs experimental data to support.