Spectral analysis for nonlinear wave forces

Applied Ocean Research 20 (1998) 199–212

Spectral analysis for nonlinear wave forces O.M. Boaghe a, S.A. Billings a, P.K. Stansby b a

Department of Automatic Control and Systems Engineering, University of Sheffield, P.O. Box-600, Mappin Street, Sheffield S1 3JD, UK b Hydrodynamics Research Group, School of Engineering, University of Manchester, Oxford Road, Manchester M13 9PL, UK Received 30 June 1997; accepted 30 March 1998

Abstract The spectral analysis of nonlinear wave forces is investigated based on the recently introduced Dynamic Morison equation. This equation accounts for nonlinearities associated with history or vortex shedding effects accurately, in contrast to the standard Morison equation which has undesirable features in this respect. An expression for the response power spectrum is initially derived in terms of the higher-order nonlinear frequency response functions and the spectrum of the input velocity. It is then shown how this can be evaluated for the Dynamic Morison model to yield an expression for the output power spectrum in terms of the coefficients of the time-domain differential equations and properties of the input. A comparison of the output spectra computed using traditional methods and the new expression of the response spectrum over several data sets is included. An interpretation of the Morison equation response spectrum is finally given. q 1998 Elsevier Science Ltd. All rights reserved. Keywords: System identification; Fatigue; Nonlinear systems; Wave forces

1. Introduction Spectral analysis is widely used in the study of wave forces and safety assessment of offshore structures. The force spectrum computation is in general based on Borgman’s linearised form of the drag term in the Morison equation (proportional to ulul where u is velocity). An important application of spectral analysis is in the determination of fatigue life and methods of prediction are well described in Ref. [1] where the importance of peak statistics and zero-crossing rates is emphasised. However, in this paper we are concerned at a more fundamental level with the validity of the Morison equation and the resulting force spectra for fixed cylinders only. Worden et al. [2] point out that for the Morison equation the linear transfer function H 1(q 1) → ` as q 1 → ` which is of course physically unrealistic, suggesting that small wave energy at high frequency produces high force energy. Force spectra have also been computed with the addition of a u 3 term to approximate better ulul. However, as well as nonlinear effects simplistically represented as drag, there are so-called history effects associated with vortex shedding [3]. Attempts have previously been made to extend the Morison equation to represent these effects, generating terms resembling those in the Duffing oscillator, producing the so-called Morison–Duffing equation [2,3]. While this equation generally improved the fit to measured data, it surprisingly proved of little value for long term

prediction. Our recent analysis suggests that the problem lies in the way the model was reconstructed and the way in which noise is handled; and noise is present in all data to a greater or lesser extent. To progress it is essential that noise on the data is adequately accommodated when the models are reconstructed. This is achieved with the NARMAX model which does not have a continuous-time form (e.g. with differentials), but rather has a discrete form with input and output terms (of any algebraic combination) from the present and previous time steps (see following section). While such a discrete form cannot be unique it does enable nonlinear frequency response functions to be generated. Swain et al. [4] used these noise-free high-order spectra to generate continuoustime force equations following the technique of Swain and Billings [18], calling the new generally valid equation the Dynamic Morison equation. This equation reproduced the NARMAX time histories extremely accurately but ironically did not fit the original data any better than the Morison equation; but this is due to the noise present in the data. The Dynamic Morison equation in fact has a simpler form than the previous Morison–Duffing equation. Reproducing force spectra via NARMAX modelling is a rather complex and highly specialist operation. In this paper a method for deriving the force spectrum directly from the terms in the Dynamic Morison equation will be presented and agreement with that obtained from NARMAX will be shown to be excellent. The scheme is applied to the five data

0141-1187/98/$ – see front matter q 1998 Elsevier Science Ltd. All rights reserved PII S 01 41 - 11 8 7( 9 8) 0 00 2 0- 0

200

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

sets used in Ref. [4], referred to as Salford [5], De Voorst and Christchurch Bay [3], where nonlinear effects range from being small (Salford and De Voorst) to substantial (Christchurch Bay). All cases relate to a fixed vertical cylinder in random waves where force and velocity data are obtained for a small spanwise element. The Salford data was obtained in a laboratory wave flume for unidirectional waves with rectangular velocity spectra pof different width; the Keulegan–Carpenter D is number (KC ¼ 2urms T=D, p  diameter) was about 5 and the Reynolds number (Re ¼ 2urms D=n, n is kinematic viscosity) about 3000, using the period T at the mid point of the spectra. The De Voorst data was obtained in a much larger outdoor flume with a Pierson–Moskowitz wave spectrum giving a Re number of 2.4 3 10 5 and a KC number of about 6; the defining period is at the peak in the spectrum. The Christchurch Bay data was obtained from an instrumented tower at a coastal site; the waves are directional, the Re number was about 5 3 10 5 and the KC number about 20. Only forces in the predominant wave direction are considered. The main advantage of the new approach is the direct link which is established between the time-domain model coefficients, the velocity input and the force spectrum. This has implications for fatigue-life analysis due to the sensitivity of the spectrum to the nonlinear model parameters. We should also mention that forces due to vortex shedding transverse to the flow direction, known as lift forces, can also be significant with frequency components higher than the wave frequency but this is not considered here. In the next section, the modelling of nonlinear forces is described in some more detail. The derivation of nonlinear spectra from general continuous-time models is then given in Section 3. Spectra are finally produced for the various data sets mentioned above.

