- Email: [email protected]
Wave predictions based on scatter diagram data. A computer program package
Advances
in Engineering
Software 23
(1995)
49-59
6 1995 Elsevier Science Limited Printed in Great Britain. All rights reserved 0965-9978/95/$09.50
0965-9978(95)00019-4
ELSEVIER
Wave predictions based on scatter diagram data. A computer program package R. Capitio National
Laboratory
of Civil Engineering
(LNEC),
Av. do Brusil, 101, 1799 Lisbon, Portugal
& R. Burrows Department
of Civil Engineering,
University qf Liverpool,
P.O. Box 147, Liverpool
L69 3BX, UK
(Received4 November 1994;revisedversion received7 April 1995;accepted27 April 1995) This paper describesa user-friendly computer program packagecontaining an approachfor the prediction of extremewavesbasedon the extrapolation of the scatter diagram of significantwave height H, and meanzero-upcrossingperiod T,. Following the validity testsof data input, the computerprogram takes into account the existenceof missingdata and calmsand includesroutines for the treatment of various theoreticaldescriptionsof H,. Sensitivity studieshave been completed to appraisethe changesin wave predictions due to missingdata consideredeither as calms or instrument failure and due to inadequaciesof measurementsystemsto record low T, values.The possibleinfluenceof currents on extreme wave height predictions has also been taken into account. For examplepurposes,the computer program is usedfor the prediction of extreme wave heights in a location (Morecambe Bay, UK) where wave data (scatter diagrams)are available. Key words: wave prediction, extremes,scatterdiagram, oceanengineering
Skew X StDevX
NOTATION A, B, C cdf, pdf
dH dHr
dT E{Hmax} H H max
HZ 1H IHr IT
N nH nHr nT R SDError
TZ Wtj
location, scaleand shapeparametersof a cdf cumulative and probability density functions incrementof H in a scatterdiagram incrementof extrapolated H in a scatter diagram incrementof T in a scatter diagram expected maximum wave height in a 3-h record period wave height, individual wave height designmaximum wave height significantwave weight lower value of H in a scatterdiagram lower value of extrapolated H in a scatterdiagram lower value of T in a scatterdiagram number of observations, samplesize, number of adopted data number of classesof H present in the scatter diagram number of extrapolated classesof H in a scatter diagram number of classesof T present in the scatter diagram return period standarddeviation error
VarX
x, y AvX Notation
skewness of a sampleof values standarddeviation of samplevalues meanzero-upcrossingwave period number of occurrencesin the class (i,j) of the scatter diagram variance of a sampleof values orthogonal axes,co-ordinates meanvalue of x. in the Computer Program
Grafyx, Graph, Utility Gamma HofP Ord PlotPos ProbX XAxis YAxis DetectCurrents DetectMissData DistParameters 49
LTS (PROGRAMS)
external units of routines Gammafunction inverseprobability distribution function Y-variate function on probability paper units plotting position functions probability distribution function X-variate function on plotting units Y-variate function on plotting units detecting possibleexistenceof currents detectingmissingvaluesin data file calculation of distribution parametersgraphically
50
Distribution DrawLine DrawTable ExtOptions GraphicParameters Graphics InputData InterExtra LeastSquares Moments Options output 1 output2 PDFandCDF PlottingPositions Probabilities Ranks ReadFile Scales ScaleYAxis Statistics TestFitting
R. Capitcio, R. Burrows
inputting data from keyboard plotting the fitted line on probability paper plotting a table of results on probability paper inputting data from keyboard inputting data from keyboard creation of graphics on screen and plotter file inputting data from keyboard calculation of R for given H method of least squares method of moments inputting data from keyboard first writing of results on results file appending results on the created results file pdf and cdf inputting data from keyboard calculations with probabilities ranking of a sample of values reading of the data file calculations of plotting scales calculation of the scale for the probability axis sample statistics testing the data fitting
INTRODUCTION For the design of an offshore or coastal structure the engineer must consider the extreme conditions of wave climate likely to occur during the design life of the structure. Extreme wave conditions are usually represented in terms of wave heights and wave periods. These parameters are commonly obtained by analyzing wave records which contain wave characteristics measured by data acquisition systems. The wave records are usually composed of lo-20min of measurement and obtained once every 3 h. In the short-term (that is, in periods of time up to, say, 3 h) the sea state intensity can be assumed to be constant and thus each wave record is regarded to be typical of the whole 3 h-period. During the wave recording period it is assumed that the statistical properties of the sea state do not vary and hence the characteristics of a uniform stochastic process can be applied to it. The short-term wave climate is well understood. The wave heights follow the Rayleigh probability distribution if one assumes a narrow-banded spectrum and considers that the short-crested random sea elevations might be modelled as a sum of a number of linear, sinusoidal, waves travelling in different directions and having random phases. Central limits theorem shows that as the number of wave components and directions tend to infinity, the probability distribution of the sea elevation tends to the Gaussian distribution and characteristics of this distribution are well known.
For the long-term wave climate analysis each wave record of lo-20min is usually reduced to two parameters, the significant wave height, Hs, and the mean zero-upcrossing period, Tz, that are supposed to represent the whole period of 3 h. These parameters can be calculated in a number of ways being at present normally found by spectral analysis of the time series of the sea surface elevation or, sometimes, by using the Tucker-Draper method.’ In addition to Hs and Tz given by the scatter diagram, another representative parameter of the sea state conditions, the expected maximum wave height - E{H,,,} - can be computed for each 3 h observation as: E{H,,,}
z 0.707Hs(log,
n)tj2 and y1= (3 x 3600)/T,
Collection of these parameters during a long period (at least 1 yr) gives to the engineer an enlarged view of the wave climate. For engineering purposes, it is useful to use the wave climate information (Hs, Tz and E{H,,,}) in order to find the long-term probability distribution of wave heights, and thus extrapolate the observations to find the conditions of design severity. This is done by use of the theory of extremes, if the conditions of application of the theory are fulfilled, or, as a last resort, by use of an empirical procedure. Both procedures use a selected cumulative probability distribution which normally belongs to the family of extremes. This is fully described in a companion article to this one.2 Extreme wave statistics can be described by several exact probability distributions. These include the extremal Fisher-Tippet distributions (Types I, II and III) for maxima and minima and the normal and lognormal distributions. Such distributions are based on firm theoretical foundations and thus data fitting should be done in a theoretical basis as well. However, the available data often does not fulfil the conditions of a random sample so that the data fitting is often empirical. Despite this, experience has shown that extrapolations made under the use of such empirical methods are reasonably accurate. Very often, an approximate representation of the wave climate is given by a bivariate histogram of Hs and Tz. The bivariate histogram, also called scatter diagram or contingency table, represents the distribution of sea state conditions during the period when wave records have been produced. For each record, the expected maximum (the most probable) wave height, E{H,,,}, can be determined if it is assumed that wave heights are Rayleigh distributed, which is a reasonable assumption. The particular case of scatter diagram information and comments on the determination of E{H,,,} are dealt with herein. Design wave conditions are usually abstracted by using procedures for the extrapolation of the significant or maximum wave heights. In particular, one empirical procedure for extrapolation of the marginal cumulative
Wave predictions
51
based on scatter diagram data
Table 1. The files required to run LTS
LTS.EXE EGAVGA.BGI CGA.BGI IBM8514.BGI LITT.CHR file.DAT
Executable File One of the Graphics Interface Files File containing Stroked Fonts The Data File (compulsorily with extension DAT)
distribution of H, is described. However, for some applications, offshore structure design for instance, it is more important to find such conditions based on the individual wave height. This can be done approximately by extrapolating the values E{H,,,} of each record or by use of Battjes’ procedure which accounts for the complete distribution of individual waves. A more detailed review of this background is given in Capitao3 and Capita0 & Burrows.* In order to use all the procedures in a logical and integrated way, a computer program able to run on personal computers has been written using the Pascal programming language. This program is briefly described below.
THE COMPUTER
PROGRAM
LTS
The computer program LTS is briefly described here in the order that the program routines flow. The package entails a PC/XT, PC/AT, PS/2 or compatible microcomputer using IBM/MS-DOS 3.2 or above. In this version, use of a 8087, 80287 or 80387 mathematics coprocessor and a minimum of 256 KByte RAM memory are necessary. The output can be directed to a line printer and/or to a plotter, so at least one of these devices should be provided. It is recommended that this package should be installed on a hard disk since the execution time can be considerably shortened. It is important as well that some free space on disk should be allowed for the creation of the results’ files. The program runs simply by executing the command LTS. However, for the graphics’ presentation on screen, graphics interface files are needed. Depending on the graphic card and monitor present in the computer, one of the following files is required: EGAVGA.BGI, HERCULES.BGI, CGA.BGI or IBM8514.BGI. Since the program uses stroked fonts for the graphics’ presentation on screen, it is necessary to use the file LITT.CHR as well. Briefly, the files required to run the computer program LTS are given in Table 1. Since the main program calls routines that belong to units of routines, the executable files containing these units are needed for the program compilation. Units are the basis of modular programming in Turbo Pascal since they are used to create libraries and to divide large programs into logically related modules.
The required units have their executable files named as: GRAPH/TPU,4 GRAFYX.TPU’ and UTILITY.TPU.3 The LTS program can be made available, at nominal charge, from the authors. Brief description
of the computer program
The computer program LTS has been written in the Pascal programming language, using Turbo Pascal (Version 6.0) environment of Borland International.4 Basically, a Pascal computer program is constituted by four major blocks:
6) The block containing the
dejinition and declara-
of units, labels, constants, types and variables of the main program. (ii) The block containing thefunctions. A function is a program part that computes and returns a value and is activated by the evaluation of function call in an expression. (iii) The block that contains the procedures. A procedure is a program part that performs a specific action, often based on a set of parameters and is activated by a procedure statement. (iv) The block containing the main program, where the organization of the program is built. tion
The computer program is herein described in the order that functions and procedures are called throughout the program run. A flow chart of the program is given below. The main program calls the following procedures: (A) InputData
This procedure permits the input of data via keyboard. The input values are: Data File Name - Name of the file that contains the data to be analyzed by the computer program. The extension must compulsorily be .DAT; Results File Name - Name of the file that will be created with the results of the program run. By default it takes the data file name with extension .RES; Plotter File Name - Name of the file that will be created for later sending to the plotter printer. This file will be created using the plotter programming language HPGL (Hewlett Packard Graphics Language). By default, the plotter file name takes the data file name with extension .PLT;
52
R. Capitiio, R. Burrows
I
Ext options
Read file I
Q
Least squ
Plot position I Distribution I
I
Draw line
I
output 2
Fig. 1. LTS flow chart. Name of the units for both axes X and Y (variate and return period) - By default they take the names ‘metre’
and ‘year’ respectively; Number of observationsper unit of return period (year)
- It is assumed by the program that, in the case of bivariate data, the observations are 3-hourly taken and hence 2920 values are expected during 1 year. For the case of the univariate data, by default, it is assumed that there is one observation per year (annual observations); Number of units (years) of observed data - By default, the number is 1.0.
(B.2) Ranks. Ranks the sample values of the univariate sample. It is assumed that the bivariate sample (scatter diagram) is already ranked in both directions. (B.3) PlottingPositions. When the required
type of analysis is univatiate, this procedure uses any of the most common plotting positions, such as i/(N + l), i/N, etc., for the assignment of probabilities to the ranked data values. For the case of bivariate analysis with grouped data, use of only the Weibull plotting position is allowed for such a purpose.
(B) ExtOptions This procedure permits the input of the following data via keyboard (the entries in this procedure are only needed for bivariate analysis): Option between run for Hs or H; Option between whether or not extrapolation of the marginal cumulative distribution is required.
This procedure calls two other procedures: (B.1) ReadFile. This procedure reads the data file. The data file can be of two different types: univariate data (can be, for instance, annual maxima of Hs or H,,,) or bivariate data (scatter diagram of Hs and Tz, for example).
(C) Distribution In this procedure, the user can select any of the available probability distributions for the fitting of the data. The available distributions (at the moment) are: Gumbel (Fisher-Tippet Type I for maxima); Weibull (Fisher-Tippet Type III for minima); In the case of the three-parameter distribution (Fisher-Tippet Type III), the location parameter A is required if graphical extrapolations are used. If the method of moments is used for the estimation of the parameters distribution, such a requirement is not necessary since all the parameters are calculated iteratively.
Wave predictions
based on scatter diagram data
(D) Probabilities According to the required type of analysis, this procedure permits the formation of the array of wave heights and mean zero up-crossing period. This procedure calls two other procedures: (D.1) Statistics. This procedure calculates some important sample statistics like the mean sample value (AwX), the Variance ( VarX), the Standard Deviation (StDewX) and the Skewness (SkewX). This procedure calls the procedure: (D.1.1) Moments. For the estimation of the parameters distribution, the method of moments is herein applied using a direct calculation for the case of FT Type I distribution and an iterative process for the case of the FT Type III distributions.
The pdf and cdf (raw and extrapolated marginal probabilities) of both significant and individual wave heights and the expected value E{H,,,} are calculated here. (0.2)
PDFandCDF.
(E) Output1 This procedure opens the results file and writes on it the statistics calculated in procedure (D.l) and the pdf/cdf calculated in (D.2). (F) GraphicParam More data is read from the keyboard if graphical procedures are required. The following entries are required for the plotting of the graphics in both screen and plotter printer: Minimum and maximum abscissa for graphics - Default values are provided: Interval between X-axis graphic units - By default is 1; Number of ticks between X-axis graphic units - By default is 10; Required number of return periods - The program calculates the return values of H or Hs corresponding to each required return period. By default, the number of required return periods is 0. (G) TestFitting To decide which distribution fits best to the sample data, calculation of the correlation factor r and standard deviation error, SDError, is done herein. These factors can only serve as an approximate indicator of the goodness of fit of the selected distribution. Depending on the selected options, this procedure might call three other procedures. (G.1) Scales. Some distributions need to have a nonlinear variate axis. Thus, a transformation of the xcoordinate axis is required for the plotting on the screen
53
or plotter printer. According to the selected distribution, this procedure scales up that axis. (G.2) LeastSquares. This procedure is called when the method of least squares is required. In order to find the best fit of a straight line through the points [H, P(H)], the standard least squares method is herein applied and thus the parameters of the straight line are found. (G.3) DistParameters. If the least squares method is required, this procedure, based on the parameters of the fitted straight line, calculates the parameters of the selected distribution.
(H) Graphics The graphic routine is designed to present, through the screen and plotter, a graph showing the straight line fitted through the data points, according to the selected distribution. This routine tests the monitor type and, by default, creates a pre-defined graphic layout. At the present, the routine permits to scale up the probability axis from 0.000001 to 0.99999999999 but this range can be extended (reduced) if required. This procedure calls three other procedures. (H.l) Scale YAxis. Graphic routine destined to scale up the probability axis using both pixel units (for screen output) or plotter units (for plotter output). (H.2) DrawLine. Graphic procedure that permits the drawing of the fitted line through the data points. (H.3) Draw Table. Graphic routine that permits the creation of a table containing the predictions of wave height, based on the fitted line, for a given return period.
(I) output2 This procedure re-opens the results file and appends on it the results of the estimation of distribution parameters, correlation coefficient, standard deviation error (graphic procedures) and return values for the required return periods. (J) Options This procedure contains a menu of options that permits the user to return to other menus, run the program again with the same data or different data, changing of graphic parameters or plotting positions (only in the case of univariate analysis). This procedure calls another procedure. (J.1) InterExtra. If required, calculation of the return period R for any given value of Hs is completed.
54
R. Capitdo, R. Burrows
Input to computer program LTS
done for both significant wave heights and individual wave heights (by using the Battjes’ method6). Samples of the input data report and tabular results arising are given in Figs 2 and 3. Plots of graphical results for both significant wave height and individual wave height, each in Gumbel and Weibull formats, are shown in Figs 4-7. An important observation in relation to the individual waveheight (H) data is that the actual points plotted represent the true distribution.2 The fitted lines, by least squares here, serve as approximations only and, in the case of the Gumbel fits, can depart markedly from the true values in the upper tail. The Gumbel distribution, therefore, gives a poor representation of the distribution of individual wave heights, whilst the Weibull distribution, with the flexibility of an extra parameter, provides acceptable fits. The sample data has been studied for the purpose of determining how missing data, considered either as calms or non-measured data, can influence the wave predictions. Furthermore, an assessment has been made to establish, for this location, whether inability of the measurement systems in recording low Tz values can affect the same wave predictions. It was observed from analysis of the results that predictions of significant wave height found by using the method of moments were slightly greater than those obtained by using the least squares method as can be seen from Table 2. The systems that exist for the measurement of wave data have limitations and the accuracy of such a measurement can vary substantially from system to system. It is important thus to be aware of the effect on wave predictions when a particular measurement system is unable to measure a certain wave characteristic. One common problem is the difficulty of measuring low values of Tz. Low values of Tz are usually correlated with low values of Hs. Consequently, on average, if one does not consider low values of T, that means the predictions of wave height will be higher than those calculated if all Tz values are used. In order to illustrate the differences in predictions, three data samples from Morecambe Bay have been considered. The first sample (all data with Tz readings beginning at 1.75s) is presented in Table 2. The second and third data samples have been created upon the first one but assuming that the measurement system was unable to record waves having Tz values less than 2.75
The computer program only needs input data from the keyboard (or default values) and a file containing the data sample. The layout of such a file depends on the type of analysis required. For an univariate sample (for instance, annual maxima), the data file is simply constituted by a column which contains the data values x1, x2, . , x,,: ( beginning of data file ) First extreme value Second extreme value Third extreme value
Xl x2 x3
( end of data file ) The data file layout for a bivariate diagram) is shown below.
nib extreme value
sample (scatter
Output of the computer program Three output files are created computer program LTS:
when running
the
(i) one file that is dumped to the screen; (ii) one text file (ASCII format) containing the results ready to be printed on a line printer (by default, with extension RES); (iii) one plotter file that can later be sent to the plotter printer (by default, with extension PLT). This file is written in HPGL programming language, so any plotter accepting the commands HPGL is suitable. Generally, most commercial plotters accept this language. APPLICATION
TO WAVE DATA
All scatter diagram data Scatter diagram data collected at Morecambe Bay in the British Isles have been considered in this study. Wave heights predictions associated with return periods of 10, 50 and 1OOyr have been calculated by the computer program LTS on both Gumbel and Weibull probability papers, fitting by using both the least squares method and method of moments. These calculations have been
( beginning of data file ) nH nT nHr w(l,l) w(2,l) ..
w(nH,l)
1H 1T 1Hr w(l,2) w(2,2)
dH dT dHr . .. .. . .. ..
41, nT) w(2,nT) I . .
w(nH,2)
. .
Number, lower value (class average) and increment of H (or Hs) Number, lower value (class average) and increment of Tz Number, lower value and increment of extrapolated H (or Hs) Occurrences of the Is’ class of H Occurrences of the Znd class of H .. w(nH,nT)
( end of data file )
Occurrences of the nth class of H
Wave predictions
based on scatter diagram data
tt*******************************************.* * LONG-TERM STATISTICS
55
+
t******t**************t*************************
******t***********t****t******tttt*tt****~~* GENERAL DATA : **t********i********t*tt********~***~*,**~********~ TITLE
: Morecambe
Descyription
: All
Data
file
Bay
Data
: MCBl.DAT
Results
file
: MCBl.RES
Platter
file
: MCBl??.PLT
Wavr
Height
unities
( H )
Wave
Period
unities
I TZ
Return
Deriod
Number
Iof
observations
L!nities
Numl~er
of
years
Number
of
.._...._
second
( T
1 .__..._
year
cause
3.00
. .._....
data
Minlmum
value
of
sample
.._...._.._
Maximum
value
of
sample
.
Value
7961
. . .._...._
Missing
Mean
year
..____..___.
per as
metry
)
observed
ohs.
__..__.._
of
2920.0
.
non-was.
799
.
0.50 6.00
.._...._
.._____._.___..____.____.
1.0724
Standard
Deviation
Variance
. . ..__.....__...____..__.._
0.5R73
Skewness
___._..__...____...___..___
1.3430
SCATTER
DIAGRAM
_._._____......._
0.7663
DATA
H'S 0.15 0.75 1.25 1.75 2.15 2.75 3.25 3.75 4.25 4.75 5.25 5.75 6.25 6.75 7.25
2.0 120 73 48 0 0 0 0 0 0 0 0 i 0 0
2.5 596 316 58 0 0 0 0 0 0 0 0 0 0 0 0
3.0 537 851 69 0 0 0 0 0 0 0 0 0 0 0 0
3.5 289 640 313 9 0 0 0 0 0 0 0 0 0 0 0
4.0 113 312 499 86 4 0 0 0 0 0 0 0 0 0 0
4.5 45 172 327 283 52 2 1 0 0 0 0 0 0 0 0
5.0 33 103 192 270 117 35 2 0 0 0 0 i 0 0
PERIOD 5.5 34 89 86 130 160 84 16 : 0 0 0 0 0 0
6.0 17 49 57 68 58 53 46 8 18 0 0 0 0 0 0
6.5 21 30 29 33 22 39 45 22 113 110 0 0 0 0
7.0
7.5
8.0
9 9 3 12 15 11 10 10 3
6 3 0 3 3 11 8 4 6
0 0 0 0
0 0 0 0
8.5
1 3 1 0 10 3 3 10 1 3 0 0" 0 0
9.0 R 3 0 1 0 2 3 0' 0 0 0 0 0
9 0 0 0 0 0 0 0 0 0 0 ; 0 0
Fig. 2. Sample input data report.
and 3.25 s, respectively. Wave predictions based on these two samples are presented on Tables 3 and 4. All the predictions were found by considering the missing data as non-measured data. For the Gumbel distribution it was found that when low values of T, are not taken into account, the wave predictions become slightly higher than those obtained when using all Tz data. However, if the data fitting is made using the Weibull distribution, the wave predictions are lower when low values of Tz are not
considered. Thus when acquisition systems are unable to read low values of T,, a slight overestimation or underestimation in wave predictions will occur if one fits the data points to the Gumbel or Weibull distribution, respectively. The effects are, however, very small, as seen from the comparison of Tables 2-4. It is now necessary to consider analysis of a scatter diagram with missing data to appraise the consequences in the context of wave predictions. Missing data can be considered as the result of calms (not measured by the
56
R. Capitdo, R. Burrows TABLE OF RESULTS ----------------------------------------------------------~ HsI Hsi Wi SWi
PDFHs
CDFHs
PDFHsX
CDFHsX
-!m!___!m!--!"o:_"ccUrl____lI1______II)_-----!~!------!~!-0.00 0.25
lR38
1838
0.75
2653
4491
1.25
1682
6173
1.75
895
7068
2.25
432
7500
2.75
238
7738
3.25
133
7871
0.50 1.00 1.50 2.00 2.50 3.00 3.50 3.75
58
7929
4.25
21
7950
4.00 4.50 4.15
8
7958
5.25
2
7960
5.75
1
7961
5.00 5.50 6.00
GUMBEL DISTRIBlJTION *ttf****f*t****ttt**ttt**~************************ SIGNIFICANT WAVE HEIGHT P(x)
= exp
I -exp
( -(
0.23088
0.23085
0.23142
0.23142
0.33325
0.56405
0.29913
0.53055
0.21128
0.77531
0.22939
0.75994
0.11242
0.88772
0.12796
0.88791
0.05426
0.94197
0.06191
0.94981
0.02990
0.97187
0.02813
0.97795
0.01671
0.98857
0.01244
0.99039
0.00729
0.99586
0.00544
0.99583
0.00264
0.99849
0.00236
0.99819
0.00100
0.99950
0.00103
0.99922
0.00025
0.99975
0.00044
0.99966
0.00013
0.99987
0.00019
0.99985
WEIBULL DISTRIBUTION *t******t*tt********t*****+tttt*t*t******************** SIGNIFICANT WAVE HEIGHT X - A J/B
)
]
P(x)
= 1 - exp
----------LEAST
SQUARES
41s Bls
= =
METHOD
0.77161 0.57908
I
METHOD
I I
Amm= Bmm=
OF
MOMENTS
r
= 0.99843 I SDError = 0.16721 I --------_------~-_~-____________________-----~~ --------
Extrapolation
of
SQUARES
Als Bls ClS
0.59751
= = =
T 10 50
100
.
Hs(LS1
.
6.73 7.66 8.06
Hs(MM)
. Hmax(LS1.
6.87 7.83 8.25
12.64 14.33 15.05
[ 12.00 i 23 1
--------
of
Extrapolation
12.90 14.64 15.39
. T . H . ____-----___------__--------------------------10 15.95 17.55 18.25
0.08321
I
Amm=
0.08321
1.09011 1.32633
I I
Bmm= cm=
1.07140 1.30186
0.99940
I
0.03518
I
Extrapolation
of
OF
Observations
.
T 10 50
.
Hs(LS)
.
6.40 7.13 7.44
Hs(MM)
. Hmax(LS).
6.50 7.26 7.58
12.05 13.38 13.93
WEIBULL DISTRIBUTION *t***tt*****t*t***tt****************~***~***,** INDIVIDUAL WAVE HEIGHT Limit of Extrapolation [ 12.00 Number of Classes I 23 1
1
Observations
METHOD
MOMENTS
--------
Hmax(MM).
100
INDIVIDUAL WAVE HEIGHT Limit of Extrapolation Number of Classes
----------_
parameters
METHOD
1
---------------
.
/ B 1°C
I
r I =
SDError Observations
I 1X - A)
Distribution
LEAST
0.72756
[ -
--------
--------
Extrapolation
-____-----__----------------------------------. T . H . 10
14.41 50 15.90 100 16.52 _----_____------___---------~~~~~~~~~--------~~
Fig. 3. Sampletabular resultsreport.
of
Hmax(MM). 12.23 13.61 14.18
]
Observations
--------
Wave predictions
based on scatter diagram data
MB -Miss. data: CALMS Gumbel distribution Significant wave heights
99.9990
t
+
+
Least sq-s -2.329 0.996
/ A= B=
/
Wave
height
(m)
Wave
height(m)
Fig. 4. Gumbel plot of significant wave heights, H,.
Fig. 6. Gumbel plot of individual wave heights, H.
acquisition system because of the very low values of Hs) or as non-measured data resulting from randomly occurring malfunction of the system. In this last case, it is important to note that the predictions can be seriously affected if occurrence of malfunction is biased to occur when the sea state is rough. In this case the data set would be unacceptable for accurate wave predictions. In order to compare wave predictions when considering missing data either as calms or non-measurements due to malfunction of the acquisition system, Tables 5 and 6 have been produced for Morecambe Bay.
From the comparison of Tables 5 and 6 it can be concluded that the decision to consider missing data as non-measured data instead of as calms produces slightly higher predictions of wave heights as might be anticipated. In other localities the occurrence of ‘calms’ will be more prevalent and the consequences of their misrepresentation will be much more significant. CONCLUSIONS Scatter diagram data have been used for predictions of wave heights by use of the computer program LTS.
MB - Miss. data: CALMS Weibull distribution Individual wave heights
Weibull distribution Significant wave heights
30.0000
t
HI101 = H[501= H(lOOl=
Least squares
Moments
0.083 1.090 I.326
0.083 I.071 1.302
6.36 7.09
6.45 7.21 7.53
7.40
2o.ml
22, .e 4 % &
Leas1 sq”areS
9o.ooooooO-
8O.OOOOOOO 7o.OOoOOO0,
Wave
height
(m)
Fig. 5. Weibull plot of significant wave heights, H,.
Wave
A= LI= c=
0.083 0.542 0.888
H [lo] =
14.39
H 150) = H[lOO]=
15.82 16.44
I 6
I II
height
i C
I J16 21
(m)
Fig. 7. Weibull plot of individual wave heights, H.
58
R. Capitfio, R. Burrows Table 2. Wave height predictions -
Morecambe Bay (all data)
R
f&m
W
Gumbel
10 50 100
Weibull
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 7.7 8.1
6.9 7.8 8.3
6.4 7.1 7.4
6.5 7.3 7.6
12.6 14.3 15.1
12.9 14.6 15.4
12.1 13.4 13.9
12.2 13.6 14.2
Table 3. Wave height predictions &I
z
Gumbel
10 50 100
- 3h (ml
Weibull
Weibull Battjes6
16.0 17.6 18.3
14.5 15.9 16.5
Morecambe Bay ( Tz > 2.75 s) ffmax - 3h W
Gumbel
Weibull
6s
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.8 7.7 8.1
7.1 8.0 8.5
6.4 7.1 7.3
6.5 7.2 7.4
12.8 14.5 15.2
13.2 15.0 15.8
12.0 13.2 13.7
12.1 13.4 14.0
Weibull
Battjes6 16.2 17.9 18.6
14.4 15.8 16.4
The number of Tz values >2.75s is 6750, compared with 7961 (all data). Table 4. Wave height predictions -
Morecambe Bay (T’ > 3.25s) fLx
&I
Gumbel
10 50 100
- 3h (ml
z
Weibull
(:I
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
7.0 7.9 8.3
7.3 8.3 8.8
6.4 7.0 7.2
6.3 6.9 7.2
13.0 14.7 15.4
13.7 15.5 16.3
12.0 13.1 13.6
11.8 13.0 13.4
Weibull
Battje@ 16.7 18.3 19.0
14.1 15.4 15.9
The number of T, values>3.25 s is 5293,comparedwith 7961(all data). LSM - Distribution parametersestimationby Least SquaresMethod. M - Distribution parametersestimationMethod of Moments. Battjes’ ( 1970)6procedurefor the prediction of individual wave heights.3
Several conclusions can be drawn from comparisons between the results obtained:
(1) Acquisition
systems that are inadequate to register low values of Tz are likely to overestimate predicted values of wave heights if the Gumbel distribution is used, whereas underestimation may occur if fitting is made using the Weibull distribution, according to the sample application made herein. (2) With missing data in the scatter diagrams, the decision to consider this either as ‘calms’ or ‘non-measurements’ leads to different predicted values of wave height. For design purposes it is important to know what has really happened with those missing data. When in doubt, it is better to consider the missing data as non-measurements since it is a conservative decision. (31 Using the ‘weighted least squares method’, a greater importance is given to the central range of points than using the method of moments and wave height predictions are lower. It is not possible to state categorically that one probability distribution is better, that is, represents better the data set, than another. The Weibull
distribution is here shown to be more flexible in its fit to a particular set of data since it contains three parameters whereas the Gumbel distribution has only two. However, in practice, the Gumbel distribution has generally been accepted as the principal extreme value distribution for the current application. It is, in fact, a distribution that generally applies very well to wave data and the estimation of parameters is straightforward. In this study, long-term variability of the characteristics of the wave climate has not been a consideration. This represents an important simplification of the problem. Appreciation must be given to the fact that scatter diagram data records covering perhaps only 1 or 2yr may not be ‘typical’ of the true long-term wave climate and some intuitive adjustment may be necessary if, for example, wind records show that the wave climate sampled may have been more, or less, severe than the average. An improvement in the quality of analysis of extreme wave data is achieved if emphasis is placed on the effort of collecting as many wave data as possible and over a period as long as is feasible. This enables the engineer to have a better perception of the true wave climate that exists at the site and the danger of selecting a
Wave predictions Table 5. Wave height predictions -
Morecambe Bay (missing data as non-measured data) f&m
- 3h Cm)
i2
Weibull
Gumbel 10 50 100
59
based on scatter diagram data
Weibull
Gumbel
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 1.7 8.1
6.9 7.8 8.3
6.4 7.1 7.4
6.5 I.3 7.6
12.6 14.3 15.1
12.9 14.6 15.4
12.1 13.4 13.9
12.2 13.6 14.2
Weibull Battjes6
16.0 17.6 18.3
14.5 15.9 16.5
Calms are 799 out of 8760. Table 6. Wave height predictions CRY4
Morecambe Bay (missing data as calms)
2
ffmax
- 3h
(ml
Gumbel 10 50 100
Weibull
(:I
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 7.6 8.0
6.8 7.8 8.2
6.4 7.1 I.4
6.5 7.2 I.5
12.5 14.2 15.0
12.8 14.5 15.3
12.0 13.3 13.9
12.2 13.5 14.1
Weibull Battjes6
15.9 17.5 18.2
14.4 15.8 16.4
Non-measurements are 799 out of 8760. LSM - Distribution parameters estimation by Least Squares Method. M - Distribution parameters estimation by Method of Moments. Battjes’ (1970)6 procedure for the prediction of individual wave heights.’
distribution that does not belong to the population will certainly be minimized. The computer program must be used carefully. It is important that the user is aware of the whole process of predicting wave characteristics. In particular, the user must know what is the reliability of the data sample that is going to be used since enormous mistakes can arise from data samples not containing a reasonable set of extreme values.
REFERENCES 1. Draper, L., Derivation of a ‘Design Wave’ from instrumental records of sea waves. Proc. Z.C.E., 1963,26, 291l 304.
2. Capitao, R. & Burrows, R., Wave predictions basedon scatter diagram data: A review. Adv. Engng Sqftwure, 1995,23, 37-47. Capitao, R., Wave predictions basedon scatter diagram data. MSc(Eng.) Thesis, Dep. of Civil Engng, Univ. of Liverpool, 1992. Borland International, Turbo Pascal -- Version 6.0. Library Reference and Programmers’ Guide. Borland International, Scotts Valley, CA. Palma, J., Bibliotecas em Pascal para Constru9Bo de Graficos Cartesianos em Microcomputador. Relatorio 3/90, LNEC, Lisboa (in Portuguese),1990. 6. Battjes,J. A., Long-term wave height distribution at seven stations around the British Isles. National Institute of Oceanography,Internal Report No. A.44. 1970.
in Engineering
Software 23
(1995)
49-59
6 1995 Elsevier Science Limited Printed in Great Britain. All rights reserved 0965-9978/95/$09.50
0965-9978(95)00019-4
ELSEVIER
Wave predictions based on scatter diagram data. A computer program package R. Capitio National
Laboratory
of Civil Engineering
(LNEC),
Av. do Brusil, 101, 1799 Lisbon, Portugal
& R. Burrows Department
of Civil Engineering,
University qf Liverpool,
P.O. Box 147, Liverpool
L69 3BX, UK
(Received4 November 1994;revisedversion received7 April 1995;accepted27 April 1995) This paper describesa user-friendly computer program packagecontaining an approachfor the prediction of extremewavesbasedon the extrapolation of the scatter diagram of significantwave height H, and meanzero-upcrossingperiod T,. Following the validity testsof data input, the computerprogram takes into account the existenceof missingdata and calmsand includesroutines for the treatment of various theoreticaldescriptionsof H,. Sensitivity studieshave been completed to appraisethe changesin wave predictions due to missingdata consideredeither as calms or instrument failure and due to inadequaciesof measurementsystemsto record low T, values.The possibleinfluenceof currents on extreme wave height predictions has also been taken into account. For examplepurposes,the computer program is usedfor the prediction of extreme wave heights in a location (Morecambe Bay, UK) where wave data (scatter diagrams)are available. Key words: wave prediction, extremes,scatterdiagram, oceanengineering
Skew X StDevX
NOTATION A, B, C cdf, pdf
dH dHr
dT E{Hmax} H H max
HZ 1H IHr IT
N nH nHr nT R SDError
TZ Wtj
location, scaleand shapeparametersof a cdf cumulative and probability density functions incrementof H in a scatterdiagram incrementof extrapolated H in a scatter diagram incrementof T in a scatter diagram expected maximum wave height in a 3-h record period wave height, individual wave height designmaximum wave height significantwave weight lower value of H in a scatterdiagram lower value of extrapolated H in a scatterdiagram lower value of T in a scatterdiagram number of observations, samplesize, number of adopted data number of classesof H present in the scatter diagram number of extrapolated classesof H in a scatter diagram number of classesof T present in the scatter diagram return period standarddeviation error
VarX
x, y AvX Notation
skewness of a sampleof values standarddeviation of samplevalues meanzero-upcrossingwave period number of occurrencesin the class (i,j) of the scatter diagram variance of a sampleof values orthogonal axes,co-ordinates meanvalue of x. in the Computer Program
Grafyx, Graph, Utility Gamma HofP Ord PlotPos ProbX XAxis YAxis DetectCurrents DetectMissData DistParameters 49
LTS (PROGRAMS)
external units of routines Gammafunction inverseprobability distribution function Y-variate function on probability paper units plotting position functions probability distribution function X-variate function on plotting units Y-variate function on plotting units detecting possibleexistenceof currents detectingmissingvaluesin data file calculation of distribution parametersgraphically
50
Distribution DrawLine DrawTable ExtOptions GraphicParameters Graphics InputData InterExtra LeastSquares Moments Options output 1 output2 PDFandCDF PlottingPositions Probabilities Ranks ReadFile Scales ScaleYAxis Statistics TestFitting
R. Capitcio, R. Burrows
inputting data from keyboard plotting the fitted line on probability paper plotting a table of results on probability paper inputting data from keyboard inputting data from keyboard creation of graphics on screen and plotter file inputting data from keyboard calculation of R for given H method of least squares method of moments inputting data from keyboard first writing of results on results file appending results on the created results file pdf and cdf inputting data from keyboard calculations with probabilities ranking of a sample of values reading of the data file calculations of plotting scales calculation of the scale for the probability axis sample statistics testing the data fitting
INTRODUCTION For the design of an offshore or coastal structure the engineer must consider the extreme conditions of wave climate likely to occur during the design life of the structure. Extreme wave conditions are usually represented in terms of wave heights and wave periods. These parameters are commonly obtained by analyzing wave records which contain wave characteristics measured by data acquisition systems. The wave records are usually composed of lo-20min of measurement and obtained once every 3 h. In the short-term (that is, in periods of time up to, say, 3 h) the sea state intensity can be assumed to be constant and thus each wave record is regarded to be typical of the whole 3 h-period. During the wave recording period it is assumed that the statistical properties of the sea state do not vary and hence the characteristics of a uniform stochastic process can be applied to it. The short-term wave climate is well understood. The wave heights follow the Rayleigh probability distribution if one assumes a narrow-banded spectrum and considers that the short-crested random sea elevations might be modelled as a sum of a number of linear, sinusoidal, waves travelling in different directions and having random phases. Central limits theorem shows that as the number of wave components and directions tend to infinity, the probability distribution of the sea elevation tends to the Gaussian distribution and characteristics of this distribution are well known.
For the long-term wave climate analysis each wave record of lo-20min is usually reduced to two parameters, the significant wave height, Hs, and the mean zero-upcrossing period, Tz, that are supposed to represent the whole period of 3 h. These parameters can be calculated in a number of ways being at present normally found by spectral analysis of the time series of the sea surface elevation or, sometimes, by using the Tucker-Draper method.’ In addition to Hs and Tz given by the scatter diagram, another representative parameter of the sea state conditions, the expected maximum wave height - E{H,,,} - can be computed for each 3 h observation as: E{H,,,}
z 0.707Hs(log,
n)tj2 and y1= (3 x 3600)/T,
Collection of these parameters during a long period (at least 1 yr) gives to the engineer an enlarged view of the wave climate. For engineering purposes, it is useful to use the wave climate information (Hs, Tz and E{H,,,}) in order to find the long-term probability distribution of wave heights, and thus extrapolate the observations to find the conditions of design severity. This is done by use of the theory of extremes, if the conditions of application of the theory are fulfilled, or, as a last resort, by use of an empirical procedure. Both procedures use a selected cumulative probability distribution which normally belongs to the family of extremes. This is fully described in a companion article to this one.2 Extreme wave statistics can be described by several exact probability distributions. These include the extremal Fisher-Tippet distributions (Types I, II and III) for maxima and minima and the normal and lognormal distributions. Such distributions are based on firm theoretical foundations and thus data fitting should be done in a theoretical basis as well. However, the available data often does not fulfil the conditions of a random sample so that the data fitting is often empirical. Despite this, experience has shown that extrapolations made under the use of such empirical methods are reasonably accurate. Very often, an approximate representation of the wave climate is given by a bivariate histogram of Hs and Tz. The bivariate histogram, also called scatter diagram or contingency table, represents the distribution of sea state conditions during the period when wave records have been produced. For each record, the expected maximum (the most probable) wave height, E{H,,,}, can be determined if it is assumed that wave heights are Rayleigh distributed, which is a reasonable assumption. The particular case of scatter diagram information and comments on the determination of E{H,,,} are dealt with herein. Design wave conditions are usually abstracted by using procedures for the extrapolation of the significant or maximum wave heights. In particular, one empirical procedure for extrapolation of the marginal cumulative
Wave predictions
51
based on scatter diagram data
Table 1. The files required to run LTS
LTS.EXE EGAVGA.BGI CGA.BGI IBM8514.BGI LITT.CHR file.DAT
Executable File One of the Graphics Interface Files File containing Stroked Fonts The Data File (compulsorily with extension DAT)
distribution of H, is described. However, for some applications, offshore structure design for instance, it is more important to find such conditions based on the individual wave height. This can be done approximately by extrapolating the values E{H,,,} of each record or by use of Battjes’ procedure which accounts for the complete distribution of individual waves. A more detailed review of this background is given in Capitao3 and Capita0 & Burrows.* In order to use all the procedures in a logical and integrated way, a computer program able to run on personal computers has been written using the Pascal programming language. This program is briefly described below.
THE COMPUTER
PROGRAM
LTS
The computer program LTS is briefly described here in the order that the program routines flow. The package entails a PC/XT, PC/AT, PS/2 or compatible microcomputer using IBM/MS-DOS 3.2 or above. In this version, use of a 8087, 80287 or 80387 mathematics coprocessor and a minimum of 256 KByte RAM memory are necessary. The output can be directed to a line printer and/or to a plotter, so at least one of these devices should be provided. It is recommended that this package should be installed on a hard disk since the execution time can be considerably shortened. It is important as well that some free space on disk should be allowed for the creation of the results’ files. The program runs simply by executing the command LTS. However, for the graphics’ presentation on screen, graphics interface files are needed. Depending on the graphic card and monitor present in the computer, one of the following files is required: EGAVGA.BGI, HERCULES.BGI, CGA.BGI or IBM8514.BGI. Since the program uses stroked fonts for the graphics’ presentation on screen, it is necessary to use the file LITT.CHR as well. Briefly, the files required to run the computer program LTS are given in Table 1. Since the main program calls routines that belong to units of routines, the executable files containing these units are needed for the program compilation. Units are the basis of modular programming in Turbo Pascal since they are used to create libraries and to divide large programs into logically related modules.
The required units have their executable files named as: GRAPH/TPU,4 GRAFYX.TPU’ and UTILITY.TPU.3 The LTS program can be made available, at nominal charge, from the authors. Brief description
of the computer program
The computer program LTS has been written in the Pascal programming language, using Turbo Pascal (Version 6.0) environment of Borland International.4 Basically, a Pascal computer program is constituted by four major blocks:
6) The block containing the
dejinition and declara-
of units, labels, constants, types and variables of the main program. (ii) The block containing thefunctions. A function is a program part that computes and returns a value and is activated by the evaluation of function call in an expression. (iii) The block that contains the procedures. A procedure is a program part that performs a specific action, often based on a set of parameters and is activated by a procedure statement. (iv) The block containing the main program, where the organization of the program is built. tion
The computer program is herein described in the order that functions and procedures are called throughout the program run. A flow chart of the program is given below. The main program calls the following procedures: (A) InputData
This procedure permits the input of data via keyboard. The input values are: Data File Name - Name of the file that contains the data to be analyzed by the computer program. The extension must compulsorily be .DAT; Results File Name - Name of the file that will be created with the results of the program run. By default it takes the data file name with extension .RES; Plotter File Name - Name of the file that will be created for later sending to the plotter printer. This file will be created using the plotter programming language HPGL (Hewlett Packard Graphics Language). By default, the plotter file name takes the data file name with extension .PLT;
52
R. Capitiio, R. Burrows
I
Ext options
Read file I
Q
Least squ
Plot position I Distribution I
I
Draw line
I
output 2
Fig. 1. LTS flow chart. Name of the units for both axes X and Y (variate and return period) - By default they take the names ‘metre’
and ‘year’ respectively; Number of observationsper unit of return period (year)
- It is assumed by the program that, in the case of bivariate data, the observations are 3-hourly taken and hence 2920 values are expected during 1 year. For the case of the univariate data, by default, it is assumed that there is one observation per year (annual observations); Number of units (years) of observed data - By default, the number is 1.0.
(B.2) Ranks. Ranks the sample values of the univariate sample. It is assumed that the bivariate sample (scatter diagram) is already ranked in both directions. (B.3) PlottingPositions. When the required
type of analysis is univatiate, this procedure uses any of the most common plotting positions, such as i/(N + l), i/N, etc., for the assignment of probabilities to the ranked data values. For the case of bivariate analysis with grouped data, use of only the Weibull plotting position is allowed for such a purpose.
(B) ExtOptions This procedure permits the input of the following data via keyboard (the entries in this procedure are only needed for bivariate analysis): Option between run for Hs or H; Option between whether or not extrapolation of the marginal cumulative distribution is required.
This procedure calls two other procedures: (B.1) ReadFile. This procedure reads the data file. The data file can be of two different types: univariate data (can be, for instance, annual maxima of Hs or H,,,) or bivariate data (scatter diagram of Hs and Tz, for example).
(C) Distribution In this procedure, the user can select any of the available probability distributions for the fitting of the data. The available distributions (at the moment) are: Gumbel (Fisher-Tippet Type I for maxima); Weibull (Fisher-Tippet Type III for minima); In the case of the three-parameter distribution (Fisher-Tippet Type III), the location parameter A is required if graphical extrapolations are used. If the method of moments is used for the estimation of the parameters distribution, such a requirement is not necessary since all the parameters are calculated iteratively.
Wave predictions
based on scatter diagram data
(D) Probabilities According to the required type of analysis, this procedure permits the formation of the array of wave heights and mean zero up-crossing period. This procedure calls two other procedures: (D.1) Statistics. This procedure calculates some important sample statistics like the mean sample value (AwX), the Variance ( VarX), the Standard Deviation (StDewX) and the Skewness (SkewX). This procedure calls the procedure: (D.1.1) Moments. For the estimation of the parameters distribution, the method of moments is herein applied using a direct calculation for the case of FT Type I distribution and an iterative process for the case of the FT Type III distributions.
The pdf and cdf (raw and extrapolated marginal probabilities) of both significant and individual wave heights and the expected value E{H,,,} are calculated here. (0.2)
PDFandCDF.
(E) Output1 This procedure opens the results file and writes on it the statistics calculated in procedure (D.l) and the pdf/cdf calculated in (D.2). (F) GraphicParam More data is read from the keyboard if graphical procedures are required. The following entries are required for the plotting of the graphics in both screen and plotter printer: Minimum and maximum abscissa for graphics - Default values are provided: Interval between X-axis graphic units - By default is 1; Number of ticks between X-axis graphic units - By default is 10; Required number of return periods - The program calculates the return values of H or Hs corresponding to each required return period. By default, the number of required return periods is 0. (G) TestFitting To decide which distribution fits best to the sample data, calculation of the correlation factor r and standard deviation error, SDError, is done herein. These factors can only serve as an approximate indicator of the goodness of fit of the selected distribution. Depending on the selected options, this procedure might call three other procedures. (G.1) Scales. Some distributions need to have a nonlinear variate axis. Thus, a transformation of the xcoordinate axis is required for the plotting on the screen
53
or plotter printer. According to the selected distribution, this procedure scales up that axis. (G.2) LeastSquares. This procedure is called when the method of least squares is required. In order to find the best fit of a straight line through the points [H, P(H)], the standard least squares method is herein applied and thus the parameters of the straight line are found. (G.3) DistParameters. If the least squares method is required, this procedure, based on the parameters of the fitted straight line, calculates the parameters of the selected distribution.
(H) Graphics The graphic routine is designed to present, through the screen and plotter, a graph showing the straight line fitted through the data points, according to the selected distribution. This routine tests the monitor type and, by default, creates a pre-defined graphic layout. At the present, the routine permits to scale up the probability axis from 0.000001 to 0.99999999999 but this range can be extended (reduced) if required. This procedure calls three other procedures. (H.l) Scale YAxis. Graphic routine destined to scale up the probability axis using both pixel units (for screen output) or plotter units (for plotter output). (H.2) DrawLine. Graphic procedure that permits the drawing of the fitted line through the data points. (H.3) Draw Table. Graphic routine that permits the creation of a table containing the predictions of wave height, based on the fitted line, for a given return period.
(I) output2 This procedure re-opens the results file and appends on it the results of the estimation of distribution parameters, correlation coefficient, standard deviation error (graphic procedures) and return values for the required return periods. (J) Options This procedure contains a menu of options that permits the user to return to other menus, run the program again with the same data or different data, changing of graphic parameters or plotting positions (only in the case of univariate analysis). This procedure calls another procedure. (J.1) InterExtra. If required, calculation of the return period R for any given value of Hs is completed.
54
R. Capitdo, R. Burrows
Input to computer program LTS
done for both significant wave heights and individual wave heights (by using the Battjes’ method6). Samples of the input data report and tabular results arising are given in Figs 2 and 3. Plots of graphical results for both significant wave height and individual wave height, each in Gumbel and Weibull formats, are shown in Figs 4-7. An important observation in relation to the individual waveheight (H) data is that the actual points plotted represent the true distribution.2 The fitted lines, by least squares here, serve as approximations only and, in the case of the Gumbel fits, can depart markedly from the true values in the upper tail. The Gumbel distribution, therefore, gives a poor representation of the distribution of individual wave heights, whilst the Weibull distribution, with the flexibility of an extra parameter, provides acceptable fits. The sample data has been studied for the purpose of determining how missing data, considered either as calms or non-measured data, can influence the wave predictions. Furthermore, an assessment has been made to establish, for this location, whether inability of the measurement systems in recording low Tz values can affect the same wave predictions. It was observed from analysis of the results that predictions of significant wave height found by using the method of moments were slightly greater than those obtained by using the least squares method as can be seen from Table 2. The systems that exist for the measurement of wave data have limitations and the accuracy of such a measurement can vary substantially from system to system. It is important thus to be aware of the effect on wave predictions when a particular measurement system is unable to measure a certain wave characteristic. One common problem is the difficulty of measuring low values of Tz. Low values of Tz are usually correlated with low values of Hs. Consequently, on average, if one does not consider low values of T, that means the predictions of wave height will be higher than those calculated if all Tz values are used. In order to illustrate the differences in predictions, three data samples from Morecambe Bay have been considered. The first sample (all data with Tz readings beginning at 1.75s) is presented in Table 2. The second and third data samples have been created upon the first one but assuming that the measurement system was unable to record waves having Tz values less than 2.75
The computer program only needs input data from the keyboard (or default values) and a file containing the data sample. The layout of such a file depends on the type of analysis required. For an univariate sample (for instance, annual maxima), the data file is simply constituted by a column which contains the data values x1, x2, . , x,,: ( beginning of data file ) First extreme value Second extreme value Third extreme value
Xl x2 x3
( end of data file ) The data file layout for a bivariate diagram) is shown below.
nib extreme value
sample (scatter
Output of the computer program Three output files are created computer program LTS:
when running
the
(i) one file that is dumped to the screen; (ii) one text file (ASCII format) containing the results ready to be printed on a line printer (by default, with extension RES); (iii) one plotter file that can later be sent to the plotter printer (by default, with extension PLT). This file is written in HPGL programming language, so any plotter accepting the commands HPGL is suitable. Generally, most commercial plotters accept this language. APPLICATION
TO WAVE DATA
All scatter diagram data Scatter diagram data collected at Morecambe Bay in the British Isles have been considered in this study. Wave heights predictions associated with return periods of 10, 50 and 1OOyr have been calculated by the computer program LTS on both Gumbel and Weibull probability papers, fitting by using both the least squares method and method of moments. These calculations have been
( beginning of data file ) nH nT nHr w(l,l) w(2,l) ..
w(nH,l)
1H 1T 1Hr w(l,2) w(2,2)
dH dT dHr . .. .. . .. ..
41, nT) w(2,nT) I . .
w(nH,2)
. .
Number, lower value (class average) and increment of H (or Hs) Number, lower value (class average) and increment of Tz Number, lower value and increment of extrapolated H (or Hs) Occurrences of the Is’ class of H Occurrences of the Znd class of H .. w(nH,nT)
( end of data file )
Occurrences of the nth class of H
Wave predictions
based on scatter diagram data
tt*******************************************.* * LONG-TERM STATISTICS
55
+
t******t**************t*************************
******t***********t****t******tttt*tt****~~* GENERAL DATA : **t********i********t*tt********~***~*,**~********~ TITLE
: Morecambe
Descyription
: All
Data
file
Bay
Data
: MCBl.DAT
Results
file
: MCBl.RES
Platter
file
: MCBl??.PLT
Wavr
Height
unities
( H )
Wave
Period
unities
I TZ
Return
Deriod
Number
Iof
observations
L!nities
Numl~er
of
years
Number
of
.._...._
second
( T
1 .__..._
year
cause
3.00
. .._....
data
Minlmum
value
of
sample
.._...._.._
Maximum
value
of
sample
.
Value
7961
. . .._...._
Missing
Mean
year
..____..___.
per as
metry
)
observed
ohs.
__..__.._
of
2920.0
.
non-was.
799
.
0.50 6.00
.._...._
.._____._.___..____.____.
1.0724
Standard
Deviation
Variance
. . ..__.....__...____..__.._
0.5R73
Skewness
___._..__...____...___..___
1.3430
SCATTER
DIAGRAM
_._._____......._
0.7663
DATA
H'S 0.15 0.75 1.25 1.75 2.15 2.75 3.25 3.75 4.25 4.75 5.25 5.75 6.25 6.75 7.25
2.0 120 73 48 0 0 0 0 0 0 0 0 i 0 0
2.5 596 316 58 0 0 0 0 0 0 0 0 0 0 0 0
3.0 537 851 69 0 0 0 0 0 0 0 0 0 0 0 0
3.5 289 640 313 9 0 0 0 0 0 0 0 0 0 0 0
4.0 113 312 499 86 4 0 0 0 0 0 0 0 0 0 0
4.5 45 172 327 283 52 2 1 0 0 0 0 0 0 0 0
5.0 33 103 192 270 117 35 2 0 0 0 0 i 0 0
PERIOD 5.5 34 89 86 130 160 84 16 : 0 0 0 0 0 0
6.0 17 49 57 68 58 53 46 8 18 0 0 0 0 0 0
6.5 21 30 29 33 22 39 45 22 113 110 0 0 0 0
7.0
7.5
8.0
9 9 3 12 15 11 10 10 3
6 3 0 3 3 11 8 4 6
0 0 0 0
0 0 0 0
8.5
1 3 1 0 10 3 3 10 1 3 0 0" 0 0
9.0 R 3 0 1 0 2 3 0' 0 0 0 0 0
9 0 0 0 0 0 0 0 0 0 0 ; 0 0
Fig. 2. Sample input data report.
and 3.25 s, respectively. Wave predictions based on these two samples are presented on Tables 3 and 4. All the predictions were found by considering the missing data as non-measured data. For the Gumbel distribution it was found that when low values of T, are not taken into account, the wave predictions become slightly higher than those obtained when using all Tz data. However, if the data fitting is made using the Weibull distribution, the wave predictions are lower when low values of Tz are not
considered. Thus when acquisition systems are unable to read low values of T,, a slight overestimation or underestimation in wave predictions will occur if one fits the data points to the Gumbel or Weibull distribution, respectively. The effects are, however, very small, as seen from the comparison of Tables 2-4. It is now necessary to consider analysis of a scatter diagram with missing data to appraise the consequences in the context of wave predictions. Missing data can be considered as the result of calms (not measured by the
56
R. Capitdo, R. Burrows TABLE OF RESULTS ----------------------------------------------------------~ HsI Hsi Wi SWi
PDFHs
CDFHs
PDFHsX
CDFHsX
-!m!___!m!--!"o:_"ccUrl____lI1______II)_-----!~!------!~!-0.00 0.25
lR38
1838
0.75
2653
4491
1.25
1682
6173
1.75
895
7068
2.25
432
7500
2.75
238
7738
3.25
133
7871
0.50 1.00 1.50 2.00 2.50 3.00 3.50 3.75
58
7929
4.25
21
7950
4.00 4.50 4.15
8
7958
5.25
2
7960
5.75
1
7961
5.00 5.50 6.00
GUMBEL DISTRIBlJTION *ttf****f*t****ttt**ttt**~************************ SIGNIFICANT WAVE HEIGHT P(x)
= exp
I -exp
( -(
0.23088
0.23085
0.23142
0.23142
0.33325
0.56405
0.29913
0.53055
0.21128
0.77531
0.22939
0.75994
0.11242
0.88772
0.12796
0.88791
0.05426
0.94197
0.06191
0.94981
0.02990
0.97187
0.02813
0.97795
0.01671
0.98857
0.01244
0.99039
0.00729
0.99586
0.00544
0.99583
0.00264
0.99849
0.00236
0.99819
0.00100
0.99950
0.00103
0.99922
0.00025
0.99975
0.00044
0.99966
0.00013
0.99987
0.00019
0.99985
WEIBULL DISTRIBUTION *t******t*tt********t*****+tttt*t*t******************** SIGNIFICANT WAVE HEIGHT X - A J/B
)
]
P(x)
= 1 - exp
----------LEAST
SQUARES
41s Bls
= =
METHOD
0.77161 0.57908
I
METHOD
I I
Amm= Bmm=
OF
MOMENTS
r
= 0.99843 I SDError = 0.16721 I --------_------~-_~-____________________-----~~ --------
Extrapolation
of
SQUARES
Als Bls ClS
0.59751
= = =
T 10 50
100
.
Hs(LS1
.
6.73 7.66 8.06
Hs(MM)
. Hmax(LS1.
6.87 7.83 8.25
12.64 14.33 15.05
[ 12.00 i 23 1
--------
of
Extrapolation
12.90 14.64 15.39
. T . H . ____-----___------__--------------------------10 15.95 17.55 18.25
0.08321
I
Amm=
0.08321
1.09011 1.32633
I I
Bmm= cm=
1.07140 1.30186
0.99940
I
0.03518
I
Extrapolation
of
OF
Observations
.
T 10 50
.
Hs(LS)
.
6.40 7.13 7.44
Hs(MM)
. Hmax(LS).
6.50 7.26 7.58
12.05 13.38 13.93
WEIBULL DISTRIBUTION *t***tt*****t*t***tt****************~***~***,** INDIVIDUAL WAVE HEIGHT Limit of Extrapolation [ 12.00 Number of Classes I 23 1
1
Observations
METHOD
MOMENTS
--------
Hmax(MM).
100
INDIVIDUAL WAVE HEIGHT Limit of Extrapolation Number of Classes
----------_
parameters
METHOD
1
---------------
.
/ B 1°C
I
r I =
SDError Observations
I 1X - A)
Distribution
LEAST
0.72756
[ -
--------
--------
Extrapolation
-____-----__----------------------------------. T . H . 10
14.41 50 15.90 100 16.52 _----_____------___---------~~~~~~~~~--------~~
Fig. 3. Sampletabular resultsreport.
of
Hmax(MM). 12.23 13.61 14.18
]
Observations
--------
Wave predictions
based on scatter diagram data
MB -Miss. data: CALMS Gumbel distribution Significant wave heights
99.9990
t
+
+
Least sq-s -2.329 0.996
/ A= B=
/
Wave
height
(m)
Wave
height(m)
Fig. 4. Gumbel plot of significant wave heights, H,.
Fig. 6. Gumbel plot of individual wave heights, H.
acquisition system because of the very low values of Hs) or as non-measured data resulting from randomly occurring malfunction of the system. In this last case, it is important to note that the predictions can be seriously affected if occurrence of malfunction is biased to occur when the sea state is rough. In this case the data set would be unacceptable for accurate wave predictions. In order to compare wave predictions when considering missing data either as calms or non-measurements due to malfunction of the acquisition system, Tables 5 and 6 have been produced for Morecambe Bay.
From the comparison of Tables 5 and 6 it can be concluded that the decision to consider missing data as non-measured data instead of as calms produces slightly higher predictions of wave heights as might be anticipated. In other localities the occurrence of ‘calms’ will be more prevalent and the consequences of their misrepresentation will be much more significant. CONCLUSIONS Scatter diagram data have been used for predictions of wave heights by use of the computer program LTS.
MB - Miss. data: CALMS Weibull distribution Individual wave heights
Weibull distribution Significant wave heights
30.0000
t
HI101 = H[501= H(lOOl=
Least squares
Moments
0.083 1.090 I.326
0.083 I.071 1.302
6.36 7.09
6.45 7.21 7.53
7.40
2o.ml
22, .e 4 % &
Leas1 sq”areS
9o.ooooooO-
8O.OOOOOOO 7o.OOoOOO0,
Wave
height
(m)
Fig. 5. Weibull plot of significant wave heights, H,.
Wave
A= LI= c=
0.083 0.542 0.888
H [lo] =
14.39
H 150) = H[lOO]=
15.82 16.44
I 6
I II
height
i C
I J16 21
(m)
Fig. 7. Weibull plot of individual wave heights, H.
58
R. Capitfio, R. Burrows Table 2. Wave height predictions -
Morecambe Bay (all data)
R
f&m
W
Gumbel
10 50 100
Weibull
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 7.7 8.1
6.9 7.8 8.3
6.4 7.1 7.4
6.5 7.3 7.6
12.6 14.3 15.1
12.9 14.6 15.4
12.1 13.4 13.9
12.2 13.6 14.2
Table 3. Wave height predictions &I
z
Gumbel
10 50 100
- 3h (ml
Weibull
Weibull Battjes6
16.0 17.6 18.3
14.5 15.9 16.5
Morecambe Bay ( Tz > 2.75 s) ffmax - 3h W
Gumbel
Weibull
6s
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.8 7.7 8.1
7.1 8.0 8.5
6.4 7.1 7.3
6.5 7.2 7.4
12.8 14.5 15.2
13.2 15.0 15.8
12.0 13.2 13.7
12.1 13.4 14.0
Weibull
Battjes6 16.2 17.9 18.6
14.4 15.8 16.4
The number of Tz values >2.75s is 6750, compared with 7961 (all data). Table 4. Wave height predictions -
Morecambe Bay (T’ > 3.25s) fLx
&I
Gumbel
10 50 100
- 3h (ml
z
Weibull
(:I
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
7.0 7.9 8.3
7.3 8.3 8.8
6.4 7.0 7.2
6.3 6.9 7.2
13.0 14.7 15.4
13.7 15.5 16.3
12.0 13.1 13.6
11.8 13.0 13.4
Weibull
Battje@ 16.7 18.3 19.0
14.1 15.4 15.9
The number of T, values>3.25 s is 5293,comparedwith 7961(all data). LSM - Distribution parametersestimationby Least SquaresMethod. M - Distribution parametersestimationMethod of Moments. Battjes’ ( 1970)6procedurefor the prediction of individual wave heights.3
Several conclusions can be drawn from comparisons between the results obtained:
(1) Acquisition
systems that are inadequate to register low values of Tz are likely to overestimate predicted values of wave heights if the Gumbel distribution is used, whereas underestimation may occur if fitting is made using the Weibull distribution, according to the sample application made herein. (2) With missing data in the scatter diagrams, the decision to consider this either as ‘calms’ or ‘non-measurements’ leads to different predicted values of wave height. For design purposes it is important to know what has really happened with those missing data. When in doubt, it is better to consider the missing data as non-measurements since it is a conservative decision. (31 Using the ‘weighted least squares method’, a greater importance is given to the central range of points than using the method of moments and wave height predictions are lower. It is not possible to state categorically that one probability distribution is better, that is, represents better the data set, than another. The Weibull
distribution is here shown to be more flexible in its fit to a particular set of data since it contains three parameters whereas the Gumbel distribution has only two. However, in practice, the Gumbel distribution has generally been accepted as the principal extreme value distribution for the current application. It is, in fact, a distribution that generally applies very well to wave data and the estimation of parameters is straightforward. In this study, long-term variability of the characteristics of the wave climate has not been a consideration. This represents an important simplification of the problem. Appreciation must be given to the fact that scatter diagram data records covering perhaps only 1 or 2yr may not be ‘typical’ of the true long-term wave climate and some intuitive adjustment may be necessary if, for example, wind records show that the wave climate sampled may have been more, or less, severe than the average. An improvement in the quality of analysis of extreme wave data is achieved if emphasis is placed on the effort of collecting as many wave data as possible and over a period as long as is feasible. This enables the engineer to have a better perception of the true wave climate that exists at the site and the danger of selecting a
Wave predictions Table 5. Wave height predictions -
Morecambe Bay (missing data as non-measured data) f&m
- 3h Cm)
i2
Weibull
Gumbel 10 50 100
59
based on scatter diagram data
Weibull
Gumbel
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 1.7 8.1
6.9 7.8 8.3
6.4 7.1 7.4
6.5 I.3 7.6
12.6 14.3 15.1
12.9 14.6 15.4
12.1 13.4 13.9
12.2 13.6 14.2
Weibull Battjes6
16.0 17.6 18.3
14.5 15.9 16.5
Calms are 799 out of 8760. Table 6. Wave height predictions CRY4
Morecambe Bay (missing data as calms)
2
ffmax
- 3h
(ml
Gumbel 10 50 100
Weibull
(:I
Gumbel
Weibull
Gumbel
LSM
M
LSM
M
LSM
M
LSM
M
6.7 7.6 8.0
6.8 7.8 8.2
6.4 7.1 I.4
6.5 7.2 I.5
12.5 14.2 15.0
12.8 14.5 15.3
12.0 13.3 13.9
12.2 13.5 14.1
Weibull Battjes6
15.9 17.5 18.2
14.4 15.8 16.4
Non-measurements are 799 out of 8760. LSM - Distribution parameters estimation by Least Squares Method. M - Distribution parameters estimation by Method of Moments. Battjes’ (1970)6 procedure for the prediction of individual wave heights.’
distribution that does not belong to the population will certainly be minimized. The computer program must be used carefully. It is important that the user is aware of the whole process of predicting wave characteristics. In particular, the user must know what is the reliability of the data sample that is going to be used since enormous mistakes can arise from data samples not containing a reasonable set of extreme values.
REFERENCES 1. Draper, L., Derivation of a ‘Design Wave’ from instrumental records of sea waves. Proc. Z.C.E., 1963,26, 291l 304.
2. Capitao, R. & Burrows, R., Wave predictions basedon scatter diagram data: A review. Adv. Engng Sqftwure, 1995,23, 37-47. Capitao, R., Wave predictions basedon scatter diagram data. MSc(Eng.) Thesis, Dep. of Civil Engng, Univ. of Liverpool, 1992. Borland International, Turbo Pascal -- Version 6.0. Library Reference and Programmers’ Guide. Borland International, Scotts Valley, CA. Palma, J., Bibliotecas em Pascal para Constru9Bo de Graficos Cartesianos em Microcomputador. Relatorio 3/90, LNEC, Lisboa (in Portuguese),1990. 6. Battjes,J. A., Long-term wave height distribution at seven stations around the British Isles. National Institute of Oceanography,Internal Report No. A.44. 1970.