Effect of the Coefficient on the Performance of A Two-Layer Boussinesq-Type Model

The coefficients embodied in a Boussinesq-type model are very important since they are determined to optimize the linear and nonlinear properties. In most conventional Boussinesq-type models, these coefficients are assigned the specific values. As for the multi-layer Boussinesq-type models with the inclusion of the vertical velocity, however, the effect of the different values of these coefficients on linear and nonlinear performances has never been investigated yet. The present study focuses on a two-layer Boussinesq-type model with the highest spatial derivatives being 2 and theoretically and numerically examines the effect of the coefficient a on model performance. Theoretical analysis show that different values for α (0.13⩽α⩽0.25) do not have great effects on the high accuracy of the linear shoaling, linear phase celerity and even third-order nonlinearity for water depth range of 0<kh⩽10 (k is wave number and h is water depth). The corresponding errors using different α values are restricted within 0.1%, 0.1% and 1% for the linear shoaling amplitude, dispersion and nonlinear harmonics, respectively. Numerical tests including regular wave shoaling over mildly varying slope from deep to shallow water, regular wave propagation over submerged bar, bichromatic wave group and focusing wave propagation over deep water are conducted. The comparison between numerical results using different values of α, experimental data and analytical solutions confirm the theoretical analysis. The flexibility and consistency of the two-layer Boussinesq-type model is therefore demonstrated theoretically and numerically.