2. Modelling nonlinear wave forces Four models are presented and discussed in more detail in this section: the Morison equation, the Morison–Duffing equation, NARMAX parametric models and Dynamic Morison models. The traditional way of modelling wave forces when drag is significant is based on the Morison equation [6] which assumes that the force is composed of a linear combination of inertia and drag forces. The Morison equation can be written as: 1 du(t) 1 du(t) þ rDCd u(t)lu(t)l ¼ Ki F(t) ¼ prD2 Cm 4 dt 2 dt þ Kd u(t)lu(t)l

ð1Þ

where F(t) is the force per unit axial length, u(t) is the instantaneous flow velocity, r is the water density and D is the cylinder diameter. The drag and inertia coefficients C d and C m depend on the characteristics of the flow. The Morison equation generally predicts the main trends in the measured data quite well, however some

characteristics of the flow are not well represented. These are revealed in a frequency domain analysis. It is not possible to directly map the Morison equation into the frequency domain using the techniques of Volterra series due to the presence of the drag term u(t)lu(t)l which has a discontinuous derivative. A necessary condition for the existence of the Volterra series is that all nonlinearities must be infinitely differentiable [7]. The drag term can however be approximated by a polynomial of the form: u(t)lu(t)l ¼ a1 u(t) þ a3 u3 (t) þ a5 u5 (t) þ …

(2)

This approximation holds under the assumption that the input signal is a zero-mean Gaussian signal whose odd order moments are identically zero [8]. The approximated Morison equation, for nonlinearities in the input up to third order, is therefore given by: du(t) þ Kd1 u(t) þ Kd3 u3 (t) (3) F(t) ¼ Ki dt The linear transfer function of the approximated Morison equation was computed by Worden et al. [2]: H1 (q1 ) ¼ Kd1 þ jq1 Ki . As mentioned in Section 1, limq1 →` lH1 (q1 )l ¼ `, when the cylinder is subject to a high frequency input wave of very small amplitude, it will experience an extremely high force. However, for most physical systems the gain of the system falls off as the frequency increases. The influence of vortex shedding was modelled in the Morison–Duffing equation by higher order and derivative terms in F(t). Additional terms proportional to F 2(t), F 3(t), (dF(t))/(dt) and (d 2F(t))/(dt 2) were included, terms also present in the Duffing nonlinear oscillator equation. After some preliminary tests on U-tube data it became apparent that the F 2(t) term could be discarded as insignificant. The model also showed a degree of improvement when the F 3(t) term was replaced by one of the form F(t)lF(t)l. The form of the Morison–Duffing equation thus became: d2 F(t) dF(t) 1 du(t) þ a3 F(t)lF(t)l þ F(t) ¼ prD2 Cm þ a2 dt 4 dt dt2 1 ð4Þ þ rDCd u(t)lu(t)l 2 Nonlinear wave force models can also be constructed using nonlinear system identification techniques based on the NARMAX approach. The NARMAX (Nonlinear Auto Regressive Moving Average with eXogenous inputs) model is an expansion of past inputs, outputs and noise terms [9].

a1

F(k) ¼ f [F(k ¹ 1), …, F(k ¹ nF ), u(k ¹ 1), …, u(k ¹ nu ), 3 e(k ¹ 1), …, e(k ¹ ne )] þ e(k)

ð5Þ

F(k), u(k) and e(k) are the discrete-time output, input and noise respectively at the time intervals k and f is a nonlinear rational or polynomial function. System identification based on this model consists of determining the model structure, followed by parameter estimation and model validation. Nonlinear wave force models were identified using the NARMAX model by Worden et al. [2].

201

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

The predictive performance of the estimated models were found to be better than the Morison equation but the models were difficult to interpret. Although identified NARMAX models are not necessarily unique for a particular data set, if they capture the underlying system characteristics all the models should exhibit the important dynamic features of the system. For example, if the NARMAX model represents the underlying system accurately it must reflect the unique linear and nonlinear frequency behaviour of the system. Swain et al. [4] developed and discussed a new model structure, the Dynamic Morison equation, based on the frequency domain uniqueness of identified NARMAX models associated with nonlinear wave forces. The procedure employed to determine the model consisted of two stages: estimation of a discrete-time (NARMAX) model from sampled input–output data and computation of the generalised frequency response functions. Continuous-time nonlinear differential equation models were then identified by curve fitting to the complex linear and higher order frequency response data using a weighted complex orthogonal estimator [4]. The continuous-time model is very similar to the Morison and Duffing–Morison equations and is referred to as the Dynamic Morison equation: d2 F(t) dF(t) 1 du(t) þ F(t) ¼ prD2 Cm þ a2 dt 4 dt dt2 1 ð6Þ þ rDCd u(t)lu(t)l 2 Depending on the different experimental conditions, the model of Eq. (6) may contain higher-order nonlinear terms of the input. For the present analysis a third-order NARMAX model was found to be sufficient to model all the relevant features of the data and hence the estimated continuous-time model incorporates nonlinear terms up to third degree: a1

d2 F(t) dF(t) du(t) þ F(t) ¼ Ki þ Kd1 u(t) þ Kd3 u(t)3 þ a2 2 dt dt dt (7) The constants and coefficients associated with the Morison (Eq. (3)) and Dynamic Morison (Eq. (7)) models were identified for the Salford Cylinder, De-Voorst and Christchurch Bay data sets by Swain et al. [4] and are summarised in Table 1. a1

The Dynamic Morison equation (Eq. (6)) has two additional terms with coefficients which may be nondimensionalised with reference to a time scale. For this purpose the wave period associated with the peak in the spectrum (or the mid point in the case of the Salford data) were chosen. Eq. (6) can then be given by gh2 Tw2

d2 F(t) dF(t) du(t) þ F(t) ¼ bm þ bd u(t)lu(t)l þ gh1 Tw 2 dt dt dt (8)

It is always desirable that coefficients of the equation should vary as little as possible from one situation to another and after experimenting algebraically, Eq. (8) can be written in the form gh0 (gh1 Tw )2

d2 F(t) dF(t) du(t) þ F(t) ¼ bm þ gh1 Tw 2 dt dt dt þ bd u(t)lu(t)l

ð9Þ

The results for the five test situations (together with T w, C d, C m, KC, Re values) are given in Table 2. The values of gh0 are very similar apart from Salford Data Set 3 which is unusual in that it has a wide band spectrum. The narrower bandwidths for Salford Data Sets 1 and 2 and for the Christchurch Bay and De-Voorst data have gh0 ¼ 0:86 6 0:04. This is a remarkably small variation considering the range of KC, Re and wave conditions (from unidirectional wave flumes to directional field conditions). It is also interesting to compare the coefficients for the Salford data with those for the De-Voorst flume. The KC values are similar (about 5) but the Reynolds numbers differ in magnitude by approximately a factor of 80. It is generally accepted that vortex shedding becomes weaker as Reynolds number increases (for a comprehensive review of such matters see Refs. [10,11]). Vortex shedding generates the history effects and consistently the value of gh1 . 0:1 for the Salford data (Re . 3 3 10 3}) and gh1 . 0:03 for the DeVoorst data (Re . 2 3 10 5}). On the other hand the Reynolds number of the De-Voorst and Christchurch Bay data are similar (within a factor of 2) but the KC values are very different. At KC ¼ 6 for the former vortex shedding is generally less strong than for KC . 20 for the latter. The corresponding values of gh1 are 0.03 and 0.06. The values of gh1 are therefore consistent with our knowledge of vortex shedding and its dependence on KC and Re.

Table 1 Coefficients for Morison (Eq. (3)) and Dynamic Morison (Eq. (7)) continuous-time models [4] Data set

Equation type

Ki

Kd1

Kd3

a1

a2

Salford Set 1

Morison Dynamic Morison Dynamic Morison Dynamic Morison Dynamic Morison Dynamic

2.32 2.14 2.12 2.08 2.14 1.90 148.67 171.14 169.66 171.94

1.56 2.09 1.39 1.99 1.72 1.49 54.33 23.00 27.05 4.81

171.68 108.12 208.40 116.44 152.47 162.32 32.02 23.25 72.22 74.00

– 0.04 – 0.04 – 0.04 – 0.27 – 0.02

– 0.22 – 0.22 – 0.19 – 0.56 – 0.16

Salford Set 2 Salford Set 3 Christchurch Bay Set De-Voorst

Morison Morison Morison Morison Morison

202

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

It is also interesting to note the values of C d and C m. C m varies little but C d . 0.4 at high Reynolds numbers is much smaller than at the lower Reynolds numbers where C d . 1.8. While this overall trend with Reynolds number is expected, C d . 0.4 is less than what is normally obtained at Re . 10 5 from the Morison equation where C d . 0.6. The results in terms of nondimensional coefficients are thus quite encouraging but many more situations need to be analysed to confirm these preliminary findings. Spectral analysis derived for the Morison equation, the NARMAX model and the Dynamic Morison equation, for all five data sets analysed will be studied in the next section. New characteristics of the Morison equation can be revealed and interpreted by using the estimated response spectrum. 3. Spectral analysis of nonlinear wave forces There are several ways in which spectral analysis can be performed. Spectral estimates can be obtained from measured input and output signals, from a predicted output signal via an identified continuous-time differential equation, or from the predicted output spectrum via the nonlinear frequency response functions. The latter case is particularly informative because an analytic expression for the power spectral density can be obtained, as a function of the nonlinear frequency response functions and the parameters of the system model. The application of spectral analysis to nonlinear wave forces is directly dependent on the type of input signal applied to the system, that is deterministic or a random signal. For deterministic input signals, the response power spectrum has been derived, for example by Worden et al. [2]. In real situations the input signal is random and a new approach has to be adopted. The response power spectrum for nonlinear wave forces is derived in the following sections, when random input signals are applied to the system.

The nth-order kernel of Eq. (11),hn ðt1 ; t 2, …, t n) is often referred to as a nonlinear impulse response of order n. The Fourier transform of h n(·) is called the nonlinear transfer function or generalised frequency response function of order n: ` ` Z Z … hn (t1 , t2 , …, tn ) exp( ¹ j Hn (q1 , q2 , …, qn ) ¼ ¹`

¹`

3 (q1 t1 þ q2 t2 þ … þ qn tn )) 3 dt1 dt2 … dtn