Introduction
The Boussinesq-type model (BT hereinafter for simplicity) is a powerful model for wave simulations in coastal and offshore areas, and its ability has been proven under various conditions. BT models have blossomed in the past 60 years since Peregrine (1967) pioneered the classical BT model for water waves over varying topographies. During this period, numerous BT models have been derived to improve both linear and nonlinear properties. Some typical reviews on BT models can be found in the literature (Madsen and Schäffer, 1999;Kirby, 2003Kirby, , 2016Madsen and Fuhr-man, 2010;Brocchini, 2013).
The core of the BT equations is shallow water equations, with some higher-order linear and nonlinear derivatives included to highlight the importance of the wave kinematic and dynamic properties in deeper water. Naturally, some coefficients exist in models. The values of these coefficients are very important to describing wave properties. In general, the Padé approximation of linear Stokes dispersion first proposed by Witting (1984) is widely used to determine the linear dispersion coefficients. As a result, these coefficients are conventionally pre-described with specific val-μ 4 μ z α μ 4 ues in most BT models. Actually the values of these coefficients depend on the approach/criterion utilized to optimize the embodied properties in BT models, therefore, are not unique. Usually, there are two ways to determine dispersion coefficients, i.e., the exact Padé type and the approximate Padé type. Madsen et al. (1991) added third-order derivatives to the classical BT model (Madsen et al., 1991;Madsen and Sørensen, 1992), which was expressed in terms of the horizontal depth-integrated velocities and surface elevation, to improve the linear dispersion only. Shortly after that, Madsen and Sørensen (1992) developed this model for the slowly varying bathymetry case, and the linear shoaling property of the BT model was highlighted. They followed Witting's method (1984) and set the dispersion coefficient B=1/15. Thus, the dispersion of the model could exactly match the Taylor expansion of the Stokes first-order dispersion relation on the order of ( =kh, where k is the wave number and h is water depth). Madsen et al.'s(1991) choice is called the exact Padé [2, 2] type. A different extended BT model was derived by Nwogu (1993), who used the velocity defined at an arbitrary water depth measured from still water level. Nwogu used =-0.531h to obtain the best fit between the linear dispersion relation of the BT model and the exact dispersion relation for a wider range of water depths. The dispersion relation did not exactly match the Taylor expansion of the Stokes first-order dispersion relation on the order of . Nwogu's choice is an approximate Padé [2, 2] type.
The exact Padé type highlights the importance of dispersion in shallow water and finite water depth, while the approximate Padé type highlights the importance of the application range, where it does not pursue accurate linear dispersion as the exact Padé type does in shallow water and finite water depth. Most BT models for coastal waves followed the abovementioned two ways to determine the dispersion relation of Stokes first-order waves. The research following the exact Padé way can be found in the literature (Beji and Nadaoka, 1996;Schäffer and Madsen, 1995;Madsen and Schäffer, 1998;Gobbi et al., 2000;Zou, 1999Zou, , 2000Zhao et al., 2004;Liu and Sun, 2005;Zou and Fang, 2008) while some researchers (Wei et al., 1995;Agnon et al, 1999;Madsen et al., 2002Madsen et al., , 2003Lynett andLiu, 2004a, 2004b;Chazel et al., 2009;Fang, 2015, 2016;Liu et al., 2018 used the approximate Padé type. For the exact Padétype BT models, the choice of the dispersion coefficients is very simple, and the results are also apparent. For the approximate Padé-type BT models, the choice of the dispersion coefficients varies greatly (Chazel et al., 2009;Galan et al., 2012;Simarro et al., 2013;Liu et al., 2018).
However, the dispersion coefficients determined only by linear dispersion often deteriorate the accuracy of the linear shoaling property since the dispersion coefficients contribute directly to the accuracy of the linear shoaling. For example, there is only one free coefficient in the model of Wei et al. (1995), this dispersion coefficient is closely related to the velocity location, and it affects not only the accuracy of dispersion but also the accuracy of the shoaling property. This also could be found in both exact Padé-type and approximate Padé-type BT models (Galan et al., 2012;Simarro et al., 2013), and they obtained some different values for the coefficients embodied in some conventional BT models (Kennedy et al., 2001;Madsen and Schäffer, 1998;Lynett and Liu, 2004a) by minimizing a globally linear error that includes celerity and shoaling errors, and the results showed that the accuracy of the linear celerity and shoaling was very sensitive to the choice of the coefficients value. To avoid the deterioration of the linear shoaling property, some artificial amending methods have been proposed in many studies (Schäffer and Madsen, 1995;Madsen and Schäffer, 1998;Zou, 1999Zou, , 2000Madsen et al., 2002Madsen et al., , 2003Madsen et al., , 2006Liu and Sun, 2005;Zou and Fang, 2008;Chazel et al., 2009;Fang, 2015, 2016;Liu et al., 2018). These methods improve the accuracy of the linear shoaling property. However, it will deteriorate the accuracy of the local velocity profiles over a varying bottom. To improve the linear shoaling property, Bingham et al. (2009) utilized different methods and new formulas, pseudo horizontal velocity potential and pseudo vertical velocity component, where the effect of linear shoaling was highlighted.
The linear shoaling properties are not only affected by the dispersion coefficients but also the shoaling coefficients. There are many methods to optimize shoaling property. Madsen and Sørensen (1992) was the first one to optimize the shoaling property in BT model, where the linear shoaling gradient was set as the target to determine the shoaling coefficient. Chen and Liu (1995) proposed a new formula to demonstrate the actual linear shoaling performance of BT model, and in this formula, the error of linear shoaling amplitude from deep water to shallow water was highlighted. Later, this formula was applied in the literatures (Galan et al., 2012;Simarro et al., 2013, Liu andFang, 2016). To avoid using the complex formula proposed by Chen and Liu (1995), Madsen et al. (2002) still used the linear shoaling gradient as a standard to determine the shoaling coefficient, and the formula divided by wave number was used to focus on the importance of shoaling performance in shallow water. Regardless of how to determine the linear dispersion or linear shoaling, the values of the coefficients embodied in the BT models are the determined ones if the criteria are set for a given application range (Galan et al, 2012;Simarro et al., 2013;Liu and Fang, 2016).
Recently,  analyzed the linear dispersion and shoaling properties of the multilayer BT model with the order of highest spatial derivatives being 2 derived by Liu et al. (2018). They did not introduce the amending methods mentioned above to improve the shoaling property. However, with a 1% tolerance error, the accurate linear shoaling property is surprisingly accurate up to k max h=15.8, 41.7 and 143.9 for the two-layer, three-layer and four-layer models, respectively. They also mentioned that the application range of the models was mainly restricted by the accuracy of the velocity profiles from the comprehensive aspects. Is the choice of the coefficients a coincidence? If this is not a coincidence, can the values of the coefficients vary within a range without apparently affecting other properties? The simulated free surface elevations almost keep unchanged when utilizing different coefficient values in numerical simulations from the two-layer BT model . We found it interesting as it goes against our previous knowledge. Can these be explained from a theoretical point view and be verified by further numerical simulations? The questions motivated us to conduct the present research. Therefore, the two-layer BT model  is considered, and the effect of the coefficient on the linear dispersion, linear shoaling and other properties is examined theoretically and numerically.
In Section 2, a brief introduction of the two-layer BT equations is given. The effect of the coefficient on linear dispersion, linear shoaling and other properties are theoretically analyzed in Section 3. In Section 4, numerical simulations of the water waves under different conditions are conducted to examine the effect of the coefficients on the numerical results. Finally, the conclusions are given in the last section. Liu et al. (2018) derived three new sets of multilayer BT models. Shortly after the derivation,  analyzed the linear and nonlinear properties of the multilayer BT model with the highest spatial derivatives being 2. The Cartesian coordinate system of the two-layer case is shown in Fig. 1, where the x-and y-axes are located on the still water plane, and the z-axis points vertically. The main assumptions employed in these derivations include a homogeneous fluid of constant density, incompressibility and irrotationality of the flow as well as the absence of viscous and surface tension effects. The resultant twolayer model embodies high accuracy in linear, nonlinear, shoaling, and kinematic properties from deep to shallow water . Specifically, with =0.172, the model is applicable up to kh ≈ 19.7 in linear dispersion, up to kh ≈ 15 in the second nonlinear property, and up to 0<kh<15.8 in the linear shoaling property within 1% error.

The two-layer BT model
The vertical 1-D two-layer BT model consists of nine equations, and the main variables in the model include the surface elevation , the horizontal velocity and the vertical velocity defined on the surface, the horizontal velocity u 10 and the vertical velocity w 10 defined on the still water level, and two groups of computational velocities (u i * , The kinematic and dynamic conditions at the free surface are (1) . (2) The expressions between the velocities on the free surface and the velocities at still water level are u η =u 10 +η (3) The expressions between the velocities at still water level and the first group of computational velocities are (5) At the interface between two layers, the velocities satisfy the conditions At the bottom, the kinematic boundary condition is where the subscription x denotes the partial derivatives. The other variables related to the water depth are According to the previous study , the coefficients must satisfy , thus, we can set , then .

Effect of the different values of the coefficient on the linear and nonlinear properties
The two-layer BT model is analyzed theoretically to examine the effect of the coefficient on its embodied characteristics with respect to six properties: linear dispersion, second-order harmonics, third-order harmonics, amplitude dispersion, linear velocity profiles, and linear shoaling amplitude. To conduct the analysis, the methods used in the literature (Liu and Fang, 2016;Liu et al., 2018) were employed, and details can be found in these two studies. The linear dispersion, linear velocity profiles and nonlinear characteristics are theoretically computed when the water depth is set to be a constant in Eqs. (1)−(9), and the results are analyzed in Subsections 3.1 to 3.3. To compute the shoaling, a mild sloping is assumed with the nonlinear terms being kicked out in Eqs. (1)−(9), and the results are computed and analyzed in Subsection 3.4. The linear velocity-dominant method (VDM) and the dispersion-dominant method (DDM) were used in Liu et al. (2018) to obtain the dispersion coefficient . At the present stage, in order to investigate the effect of the coefficient on the linear dispersion, we do not use any method to determine the unique coefficient. The error (100(c b− c s )/c s (%)) of the linear phase celerity between the present model and analytical solutions is plotted in Fig. 2, in which the range of the coefficient is from 0 to 0.5, and the error figure is symmetric at symmetric line =0.25. As shown in this figure, with a 0.1% tolerance error in the range of 0<kh<12, the coefficient roughly varies from 0.13 to 0.37, and with a 1% tolerance error in the range of 0<kh<18, the coefficient roughly varies from 0.065 to 0.435. To further demonstrate the phase celerity curve, the results are plotted in Fig. 3 when equals 0.13, 0.172 and 0.25. The phase celerity is the most accurate within the range of 0<kh<10 when =0.25. However, it should be noted that phase celerity is also accurate (0.1% error) when =0.13 or =0.172. α Therefore, with respect to linear dispersion, we recommend choosing the coefficient from 0.13 to 0.25. The velocity errors for F u (kh) and F w (kh) proposed by Madsen et al. (2003) and later used in Liu and Fang (2016) and Liu et al. (2018 are also employed in the present study to examine the effect of , and the formulas are written as: where u s (0) and w s (0) are the maximum horizontal and vertical velocities at the still water level, respectively. = 0.13, the maximum application water depths are 9.8 and 4 at a 2% tolerance error and 1% tolerance error, respectively. When = 0.172, the maximum application water depths are 7.7 and 4.5 at a 2% tolerance error and 1% tolerance error, respectively. When = 0.25, the maximum application water depths are 5 and 3.8 at a 2% tolerance error and 1% tolerance error, respectively. Thus, when the value of is taken from 0.13 to 0.25, the maximum application water depth is 5 at a 2% tolerance error and 3.8 at a 1% tol-   α erance error. When the value of is taken from 0.13 to 0.172, the maximum application water depth is 7.7 at a 2% tolerance error and 4 at a 1% tolerance error. The error of the second-order harmonics, third-order harmonics and amplitude dispersion between the present model and analytical solutions are plotted in Fig. 8, Fig. 9 and Fig. 10, where the value of varies from 0.1 to 0.25. And the effects of on second-order harmonics, third-order harmonics and amplitude dispersion are clearly shown in these figures. Especially, the detailed corresponding results for =0.13, 0.172 and 0.25 are plotted in Fig. 11. With the increase in the coefficient from 0.13 to 0.25, the tendency is to be more accurate for the nonlinear property for kh ranges from 0 to 10. As shown in this figure, with respect to nonlinear performance, the application range of the model is restricted by the accuracy of the third-order harmonics. Therefore, the application range of the model is described with respect to third-order harmonics. At a 1% tolerance error, the maximum application range is kh=14.65, 13.57 and 12.64 for equal to 0.13, 0.172 and 0.25, respectively. Thus, with a 1% error for 0<kh<10, the accurate nonlinear property can α be found for any values of in the range of 0.13−0.25. It is also noticed that the accuracy of the nonlinear property decreases slightly when compared against the linear dispersion mentioned above. However, with respect to the nonlinear performance, this model is much more accurate than any conventional BT model. Since the coefficient can be chosen in the range of 0.13−0.25, the model is also more attractive.

Effect of on linear shoaling
The shoaling coefficients are often determined by matching the shoaling gradients between the models and analytical solutions (Madsen and Schäffer, 1998;Zou, 1999Zou, , 2000Liu and Sun, 2005). However, this may decrease the accuracy of the shoaling amplitude (Chen and Liu, 1995;Galan et al., 2012;Simarro et al., 2013;Liu and Fang, 2016). To avoid this problem, the following expression is usually employed to minimize the errors between the shoal-      40 SUN Jia-wen et al. China Ocean Eng., 2021, Vol. 35, No. 1, P. 36-47 ing amplitudes of the BT model and the analytical solution (A s ) (Chen and Liu, 1995;Simarro et al., 2013;Liu and Fang, 2016).
, and k 0 is the wavenumber in deep water. Through using the following formulation to measure error, the effect of the coefficient on linear shoaling performance is examined. (12) Fig. 12 shows the linear shoaling gradient of the model and the error of shoaling amplitude from the model and the analytical solution. The shoaling gradient varies little for 0<k 0 h<10 for these three coefficients. At a 1% tolerance error, the maximum application water depth is k 0 h =14.35 for the three coefficients, and at a 0.1% tolerance error, the maximum application water depth is k 0 h =10 for the three coefficients. Thus, for a given value ranging from 0.13 to 0.25, the accurate linear shoaling amplitude can be found for the application range 0<k 0 h<10. For a given range 0<kh≤10, when any value of is in the range of 0.13≤ ≤0.25, the error of both the dispersion and linear shoaling is smaller than 0.1% when 0.13≤ ≤ 0.25, and the error of nonlinear property (restricted by thirdorder harmonics) is smaller than 1%. As stated in , the accuracy of the linear velocity profile is the main restriction of the two-layer model, and this still holds true for 0.13≤ ≤0.25. According to the analysis in Subsection 3.2, we can see that =0.13 might be a good choice (2% tolerance error) for 0<kh≤10 and =0.172 might be a good choice (1% tolerance error) for the range of 0<kh≤ α 4.5. Therefore, it is not unique to choose the values according to the conditions. However, notice that the velocity profiles are mainly in the range of 0<kh≤3.14 near coastal areas, the velocity profiles make no greater difference for 0.13≤ ≤0.25 since the accuracy of the velocity profiles is more than 99%.
If the evolution of surface elevation is the only target for the two-layer BT model, any value of the coefficient can be chosen in the range of [0.13, 0.25]. Therefore, it is more flexible to choose the value of the coefficient. The high accuracy enables the model to be more suitable for wave Bragg reflection over a more complex bottom case (Madsen et al., 2006), in which they used fixed values to enhance the shoaling property.

Numerical verifications
α α  established multilayer BT numerical models with the highest spatial derivatives being 2. Detailed numerical implantation can be found in . The only change in this two-layer BT is that the coefficient is within 0.13≤ ≤0.25. Without loss of generality, all the following numerical simulations are conducted with =0.13, 0.172 and 0.25 in the model. 4.1 Linear wave shoaling from deep water to shallow water (0.25≤kh≤10) To investigate the linear shoaling performance of the two-layer model with different values, numerical experiments of a linear wave propagating over a slowly varying topography are carried out. The linear waves with T=6.4 s propagate from a deep water depth of h 1 =101.782 m (kh=10) to a shallow water depth of h 2 =0.623 m (kh=0.25), and the connection between the two water depths is a smooth bottom profile (Bingham and Zhang, 2007;Benoit et al., 2017;, which can be written as: where the length l=3 180 m, and the starting location at x b =270 m is used in a numerical model. The total distance  42 SUN Jia-wen et al. China Ocean Eng., 2021, Vol. 35, No. 1, P. 36-47 ution over a submerged bar was conducted by Luth (1994). The experimental data are widely used for the validations of the various numerical models. The case of regular wave period T=2.02 s (kh=0.673) with wave height H=0.02 is used in the present paper. According to previous studies and our experience in using BT numerical models, accurately simulating the wave evolution requires a model to possess the maximum range kh=6 at 1% tolerance error. As shown, the linear and nonlinear performance of the two-layer BT model with different values ranging from 0.13 to 0.25 is suitable for this case. between the numerical results and the experimental data is found, which indicates that the effect of the coefficient on the numerical results is small. The velocities on the surface elevations in such a case have never been shown numerically and experimentally in any literature. Therefore, we plot the time history of the computed horizontal velocity and vertical velocity on the surface elevations, as shown in Fig. 15 and Fig. 16. It can be seen that there is no apparent difference among the computed time history of the horizontal velocity and vertical velocity on the surface elevation with different values of =0.13, 0.172 and 0.25, which also demonstrates that the effect of the coefficient value on numerical results is small. Overall, the value of coefficient having no great effect in the two-layer BT model is numerically validated. Stansberg (1993) carried out an experimental study of  (1996) carried out an experimental study of a focusing wave group evolution over a constant water depth. In the present simulations, the linear incident wave at the left boundary is used like that in . A small error appeared in the numerical codes of avoiding numerical filtering at the first simulation point. Thus, this small error is corrected before conducting numerical simulations. The simulations on the D55 case (kh=2. 026-4.403, Baldock et al., 1996)   α file at the largest wave crest location from the two layers with different values of =0.13, 0.172 and 0.25 and the experimental data from Baldock et al. (1996) are made, as shown in Fig. 18. Overall, both the surface elevation and the horizontal velocity profile at the focused point calculated from the numerical model with different values of the coefficient are in good agreement with the experimental data, α α demonstrating that the two-layer model indeed possesses a high accuracy in linear and nonlinear properties even for the different values taken from 0.13≤ ≤0.25. A larger surface elevation at the focused position is found being compared with that in the previous study in . This difference is mainly caused by incorrect numerical filtering, as noted in the abovementioned content. Overall, the effect of the different values of 0.13≤ ≤0.25 is small.

Conclusions α α
Based on the two-layer BT model (Liu et al, 2018), the effect of the different values of the coefficient embodied in the model on the linear and nonlinear properties is investigated theoretically and numerically. An interesting result is that the different values for (0.13≤ ≤0.25) do not have great effects on the high accuracy of the linear shoaling, linear phase celerity and even third-order nonlinearity for the range of 0<kh≤10, which was demonstrated theoretically and partially demonstrated numerically. The main detailed conclusions are drawn as follows. α (1) From a theoretical perspective, when 0.13≤ ≤0.25, for the application water depth range of 0<kh≤10, the error of both the linear dispersion and linear shoaling amplitude can be controlled within 0.1%, and the error of second-order harmonics, third-order harmonics and amplitude dispersion can be controlled within 1%.  (2) From a theoretical perspective, when 0.13≤ ≤ 0.172, the error of the linear velocity profiles can be controlled within 2% for the application water depth range of 0<kh≤7.8. When 0.13≤ ≤0.25, the error of the linear velocity profiles can be controlled within 2% for the application water depth range of 0<kh≤4.5. α α α (3) From a numerical aspect, for 0.25<kh≤10, the high accuracy of the linear shoaling in the two-layer model is demonstrated for 0.13≤ ≤0.25. In the other numerical experiments provided in the present study, the effect of different values of (0.13≤ ≤0.25) on the numerical results is small.