ð12Þ

The power spectral density for a nonlinear system was derived by Rugh [13], with a real, stationary, zero-mean and Gaussian input and is applied below to the nonlinear wave forces case. For a stationary random input signal the output power spectral density can be expressed [13] as: 1 SFF (q) ¼ (2p)2n ¹ 1

` Z

…

¹`

` Z

(2n) SFF (g1 , …, g2n )d

¹`

3 (q ¹ g1 ¹ … ¹ gn ) dg1 … dg2n ¼ ¼ 3

` Z

…

¹`

` Z

1 (2p)2n ¹ 1

Hn (g1 , …, gn )Hn (gn þ 1 , …, g2n )

¹`

3 S(2n) uu (g1 , …, g2n )d(q ¹ g1 ¹ …¹ gn )dg1 … dg2n ð13Þ To achieve further simplification it is assumed that the input is a real, stationary random process with zero-mean and a Gaussian distribution. In this case the partial output spectrum is SFn Fm ¼ 0, when n þ m is odd, and when n þ m is even [13], is given by: SFn Fm (q) ¼

3.1. Response power spectrum estimation for nonlinear systems — an overview

1 (2p)(n þ m ¹ 2)=2

` X Z p ¹`

…

` Z

Hn (g1 , …, gn )Hm

¹`

3 (gn þ 1 , …, gn þ m )d(q ¹ g1 ¹ … ¹ gn )

nY þm

Suu

j, k

Since all the models discussed above are nonlinear, methods of analysing nonlinear models in the frequency domain will be briefly discussed below. The traditional description of nonlinear systems has been based on the Volterra series. Volterra [12] showed that every continuous functional F(u) in the field of continuous functionals can be represented by the Volterra expansion: ` X Fn [u] (10) F(u) ¼ n¼0

in which the index n is the degree of the functional and F n[·] is a regular homogeneous functional of the form: ` ` Z Z … hn (t1 , t2 , …, tn )u(t ¹ t1 ) Fn [t] ¼ ¹`

¹`

u(t ¹ t2 )…u(t ¹ tn ) dt1 dt2 … dtn

ð11Þ

X

3 (gj )d(gj þ gk ) dg1 … dgn þ m

ð14Þ

where p is a sum of products over a set of n/2 (unordered) pairs of integers from 1, 2, …, n. 3.2. Response power spectrum estimation for wave forces Both the Morison (Eq. (1)) and Dynamic Morison (Eq. (6)) equations contain the drag term u(t)lu(t)l, which can be approximated in certain conditions by an infinite polynomial series [8]. Swain et al. [4] found in the analysis of five data sets that a third-order nonlinearity was sufficient to model all the relevant features. In this section, the general formula of Eq. (14) is applied to the Morison (Eq. (3)) and Dynamic Morison (Eq. (7)) equations, in the case where nonlinearities up to the third-order are considered. In the general case, where nth-order nonlinearities need to be

203

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

considered, the analysis would follow the same steps and additional higher order terms would be included. The response spectrum formula of Eq. (14) for wave forces applied to the Morison (Eq. (3)) and Dynamic Morison (Eq. (7)) equations is based on the nonlinear frequency response functions up to third order. By applying Eq. (14) for a thirdorder nonlinear system, the output power spectrum is given as: ` Z 3 H3 SFF (q) ¼ H1 (q)H1 ( ¹ q)Suu (q) þ H1 (q)Suu (q) 2p ¹`

3 ( ¹ q, g, ¹ g)Suu (g) dg þ 3

` Z ¹`

þ

1 p

` Z

3 H ( ¹ q)Suu (q) 2p 1

1 H3 (q, g, ¹ g)Suu (g) dg þ d(q) 2p

H2 (q1 , q2 ) ¼ 0 H3 (q1 , q2 , q3 )

` ` Z Z

H2 ¹` ¹`

3 (g1 ,¹ g1 )H2 (g2 , ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2 þ H2 (q ¹ g, g)H2 ( ¹ q þ g, ¹ g)Suu (g)Suu (q ¹ g)dg

¹` `

6 þ (2p)2

Z Z

H3 (q ¹ g1 ¹ g2 , g1 , g2 )H3 ( ¹ q þ g1 þ g2 ,

¹` ¹`

¹ g1 , ¹ g2 )Suu (q ¹ g1 ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2 þ 3 Suu (q)

9 (2p)2

H3 (q, g1 , ¹ g1 )H3 ( ¹ q, g2 , ¹ g2 )Suu (g1 )Suu

¹` ¹`

3 (g2 ) dg1 dg2 ð15Þ This expression can be simplified if it is assumed that H 2(q 1,q 2) ¼ 0 for nonlinear wave forces and if the symmetry properties available for the input spectrum S u u and the transfer functions are considered. Since the properties of spectral conjugation hold for the nonlinear transfer function Hnp (q1 , …, qn ) ¼ Hn ( ¹ q1 , …, ¹ qn )

(16)

Eq. (15) then reduces to: 3 SFF (q) ¼ lH1 (q)l Suu (q) þ H1 (q)Suu (q) 2p 2

` Z

6 (2p)2

H3 ( ¹ q, g,

¹`

3 ¹ g)Suu (g) dg þ H1 ( ¹ q)Suu (q) 2p

þ

Kd3 ð18Þ 1 ¹ a2 (jq1 þ jq2 þ jq3 ) ¹ a1 (jq1 þ jq2 þ jq3 )2 The transfer function expressions (Eq. (18)) can now be substituted into Eq. (17) and the response spectrum determined as: ¼

SFF (q) ¼ Suu (q)

`

` ` Z Z

density, the first- and third-order transfer functions H 1(q) and H 3(q 1,q 2,q 3) should be determined and substituted in Eq. (17). Consider the Dynamic Morison model (Eq. (7) ) developed by Swain et al. [4]. By applying the probing method developed by Peyton-Jones and Billings [14] to the continuous-time differential equation (Eq. (7)), the nonlinear frequency response functions can be expressed explicitly in terms of the continuous-time model coefficients as: Kd1 þ Ki (jq) H1 (q) ¼ 1 ¹ a2 (jq) ¹ a1 (jq)2

` Z

H3 ¹`

3 (q, g, ¹ g)Suu (g) dg þ ` ` Z Z lH3 (q ¹ g1 ¹ g2 , g1 , g2 )l2 Suu (q ¹ g1 ¹ g2 ) ¹` ¹`

9 Suu (q) 3 Suu (g1 )Suu (g2 ) dg1 dg2 þ (2p)2

` ` Z Z

H3 (q, g1 ,

¹` ¹`

ð17Þ ¹ g1 )H3 ( ¹ q, g2 , ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2 Having determined an expression for the power spectral

lKd1 þ Ki (jq)l2 l1 ¹ a02 (jq) ¹ a1 (jq) l

2 2

Kd1 Kd3 3 (1 þ a1 q2 )2 þ (a2 q)2

` Z

þ

6 S (q) 2p uu

Suu (g) dg þ

¹`

6 (2p)2

` Z Z 3 Suu (q ¹ g1 ¹ g2 )Suu (g1 ) l1 ¹ a2 jq þ a1 q2 l2¹ ` ¹ `

lKd3 l2

`

3 Suu (g2 ) dg1 dg2 þ þ

3

lKd3 l2 l1 ¹ a2 jq þ a1 q2 l2

0 @

9 Suu (q) (2p)2 ` Z

12

Suu (g1 ) dg1 A

ð19Þ

¹`

It should be emphasised that the response spectrum is computed based on the nth-order input spectrum, which has been expressed as a function of the first-order input spectrum S u u(q) in Eq. (14), when the input is a real, stationary, zero-mean and Gaussian signal. The significance of the expressions in Eqs. (17) and (18) is that, for the first time, the power spectral density of nonlinear wave force models has been expressed in a form which explicitly shows how the power spectral density (PSD) depends on the model coefficients, in the case of a real, stationary, zero-mean and Gaussian input signal. These results provide significant insight into the operation of these systems and it is relatively easy to apply a similar analysis to other nonlinear systems. The results which were obtained by applying the probing method to the continuous-time differential equations, which define the system dynamics, are given in Eqs. (17) and (18). However the estimation procedure can only be applied for continuous-time differential equations, with no discontinuities, that is, in the case of nonlinear wave forces, only for a polynomial approximation of the drag term u(t)lu(t)l.

204

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

Table 2 Summary of results for fixed cylinder Model

Tw

gh0

gh1

Cm

Cd

KC

Re

Salford Set 1 Salford Set 2 Salford Set 3 Christchurch Bay Set De-Voorst Set

2.164 2.2109 2.2232 9.090

0.90 0.83 1.26 0.88

0.103 0.0997 0.0869 0.062

1.89 1.84 1.67 1.77

1.82 1.75 1.74 0.40

4.74 4.95 4.46 20.5

3.16 3 10 3 3.23 3 10 3 2.90 3 10 3 5.21 3 10 5

5.5555

0.87

0.0289

1.77

0.41

5.80

2.4 3 10 5

4. Spectral analysis performed on a fixed cylinder The data sets consist of velocities and forces for three different rectangular wave spectra for a fixed cylinder, obtained from the University of Salford, forces and velocities measured on the Christchurch Bay Tower and a data set from the De Voorst facility at Delft Hydraulics. In this section spectral analysis will be performed on the three model types: the Morison equation, the NARMAX parametric model and the Dynamic Morison equation. In each case the spectral estimates will be computed using two approaches. Both are based on the traditional methods of computing the spectral density but the first operates on the measured output signal and the second on the output produced by simulating the identified differential equation. The response power spectrum will be estimated for Morison and Dynamic Morison equations and a comparison will be made. 4.1. Spectral analysis of the time-domain models The power spectral density was computed for the measured output signal (the wave force) and for the output signal computed by simulating the identified continuous-time differential equations presented in Table 1. In all cases the input was a narrow band random stationary signal. The Welch nonparametric method [15] was applied to the time sequences in order to find the power spectral density. The results obtained for all five data sets are presented in Figs 1–5. The original input and output and the smooth output spectra are represented in parts (a)–(c) of Fig. 1Figs 2–5, where the dashed lines are 95% confidence bands. The confidence bands represent a precise measure of the quality of the spectrum estimate and are computed using the method presented in [16]. The dashed lines in the parts (d)–(f) are the spectra computed for NARMAX, Morison and Dynamic Morison models, while the solid lines in the same figures represent the smoothed original output spectrum. The signals are presented on a full scale in parts (a)–(c) of Figs 1–5 and only on the range relevant for the cubic nonlinearity in parts (d)–(f). In order to investigate the differences between these spectra, a nonlinear smoothing operation was applied to the high frequency bands of the original output signal. This has the effect of reducing the component of the signal which is produced by white noise. In all subsequent data sets the output

signal will be smoothed to remove the noise and both spectra from the original and smoothed signals will be presented. In parts (d)–(f) of Figs 1–5 the spectrum of the smooth measured output (parts (c)) is used as a reference and it should be compared with spectra for all three models (parts (d)–(f)). By analysing the graphs, it can be observed that the NARMAX (parts (d)) and the Dynamic Morison (parts (f)) models are almost identical in terms of spectrum shape. This fact is a new confirmation of the Dynamic Morison identification procedure, presented by Swain et al. [4]. The spectrum of the smoothed output shown in parts (c) of Figs 1–5 compares well with the spectrum produced by the noise free outputs generated from the NARMAX and Dynamic Morison equation models. The slight difference between the original output spectrum and the NARMAX estimated output spectrum corresponds to the noise model detected by the orthogonal least-squares algorithm. The NARMAX output in parts (d) is noise free, while the original output in parts (b) includes coloured noise effects. The NARMAX model (Eq. (5)) was validated using higher-order correlation tests. The correlation tests hold if the prediction error sequence e(k) in Eq. (5) is reduced to an unpredictable sequence. This will be achieved and the estimates will be unbiased if the correlated noise is correctly accommodated in the identified noise model [17]. The Morison spectrum is also different from the measured output spectrum. The difference is caused by the Morison model which cannot capture the dynamics of the nonlinear system. As discussed by Swain et al. [4], the Morison model is likely to be inappropriate at high frequency components of the input signal, where the Morison model introduces additional spectral components which are not present in the smoothed output spectrum and which appear to be spurious. This provides further evidence that the NARMAX and Dynamic Morison models are a more appropriate fit to the underlying system. 4.2. Spectral analysis via the frequency transfer functions The response spectra can also be computed by using the results in Section 3.2. This consists of mapping the Dynamic Morison differential equation models from Table 1 directly into the frequency domain using Eq. (19). The result of this analysis for all the data sets is summarised in Fig. 6, where the spectra obtained analytically from Eq. (14) (dash–dot lines) are compared to the spectra

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

obtained by simulating the Dynamic Morison differential equation in the time domain (solid lines). A measure of comparison between the power spectral density curves can be provided by computing the R and NMSE values defined below. NMSE and R values have been computed for all the data sets and the results are presented in Table 3. The NMSE is given by: q X (SFˆ Fˆ (k) ¹ SFF (k))2 X (20) NMSE ¼ (SFF (k) ¹ S¯ FF (k))2 where SFˆ Fˆ (k) is the spectrum computed by simulating one of the three models presented in Table 1, S F F(k) is the spectrum computed from the measured output and S¯ FF (k) is the mean value of the spectrum for the measured output.

205

The average power P F in the signal spectra can also be considered. 1X SFF (21) PF ¼ E[SFF ] ¼ N The normalised difference between the average power of the measured output and the computed output can then be computed as: S ¹ SFˆ Fˆ (22) R ¼ FF SFF The accuracy of the estimate S F F is directly dependent on the input spectrum estimate S u u. The estimate S u u was obtained for 8192 samples of the measured input and the Welch method of spectrum estimation assures an unbiased and very accurate estimate.

Fig. 1. Power spectral density computed for the Salford Data Set 1: (a) measured input data; (b) measured output data; (c) smooth output data; (d) NARMAX (dashed) and smooth output (solid) data; (e) Morison (dashed) and smooth output (solid) data; (f) Dynamic Morison (dashed) and smooth output (solid) data.

206

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

The processing time for the method based on Eq. (19) is determined by the third-order convolution of the input spectrum S u u in Eq. (19) which requires an exponential computation time. Therefore a trade-off between accuracy and processing time was adopted. However, a small number of samples (128 in Fig. 6) should not affect the accuracy of the estimate S F F, providing the input spectrum estimate S u u is sufficiently accurate. The advantage of using Eq. (17) or Eq. (19) is that the dependence of the spectral density on the model structure and coefficients is transparent. It is therefore very easy to visualise what a change in the model terms, the nonlinearity, the coefficients or the input spectrum will have on the

power spectrum and to assess the sensitivity to these effects.

4.3. Spectral analysis for the Morison equation The analytic expression for the response spectrum (Eq. (17)) allows an easy and simple interpretation of the wave force power spectrum for both the Morison (Eq. (3)) and Dynamic Morison equations (Eq. (7)) with nonlinearities up to the third order. The response power spectrum for the Dynamic Morison (Eq. (7)) equation has been derived in Eq. (19) and this can

Fig. 2. Power spectral density computed for the Salford Data Set 2: (a) measured input data; (b) measured output data; (c) smooth output data; (d) NARMAX (dashed) and smooth output (solid) data; (e) Morison (dashed) and smooth output (solid) data; (f) Dynamic Morison (dashed) and smooth output (solid) data.

207

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

be re-expressed as:  6 SFF (q) ¼ Suu (q) lH1 (q)l2 þ Real(H1 (q))Real(H3 (q, 0, 0)) 2p 3



` Z

Suu (g) dg þ þ

¹`

0 3@



` Z ¹`

3 2p

2

Suu (g) dgA ÿ þ

3

lH3 (q, 0, 0)l2

12 6 lH3 (q, 0, 0)l2 (2p)2

3 Suu(q ¹ g1 ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2

2 3k 6 SFF (q) ¼ Suu (q) H1 (q)þ H3 (q, 0, 0) þ lH (q, 0, 0)l2 (2p)2 3 2p ` ` Z Z

Suu (q ¹ g1 ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2

¹` ¹` ` ` Z Z ¹` ¹`

Z`

ð23Þ

Suu (g) dg. where the constant k ¼ ¹` This equation clearly shows the contribution that the linear and nonlinear terms in the model make to the response power spectrum. A similar expression for the response

Fig. 3. Power spectral density computed for the Salford Data Set 3: (a) measured input data; (b) measured output data; (c) smooth output data; (d) NARMAX (dashed) and smooth output (solid) data; (e) Morison (dashed) and smooth output (solid) data; (f) Dynamic Morison (dashed) and smooth output (solid) data.

208

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

power spectrum can be found for the Morison (Eq. (3)) equation, considering the frequency response functions derived for Eq. (3): H1 (q1 ) ¼ Kd1 þ jq1 Ki

3

H2 (q1 , q2 ) ¼ 0 H3 (q1 , q2 , q3 ) ¼ Kd3

2 3k 6 SFF (q) ¼ Suu (q) Kd1 þ jqKi þ Kd3 þ lKd3 l2 2p (2p)2 ` ` Z Z

Suu(q ¹ g1 ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2

¹` ¹`

(24)

In this case, the third-order transfer function H 3(q 1, q 2, q 3) reduces to a constant K d 3. In this case, the response power spectrum (Eq. (23)) can be simplified to:

¼ Suu (q)lH1 (q)l þ k1 Suu (q) þ k2 2

` ` Z Z

Suu (q ¹ g1 ¹ g2 )Suu

¹` ¹`

3 (g1 )Suu (g2 ) dg1 dg2

ð25Þ

Fig. 4. Power spectral density computed for the Christchurch Bay Data Set: (a) measured input data; (b) measured output data; (c) smooth output data; (d) NARMAX (dashed) and smooth output (solid) data; (e) Morison (dashed) and smooth output (solid) data; (f) Dynamic Morison (dashed) and smooth output (solid) data.

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

where k1 ¼

k2 ¼

6 K K 2p d1 d3

` Z ¹`

0 Suu (g) dg þ @

3 K 2p d3

` Z

12 Suu (g) dgA

¹`

6 2 Kd3 (2p)2

The first term term1 ¼ Suu (q)lH1 (q)l2 in Eq. (25) represents the contribution of the linear terms in Eq. (3). The second term term2 ¼ k1 Suu (q) is a fraction of the input spectrum and Z` Z` the third term Suu (q ¹ g1 ¹ g2 )Suu (g1 )Suu (g2 ) dg1 dg2 term3 ¼ k2 ¹` ¹` is the third-order convolution applied to the input spectrum respectively. These terms were estimated for the Salford Data Set 1 and are illustrated in Fig. 7.

209

Fig. 7 shows that the main contribution to the response spectrum is generated by the linear terms, especially at high frequencies. Any linear frequency function H 1(q), for which limq→` lH1 (q)l ¼ `, will produce a response spectrum with the same behaviour: limq→` SFF (q) ¼ `. For example, when a high frequency input wave with a very small amplitude (Fig. 8(a)) is applied to the Morison model, the force spectrum increases as shown in Fig. 8(b). This can have critical consequences on the fatigue design for offshore structures in which the stress is directly proportional to the wave force amplitude. Therefore the stress S is expected to increase infinitely limq→` S ¼ ` and therefore the number of cycles before failure limq→` N ¼ 0. Consequently the offshore structure will tend to be over engineered, were nonlinear effects are significant.

Fig. 5. Power spectral density computed for the De-Voorst Data Set: (a) measured input data; (b) measured output data; (c) smooth output data; (d) NARMAX (dashed) and smooth output (solid) data; (e) Morison (dashed) and smooth output (solid) data; (f) Dynamic Morison (dashed) and smooth output (solid) data.

210

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

Fig. 6. Power spectral density estimated for: (a) Salford Data Set 1; (b) Salford Data Set 2; (c) Salford Data Set 3; (d) Christchurch Bay Data Set; (e) De-Voorst Data Set [1ex] The solid lines represent the spectra from the Dynamic Morison equation calculated by standard methods, the dashed lines show the estimates obtained using equation (19).

For the Dynamic Morison model (Fig. 8(c)) the spectrum is not dominated by the linear transfer function at high frequencies and therefore it is expected to provide an improved and more reliable estimate of the fatigue life for offshore structures.

5. Conclusions An expression for the output spectrum of nonlinear wave forces has been derived based on the higher order frequency response functions and properties of the input. It was shown

Fig. 7. The force power spectrum for the Morison equation — Data Set 1 (Eq. (25)).

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

211

Fig. 8. The force power spectrum for (a) high frequency input applied to (b) Morison and (c) Dynamic Morison equations — Data Set 1

how this can be evaluated for the Dynamic Morison and Morison equations to reveal the explicit dependence of the spectra on the coefficients in the nonlinear differential equation. An analysis and interpretation of the Morison equation response spectrum was also presented. Similar results can easily be derived for other nonlinear differential equations or difference (NARMAX) equation models by evaluating the higher order frequency response functions for these and substituting them into the expression for the output spectrum. These results provide a unique insight by explicitly showing how the time-domain coefficients influence the output spectra and vice versa. Application of the results was demonstrated using five data sets from Salford, De-Voorst and Christchurch Bay. It is realised that it is not straightforward to directly incorporate these theoretical improvements into existing design procedure where specialist practices have been evolved in relation to wave forces (see, for example, Ref. [1]). Also, while the values of the nondimensional coefficients in the Dynamic Morison equation are consistent with their Table 3 R and NMSE values for the estimated spectra using Equation (19) and the Dynamic Morison model (Eq. (7)) Data Set

R

NMSE

Salford Data Set 1 Salford Data Set 2 Salford Data Set 3 Christchurch Bay Data Set De-Voorst Data Set

9.98% 12.66% 6.03% 9.48% 3.91%

22.77% 14.74% 16.42% 16.18% 28.86%

physical origin, their magnitudes for wide parametric ranges require further evaluation (of high-quality data sets). However as offshore structures move to even-deeper waters, high-frequency forcing assumes increasing significance and methods of analysis, like those presented above, become more important.

Acknowledgements Thanks are due to Delft Hydraulics for access to wave force data obtained in the Delta flume in the De Voorst facility. This data was originally provided by M. Birkenshaw of the Health and Safety executive and Prof. Bearman of Imperial College.

References [1] Lin YK. Probabilistic Theory of Structural Dynamics. New York: McGraw-Hill, 1967. [2] Worden K, Stansby PK, Tomlinson GR, Billings SA. Identification of nonlinear wave forces. Journal of Fluids and Structures 1994;8:19–71. [3] Stansby PK, Worden K, Billings SA, Tomlinson GR. Improved wave force classification using system identification. Applied Ocean Research 1992;14:107–118. [4] Swain AK, Billings SA, Stansby PK, Baker M. Accurate prediction of nonlinear wave forces: Part I (fixed cylinder), Journal of Mechanical Systems and Signal Processing, 1998, to be published. [5] Baker M. Wave loading on a small diameter flexibly mounted cylinder in random waves, Ph.D. Dissertation, University of Salford, 1994. [6] Morison JR, O’Brien MP, Johnson JW, Schaf SA. The force exerted by surface waves on piles. Petroleum Transactions 1950;189:189–202.

212

O.M. Boaghe et al. / Applied Ocean Research 20 (1998) 199–212

[7] Palm G, Poggio T. The Volterra representation and the Wiener expansion: validity and pitfalls. SIAM Journal on Applied Mathematics 1987;33(2):195–216. [8] Bendat JS, Piersol AG. Decomposition of wave forces into linear and non-linear components. Journal of Sound and Vibration 1986;106(3):391–408. [9] Leontaritis IJ, Billings SA. Input–output parametric models for nonlinear systems; Part 1: deterministic nonlinear systems. International Journal of Control 1985;41:303–328. [10] Sarpkaya T, Isaacson M. Mechanics of Wave Forces on Offshore Structures, London: Van Nostrand Reinhold Co., 1981. [11] Stansby P, Isaacson M. Recent developments in offshore hydrodynamics, a workshop report. Applied Ocean Research 1987;9:118–127. [12] Volterra V. Theory of Functionals and of Integral and Integrodifferential Equations. New York: Dover, 1959. [13] Rugh WJ. Nonlinear Systems Theory — The Volterra/Wiener

[14]

[15] [16]

[17]

[18]

approach. London: The Johns Hopkins Press Ltd, 1981, pp. 213– 223. Peyton-Jones JC, Billings SA. Mapping non-linear integro-differential equations into the frequency domain. International Journal of Control 1990;52(4):863–879. Proakis JG, Rader CM, Ling F, Nikias CL. Advanced Digital Signal Processing. New York: Macmillan, 1992. Wellstead PE, Spectral Analysis and its Applications, Lecture Notes in Control and Information Sciences: Signal Processing for Control, Berlin, Springer-Verlag, 1986, pp. 210–245. Korenberg M, Billings SA, Liu YP, McIloy PI. Orthogonal parameter estimation algorithm for non-linear stochastic systems. International Journal of Control 1988;48(1):193–210. Swain AK, Billings, SA. Weighted complex orthogonal estimation for identifying linear and nonlinear continuous time models from generalised frequency response functions. Journal of Mechanical Systems and Signal Processing 1998, to be published.