- Email: [email protected]
Investigation of model and parameter uncertainty in water quality models using a random walk method
Journal of Marine Systems 28 Ž2001. 269–279 www.elsevier.nlrlocaterjmarsys
Investigation of model and parameter uncertainty in water quality models using a random walk method A.M. Riddle ) Brixham EnÕironmental Laboratory, AstraZeneca UK Limited, Freshwater Quarry, Brixham, DeÕon, TQ5 8BA, UK Received 3 September 2000; accepted 12 February 2001
Abstract A mathematical model to predict the effect of chemical spills in the Forth estuary in Scotland has been in use for many years. The model, based on the random walk method, predicts chemical concentrations in the estuary waters and estimates the elapsed time before the dilution is sufficient to render the spill harmless Žmaking use of a toxicity measure such as the LC50 or a water quality standard.. The model gives a deterministic result without any estimate of the uncertainty. Field studies using tracer dyes to measure the horizontal and vertical mixing rates in the estuary show that these rates vary over time. The literature on turbulent diffusion includes modelling applications using different parameterisations of the mixing process. This paper investigates the uncertainties in predicted concentrations due to model parameterisation of horizontal mixing and due to the variability in the measured mixing rates determined from surveys in the estuary. Estimates of the range of concentrations for a specific spill scenario are presented. The study shows that model formulation and parameter uncertainty are both important factors in estimating the uncertainty in model predictions. The uncertainty caused by the variations with time found in the measured mixing rates is found to be of similar magnitude to the differences in concentration resulting from using three different methods for modelling the horizontal mixing in the estuary. Uncertainties associated with model formulation could be reduced if a small number of longer timescale Že.g. 24 h. dispersion experiments were available. In addition, further data from short-term Ž; 3 h. dispersion experiments would give a better understanding of the distribution of mixing coefficients and how the mixing relates to other parameters such as tidal range and wind speed and direction. q 2001 Elsevier Science B.V. All rights reserved. Keywords: Uncertainty; Random walk; Hydrodynamic; Model; Mixing
1. Introduction A model for predicting the effect of a chemical spill in the Forth estuary in Scotland was set up in 1990 ŽRiddle, 1991.. The model covered the lower estuary from the Forth bridges in the east to Alloa in )
Tel.: q44-1803-882-882; fax: q44-1803-882-974. E-mail address: [email protected] ŽA.M. Riddle..
the west ŽFig. 1., a distance of approximately 25 km. In this area, the estuary contracts rapidly from a width of approximately 3.5 km and depth of 20 m in the east to a width of 0.25 km and depth of 2.5 m at Alloa. The estuary is macro-tidal with a mean neap and spring tidal ranges of 2.5 and 5 m, respectively, and peak tidal currents of 1.4 m sy1 in the main channel off Grangemouth. The estuary is generally well mixed, however, a small degree of salinity stratification Žsurface to bottom salinity difference
0924-7963r01r$ - see front matter q 2001 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 4 - 7 9 6 3 Ž 0 1 . 0 0 0 2 7 - 6
270
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Fig. 1. Lower Forth estuary, Scotland.
; 4 kg my3 . forms during the flood tide on neap tidal ranges. A more detailed description of the bathymetry, hydrodynamics and mixing characteristics of the estuary is given by Webb and Metcalf Ž1987.. The estuary model was based on diffusion of a Gaussian patch method to simulate the concentrations in the estuary, using interpolated field measurements of tidal current on a 200 = 200-m grid to simulate the transport processes. This original model has been enhanced to include a hydrodynamic simulation of the tidal currents using a 50 = 50-m grid and a three-dimensional random walk model to determine the dispersion ŽRiddle,
1998.. The hydrodynamic model is based on the finite difference solution of the depth averaged Navier–Stokes equations ŽLeendertse, 1967., with modifications to allow for flooding and drying of intertidal mudflats ŽFalconer and Chen, 1991.. The hydrodynamic model has been calibrated against a range of field data from the estuary including data from recording current meters, drogue tracking and dye movement from a range of tracer experiments ŽRiddle, 1998.. The water quality model uses the hydrodynamic data, stored at half hourly intervals throughout a tidal cycle together with a vertical current profile in a
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
three-dimensional particle tracking random walk approach ŽWebb, 1982; Van Dam and Louwersheimer, 1992; Van Dam and Geurtz, 1994. to simulate the chemical concentrations in the estuary ŽRiddle, 1998.. The calibration of dispersion in the water quality model was based on mixing coefficients determined from dye tracer experiments; the information consists of measurements of the dye patch size Žhorizontal and vertical. and concentration, taken over a relatively short timescale Ž6 h or less.. The model has subsequently been used for simulating effects that last for several tidal cycles. The horizontal scale of the model is small enough to resolve the current field and shear to produce chaotic motions such as described by Ridderinkof and Zimmerman Ž1992. for the Wadden Sea. The aim of this paper is to investigate the uncertainties in predicting concentrations in the estuary over a 48-h period resulting from a chemical discharge. Two aspects of uncertainty in predicted concentrations are considered in the paper:
271
obtained, which is comparable with the values obtained by Van Veen Ž1938. and employed by Van Dam and Louwersheimer Ž1992.. The use of a full three-dimensional hydrodynamic model is planned for future work, which should provide an improved simulation, particularly allowing for changes in flow direction through the water column in areas of rapidly changing water depth. The random walk model uses the appropriate tidal data from the hydrodynamic model at each timestep and represents the effluent discharge by placing a fixed number of AparticlesB at the outfall position on each timestep. These particles are moved on each subsequent timestep according to the following method ŽWebb, 1982.:
ž / .ž /
X new s Xold q U Ž p q 1 . Ynew s Yold q V Ž p q 1
p
Zold d
Zold d
dtqfx a
p
dtqfy b
Znew s Zold q f z g 1. The differences resulting from applying three different theoretical formulations, used by other researchers to model horizontal dispersion. 2. The differences arising from the uncertainty in the horizontal and vertical mixing rates and tidal current speeds, as shown by field measurements. The model theory is presented in Section 2, followed by a description of the uncertainty analysis of model formulation and parameter variability in Section 3. The results of the modelling study are discussed at the end of Section 3, and conclusions are presented in Section 4 of the paper.
2. Model theory
where X, Y and Z represent the position of a particle and subscripts ‘old’ and ‘new’ represent the positions at the start and end of a model timestep; d t is the timestep length, U and V are the depth averaged horizontal velocity components from the hydrodynamic model and d the local water depth. The functions f x , f y and f z define the mixing process and a , b and g are random numbers from a standard normal distribution. The vertical velocity is assumed to be zero. Each particle represents a fixed mass of effluent, and it is assumed that no degradation takes place. A constant diffusion coefficient, Fickian, model Žtheory K h ., has been used to parameterise the horizontal and vertical diffusion ŽBowden, 1983; Webb, 1982.:
(
f x s f y s 2 K hd t The three-dimensional random walk model uses tidal flow data from a two-dimensional depth averaged hydrodynamic model together with a prescribed vertical current profile in the direction of the flow ŽVan Dam and Louwersheimer, 1992; Van Dam and Geurtz, 1994.. The vertical current profile is expressed as a power law function and calibrated with local field data; an exponent Ž p . of 0.192 has been
(
fz s 2 Kzd t
Ž 1. Ž 2.
where K h and K z are the horizontal and vertical mixing coefficients Žm2 sy1 ., obtained from field studies. Coefficient K h represents the effect of turbulent mixing, since the shear is modelled directly by the spatial distribution in the hydrodynamic data and by the vertical current profile.
272
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
The concentration distribution of material in the estuary is estimated using a counting grid Žnot necessarily the same as the grid used by the hydrodynamic model.. The number of particles in a grid square
over a depth interval from the water surface down to a specified depth Žin this case a 1-m layer was used. is counted, giving the mass of material in a known volume, and therefore the concentration. An example model run has been carried out for a hypothetical spill of 25 kg over 0.5 h of a material with a safe concentration of 0.005 mg ly1 . Fig. 2a–c shows the movement and dilution of the material at 6, 12 and 24 h from discharge, whilst Fig. 2d shows the 95 percentile concentration of the patch plotted against time. The concentration is predicted to become safe after 20 h.
3. Uncertainty in model predictions 3.1. Mixing theory
Fig. 2. Predictions from the model for a discharge initially at position =; Ža. 6 h, Žb. 12 h and Žc. 24 h after discharge. Žd. Shows the time change in the predicted concentration Žthe dotted line represents the safe concentration..
The literature on modelling contaminant dispersion in estuaries and coastal areas shows a number of methods for parameterising the horizontal mixing. The application of these different methods can lead to a range of predicted concentrations and therefore to a range in the predicted time for concentrations to reach a safe level. The aim of this section is to quantify the differences resulting from using three different approaches to modelling the horizontal dispersion of a chemical spill in the Forth estuary, and to determine whether the formulations lead to a practical difference in the predictions. Three different theoretical formulations have been considered for the parameterisation of the horizontal mixing in the model; these are the following. Ø Fickian diffusion parameterised by a constant mixing coefficient as shown by Eq. Ž1. above; this is referred to as theory K h . Ø A mixing velocity formulation originally described by Joseph and Sender Ž1958. and used in an equivalent form by Elliott and Wallace Ž1989. and Rye et al., Ž1998.. This is parameterised as a constant diffusion velocity, Vm in Eq. Ž3. Žreferred to as theory V.. This method increases the horizontal mixing with time and simulates the effect of larger eddies contributing to the patch spread as the size of the patch increases. f x s f y s Vm 'td t
Ž 3.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Ø A method relating the mixing rate to the tidal current speed is based on the local application of the streamflow dispersion equation of Elder Ž1959. and is parameterised by a dispersion constant, k, in Eq. Ž4. Žsee Fischer et al., 1979, p. 125.. An example of the application of this formulation is given by Falconer Ž1986. in a two-dimensional depth averaged model for pollutant dispersion in Poole harbour England Žreferred to as theory U.: f x s f y s '2 kUHd t
Ž 4.
where U is the tidal velocity and H is the local water depth. This formulation specifies a spatially varying diffusion which can lead to the accumulation of particles in low diffusion areas; to compensate for this, a drift term dependent on the variation in the diffusion coefficient ŽHunter et al. 1993. has been incorporated into the model. The vertical mixing is parameterised by a constant mixing coefficient, K z , which simulates the vertical spread of material using Fickian diffusion ŽEq. Ž2.. as described by Lewis Ž1997.. Measurements of the vertical spread of dye patches in many areas around the UK ŽRiddle and Lewis, 2000. have shown applicability of this approach. In addition, the data for the Forth Žsee Section 3.2. indicate that mixing over depth will be complete in less than 4 h, which is a much shorter period than used for the model simulations. To calibrate each model results from a dye dispersion experiment, carried out in the Forth estuary on 13 July 1988 over a period of 6 h on the flood and ebb tides, were used. The data collected during the experiment consisted of horizontal patch size Žlength and width., patch depth and dye concentration. Beyond a time of 4 h the patch size data became unreliable due to difficulties in detecting the edges of the patch as it became dilute. The patch limits were determined using a submersible fluorometer with a detection limit of 0.001 mg ly1 . This has an influence of the fraction of the dye patch determined by the measurement, but in all cases the peak concentration was greater than a factor of 100 times the detection limit indicating that at least 99.75% of the patch had been measured Žassuming a Gaussian distribution of concentration in the horizontal direction..
273
Concentration data was only collected for the first 1.3 h of the experiment. The random walk model was used to simulate the dye dispersion in three different ways using the theoretical parameterisations for horizontal mixing represented by Eqs. Ž1., Ž3. and Ž4.. Eq. Ž2. was used to describe the vertical diffusion in each instance. All of the methods can adequately represent the measured dye patch spread and dilution. Fig. 3 shows the comparison of the methods against the measured data for dye patch size, expressed as standard deviation Ždefined as the square root of the patch length times the patch width. and the dye concentration. Each model gives a good prediction of the patch size over the 4-h timescale and the concentration over the shorter timescale of 1.3 h. For a model simulation that exceeds the length of the calibration period it might be expected that the three methods would give different concentration predictions. This has been investigated by running the calibrated models, using each of the mixing formulations, to simulate the rate of dilution of a chemical spill in the estuary. The conditions outlined in Section 2 have been used to predict concentrations over a 48-h period. Fig. 4 shows the three different predictions of concentration, which vary by up to a factor of 10 times Žat time 20 h after the spill.. By the end of four tides, the predictions have converged, because the effluent has mixed fully across the estuary cross-section area. The time required for the concentration to reach the safe concentration of 0.005 mg ly1 is 10, 22 and 36 h ŽFig. 4. for the diffusion velocity Ž V ., constant diffusion coefficient Ž K h . and tidal current ŽU . methods, respectively. Thus, the uncertainty in the duration of the incident is a factor of 3.6 times.
3.2. Mixing parameter Õalues In addition to uncertainty due to the theoretical modelling approach, the parameters used in simulation may be variable; for example, the tidal currents and mixing coefficients for the Forth estuary vary, as shown by the results of repeated studies in the estuary. This can lead to a significant variability in the predicted transport and dilution in the estuary
274
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Fig. 3. Comparison of predicted and measured values for Ža. standard deviation of the horizontal patch size and Žb. peak dye concentration Žmg ly1 ..
from day to day. Data from seven different dye experiments ŽTable 1. show that the horizontal dis-
persion coefficient Ž K h . varies by a factor of 17 and vertical mixing coefficient Ž K z . by a factor of 16.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
275
Fig. 4. Predicted concentration Žmg ly1 . in the estuary from the three mixing models.
These data for horizontal and vertical mixing are shown in Fig. 5 as cumulative frequency plots against logŽ K h . and logŽ K z . ŽFig. 5a,b.. The quantity of data is small but the plots indicate a linear relationship between cumulative frequency and the log of the mixing coefficient. A positive lower bound would be expected for both parameters and minimum values have been estimated using linear regression fits to the data giving K hmin s 0.024 m2 sy1 and K zmin s 0.0002 m2 sy1 . These cumulative frequency plots correspond to an inverse relationship between the frequency of occurrence and the mixing coefficient. The random walk model has been used to simulate the effect of this variability on the spill dilution
and for this simulation the mixing has been simulated using the constant diffusion coefficient Ž K h . model. In the simulation, each model particle has been tagged with a horizontal and vertical mixing coefficient taken from the inverse distributions defined above. This is achieved as follows: LK h s LK hmin q a Ž LK hmax y LK hmin . LK z s LK zmin q b Ž LK zmax y LK zmin . where LK represents the log value of the mixing coefficient and a and b are random numbers from a uniform distribution between 0 and 1. For these
Table 1 Dispersion coefficients in the Forth estuary—from dye dispersion experiments Date
Time ŽBST.
Tide
Duration Žh.
Horizontal dispersion coefficient Ž K h m2 sy1 .
Vertical dispersion coefficient Ž K z m2 sy1 .
7r7r1986 7r7r1986 8r7r1986 9r7r1986 21r7r1986 13r7r1988 14r7r1988
0504 1127 0600 1212 1215 1120 0500
Ebb Flood Ebb Flood High water slack FloodrEbb Ebb
0.7 1.2 1.2 1.4 1.1 5.5 3.0
– – 0.18 0.05 0.16 0.75 0.85
0.0029 0.0008 0.0014 0.0009 – 0.0047 0.0003
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
276
Fig. 5. Cumulative frequency distributions Ž%. of the measured horizontal and vertical mixing values plotted against the logarithm of the mixing coefficient.
calculations, the regression equations have been used to compute values of 1.17 m2 sy1 for K hmax and 0.0047 m2 sy1 for K zmax . Repeated measurement of the tidal current at five positions in the estuary over 2 days on spring tides Ž25r26 June 1979. and 2 days on neap tides Ž2r3 July 1979. showed a variability of approximately 10% in the peak current amplitude in the model domain. This uncertainty in the tidal currents is modelled using the factor, C, to multiply the tidal currents: C s Cm q g Cs
Cm s 1
Cs s 0.05
where g is a random number from a standard normal distribution. Each particle has a unique values for the factor C as well as K h and K z throughout the simulation. At selected times, the information for
each particle was written to a file together with the particle position information. These data are then analysed using a postprocessing program to give concentration information for different ranges of the input parameters. Fig. 6 illustrates this process. The model output data, for a given time after the spill, can be represented in a three-dimensional space with axes K h , K z and U. This space can be divided up into elements; the particles within each element represent the chemical concentration and location for the given range Žsmall. of the parameter values. In Fig. 6, the right-hand element represents higher mixing conditions and gives rise to a larger but less concentrated patch of chemical than the left-hand element which represents lower mixing conditions. Taking the 95 percentile concentration from each element in the parameter space yields the overall concentration distribution from the uncertainty analysis. Fig. 7 shows the maximum and minimum predicted concentrations over the first 30 h of the simulation, together with the predicted concentration values obtained from a deterministic model run using the standard conservative values for the parameters. The largest range between the maximum and minimum predicted concentrations occurred 15 h after the spill. The predicted time to reach the safe concentration varied from 3 to 24 h. The worst prediction of time to reach the safe concentration is only 20% greater than the deterministic prediction, whereas the best prediction is a factor of seven times less, due to the fact that the initial choice of model parameters was made to estimate low mixing conditions which represent the worst case for the dispersion of a spill. 3.3. Discussion It is expected that different assumptions for modelling horizontal mixing in an estuary and variability in the mixing parameters would lead to differences in the predicted concentrations. The simulations reported have quantified these differences for a practical case in the Forth estuary. It is interesting to see that the parameter uncertainty could lead to extremely rapid dilution with the safe limit being achieved within 3 h due to the highly dispersive conditions that occur in the estuary
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
277
Fig. 6. Representation of the model output in three-dimensional parameter space, UrK h rK z .
at certain times. At the other extreme, the model representing mixing as a function of the local tidal current leads to the most conservative prediction of
36 h to reach the safe concentration. Variability in the mixing rates in the lower Forth estuary has been established from the dye experiments carried out
Fig. 7. Maximum and minimum predicted concentrations throughout the spill Žsolid lines., results of the deterministic predictions Ždashed line. and the safe concentration Ždotted line..
278
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
between 1986 and 1988. It would therefore be more productive to carry out a long-term dispersion experiment Ž24 h or greater. to determine the actual course of dilution over a longer timescale. This could then be used to evaluate the different approaches to modelling the horizontal mixing. For planning dispersion experiments in a new area, the results of these simulations indicate that: 1. A number of short-term experiments should be carried out to estimate the variability in the mixing parameters for the region. 2. The addition of a longer term tracer experiment would provide valuable extra information to determine exactly how the mixing should be modelled over a longer time period.
4. Conclusions A random walk tidal model of the Forth estuary in Scotland has been used to investigate the uncertainty in predicted estuarine concentrations resulting from a discharge of material to the estuary. Two aspects of uncertainty are investigated: that due to the theoretical formulation used to model the horizontal mixing and that due to variations with time in the mixing coefficients. Different methods are reported in the literature for modelling the horizontal mixing in estuaries and coastal waters. It might be expected that different modelling approaches would lead to differences in predicted concentrations resulting from a chemical spill into an estuary. In this paper, simulations using three different formulations for the modelling of horizontal mixing are reported in an attempt to quantify the likely effects in a real estuary. In each case, the model has been calibrated against data from a 6-h dye dispersion experiment and has then been used to model the spread of a chemical in the estuary over a period of 48 h. Different rates of dilution were predicted but by 48 h the concentrations had once again converged because the patch had become mixed across the whole cross-section of the estuary. A maximum difference in the predicted concentrations was found to be 10 times and occurred 20 h after the
discharge. This uncertainty in the predicted concentration leads to an uncertainty in the time before a safe concentration is reached which, in this example, ranged from 10 to 36 h. A series of seven dye dispersion experiments carried out in the estuary between 1986 and 1988 has provided information on the variability of the mixing, and data from fixed station surveys at five positions in the estuary over a period of 2 days on spring tides and 2 days on neap tides ŽJunerJuly 1979. have given information on the variability of the tidal current. The random walk model was used to simulate the uncertainty in the predicted concentrations resulting from the variability in the measured mixing rates and tidal currents, for the same spill conditions considered above. Post-processing of the model output showed the distribution of possible concentrations for a range of times from discharge. The predicted time needed for the spill to be diluted to below the safe concentration has a range from 3 to 24 h. This study shows that model formulation and parameter uncertainty are both important factors in estimating the uncertainty in model predictions. The uncertainty caused by the variations with time found in the measured mixing rates was found to mainly improve predicted dispersion due to the conservative choice of mixing coefficients for the operational model. Uncertainties associated with model formulation lead to the conclusion that it could take longer for a chemical spill to be sufficiently diluted to be safe. However, it is possible this could be reduced if a small number of longer timescale Že.g. 24 h. dispersion experiments were available. In addition, further data from short-term Ž; 3 h. dispersion experiments would give a better understanding of the distribution of mixing coefficients and how the mixing relates to other parameters such as tidal range and wind speed and direction.
Acknowledgements I would like to thank Avecia Ltd., Grangemouth, for permission to use the Forth estuary model and the two referees for their helpful comments.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
References Bowden, K.F., 1983. Physical Oceanography of Coastal Waters. Ellis Horwood, Chichester, England. Elder, J.W., 1959. The dispersion of marked fluid in turbulent shear flow. J. Fluid Mech. 5 Žpart 4., 544–560. Elliott, A.J., Wallace, D.C., 1989. Dispersion of surface plumes in the southern North Sea. DT. Hydrogr. Z. 39, 114–137. Falconer, R.A., 1986. A two-dimensional mathematical model study of the nitrate levels in an inland natural basin. Water Quality Modelling in the Inland Natural Environment. Proceedings of the International Conference, Bournmouth, England, 10–13 June. Falconer, R.A., Chen, Y., 1991. An improved representation of flooding and drying and wind stress effects in a two-dimensional tidal numerical model. Proc. Inst. Civ. Eng. 91 ŽPart 2., 659–678, Dec. 1991. Fischer, H.B., List, E.J., Koh, R.C.Y., Imberger, J., Brooks, N.H., 1979. Mixing in Inland and Coastal Waters. Academic Press, New York. Hunter, J.R., Craig, P.D., Phillips, H.E., 1993. On the use of random walk models with spatially variable diffusivity. J. Comput. Phys. 106, 366–376. Joseph, J., Sender, H., 1958. Uber die horizontale Diffusion im Meere. Dt. Hydrogr. Z. 11, 49–77. Leendertse, J.J., 1967. Aspects of a computational model for long-period water-wave propagation. Memorandum RN-5294PR ŽMay 1967.. United States Air Force Project Rand, The Rand Corporation, Santa Monica, CA. Lewis, R.E., 1997. Dispersion in Estuaries and Coastal Waters. Wiley, Chichester, 312 pp. Ridderinkof, H., Zimmerman, J.T.F., 1992. Chaotic stirring in a tidal system. Science 258, 1107–1111.
279
Riddle, A.M., 1991. Spill model for the protection of the Forth estuary. In: Farmer, D.G., Rycroft, M.J. ŽEds.., Computer Modelling in the Environmental Sciences. Clarendon Press, Oxford, pp. 49–64. Riddle, A.M., 1998. Modelling tidal flows in the Forth estuary and over the intertidal mudflats at Grangemouth. In: McLusky, D.S. ŽEd.., Coastal Zone Topics: Process, Ecology and Management: 3. The Estuaries of Central Scotland. Joint Nature Conservation Committee and the Estuarine and Coastal Sciences Association, pp. 11–18. Riddle, A.M., Lewis, R.E., 2000. Dispersion experiments in U.K. coastal waters. Estuarine Coastal Shelf Sci. 51, 243–254. Rye, H., Reed, M., Ekrol, N., 1998. Sensitivity analysis and simulation of dispersed oil concentrations in the North Sea with the PROVANN model. Environ. Modell. Software 13, 421–429. Van Dam, G.C., Geurtz, M.J.J.M., 1994. Model particles as representatives of particles in nature. Neth. J. Aquat. Ecol. 28 Ž3–4., 317–328. Van Dam, G.C., Louwersheimer, R.A., 1992. A three-dimensional transport model for dissolved and suspended matter in estuaries and coastal seas. In: Prandle, D. ŽEd.., Dynamics and Exchanges in Estuaries and the Coastal Zone, vol. 40, Coastal and Estuarine Studies, Washington D.C., pp. 481–506. Van Veen, J., 1938. Water movements in the Straits of Dover. J. Cons., Cons. Perm. Int. Explor. Mer. 13, 7–38. Webb, A.J., 1982. A random walk model of the dispersion of caesium-137 in the Irish sea. Thesis for MSc, University of Wales, 1982. Webb, A.J., Metcalf, A.P., 1987. Physical aspects, water movements and modelling studies of the Forth estuary. The Natural Environment of the Estuary and Firth of Forth: Proceedings of the Royal Society of Edinburgh, vol. 93B, pp. 259–272.
Investigation of model and parameter uncertainty in water quality models using a random walk method A.M. Riddle ) Brixham EnÕironmental Laboratory, AstraZeneca UK Limited, Freshwater Quarry, Brixham, DeÕon, TQ5 8BA, UK Received 3 September 2000; accepted 12 February 2001
Abstract A mathematical model to predict the effect of chemical spills in the Forth estuary in Scotland has been in use for many years. The model, based on the random walk method, predicts chemical concentrations in the estuary waters and estimates the elapsed time before the dilution is sufficient to render the spill harmless Žmaking use of a toxicity measure such as the LC50 or a water quality standard.. The model gives a deterministic result without any estimate of the uncertainty. Field studies using tracer dyes to measure the horizontal and vertical mixing rates in the estuary show that these rates vary over time. The literature on turbulent diffusion includes modelling applications using different parameterisations of the mixing process. This paper investigates the uncertainties in predicted concentrations due to model parameterisation of horizontal mixing and due to the variability in the measured mixing rates determined from surveys in the estuary. Estimates of the range of concentrations for a specific spill scenario are presented. The study shows that model formulation and parameter uncertainty are both important factors in estimating the uncertainty in model predictions. The uncertainty caused by the variations with time found in the measured mixing rates is found to be of similar magnitude to the differences in concentration resulting from using three different methods for modelling the horizontal mixing in the estuary. Uncertainties associated with model formulation could be reduced if a small number of longer timescale Že.g. 24 h. dispersion experiments were available. In addition, further data from short-term Ž; 3 h. dispersion experiments would give a better understanding of the distribution of mixing coefficients and how the mixing relates to other parameters such as tidal range and wind speed and direction. q 2001 Elsevier Science B.V. All rights reserved. Keywords: Uncertainty; Random walk; Hydrodynamic; Model; Mixing
1. Introduction A model for predicting the effect of a chemical spill in the Forth estuary in Scotland was set up in 1990 ŽRiddle, 1991.. The model covered the lower estuary from the Forth bridges in the east to Alloa in )
Tel.: q44-1803-882-882; fax: q44-1803-882-974. E-mail address: [email protected] ŽA.M. Riddle..
the west ŽFig. 1., a distance of approximately 25 km. In this area, the estuary contracts rapidly from a width of approximately 3.5 km and depth of 20 m in the east to a width of 0.25 km and depth of 2.5 m at Alloa. The estuary is macro-tidal with a mean neap and spring tidal ranges of 2.5 and 5 m, respectively, and peak tidal currents of 1.4 m sy1 in the main channel off Grangemouth. The estuary is generally well mixed, however, a small degree of salinity stratification Žsurface to bottom salinity difference
0924-7963r01r$ - see front matter q 2001 Elsevier Science B.V. All rights reserved. PII: S 0 9 2 4 - 7 9 6 3 Ž 0 1 . 0 0 0 2 7 - 6
270
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Fig. 1. Lower Forth estuary, Scotland.
; 4 kg my3 . forms during the flood tide on neap tidal ranges. A more detailed description of the bathymetry, hydrodynamics and mixing characteristics of the estuary is given by Webb and Metcalf Ž1987.. The estuary model was based on diffusion of a Gaussian patch method to simulate the concentrations in the estuary, using interpolated field measurements of tidal current on a 200 = 200-m grid to simulate the transport processes. This original model has been enhanced to include a hydrodynamic simulation of the tidal currents using a 50 = 50-m grid and a three-dimensional random walk model to determine the dispersion ŽRiddle,
1998.. The hydrodynamic model is based on the finite difference solution of the depth averaged Navier–Stokes equations ŽLeendertse, 1967., with modifications to allow for flooding and drying of intertidal mudflats ŽFalconer and Chen, 1991.. The hydrodynamic model has been calibrated against a range of field data from the estuary including data from recording current meters, drogue tracking and dye movement from a range of tracer experiments ŽRiddle, 1998.. The water quality model uses the hydrodynamic data, stored at half hourly intervals throughout a tidal cycle together with a vertical current profile in a
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
three-dimensional particle tracking random walk approach ŽWebb, 1982; Van Dam and Louwersheimer, 1992; Van Dam and Geurtz, 1994. to simulate the chemical concentrations in the estuary ŽRiddle, 1998.. The calibration of dispersion in the water quality model was based on mixing coefficients determined from dye tracer experiments; the information consists of measurements of the dye patch size Žhorizontal and vertical. and concentration, taken over a relatively short timescale Ž6 h or less.. The model has subsequently been used for simulating effects that last for several tidal cycles. The horizontal scale of the model is small enough to resolve the current field and shear to produce chaotic motions such as described by Ridderinkof and Zimmerman Ž1992. for the Wadden Sea. The aim of this paper is to investigate the uncertainties in predicting concentrations in the estuary over a 48-h period resulting from a chemical discharge. Two aspects of uncertainty in predicted concentrations are considered in the paper:
271
obtained, which is comparable with the values obtained by Van Veen Ž1938. and employed by Van Dam and Louwersheimer Ž1992.. The use of a full three-dimensional hydrodynamic model is planned for future work, which should provide an improved simulation, particularly allowing for changes in flow direction through the water column in areas of rapidly changing water depth. The random walk model uses the appropriate tidal data from the hydrodynamic model at each timestep and represents the effluent discharge by placing a fixed number of AparticlesB at the outfall position on each timestep. These particles are moved on each subsequent timestep according to the following method ŽWebb, 1982.:
ž / .ž /
X new s Xold q U Ž p q 1 . Ynew s Yold q V Ž p q 1
p
Zold d
Zold d
dtqfx a
p
dtqfy b
Znew s Zold q f z g 1. The differences resulting from applying three different theoretical formulations, used by other researchers to model horizontal dispersion. 2. The differences arising from the uncertainty in the horizontal and vertical mixing rates and tidal current speeds, as shown by field measurements. The model theory is presented in Section 2, followed by a description of the uncertainty analysis of model formulation and parameter variability in Section 3. The results of the modelling study are discussed at the end of Section 3, and conclusions are presented in Section 4 of the paper.
2. Model theory
where X, Y and Z represent the position of a particle and subscripts ‘old’ and ‘new’ represent the positions at the start and end of a model timestep; d t is the timestep length, U and V are the depth averaged horizontal velocity components from the hydrodynamic model and d the local water depth. The functions f x , f y and f z define the mixing process and a , b and g are random numbers from a standard normal distribution. The vertical velocity is assumed to be zero. Each particle represents a fixed mass of effluent, and it is assumed that no degradation takes place. A constant diffusion coefficient, Fickian, model Žtheory K h ., has been used to parameterise the horizontal and vertical diffusion ŽBowden, 1983; Webb, 1982.:
(
f x s f y s 2 K hd t The three-dimensional random walk model uses tidal flow data from a two-dimensional depth averaged hydrodynamic model together with a prescribed vertical current profile in the direction of the flow ŽVan Dam and Louwersheimer, 1992; Van Dam and Geurtz, 1994.. The vertical current profile is expressed as a power law function and calibrated with local field data; an exponent Ž p . of 0.192 has been
(
fz s 2 Kzd t
Ž 1. Ž 2.
where K h and K z are the horizontal and vertical mixing coefficients Žm2 sy1 ., obtained from field studies. Coefficient K h represents the effect of turbulent mixing, since the shear is modelled directly by the spatial distribution in the hydrodynamic data and by the vertical current profile.
272
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
The concentration distribution of material in the estuary is estimated using a counting grid Žnot necessarily the same as the grid used by the hydrodynamic model.. The number of particles in a grid square
over a depth interval from the water surface down to a specified depth Žin this case a 1-m layer was used. is counted, giving the mass of material in a known volume, and therefore the concentration. An example model run has been carried out for a hypothetical spill of 25 kg over 0.5 h of a material with a safe concentration of 0.005 mg ly1 . Fig. 2a–c shows the movement and dilution of the material at 6, 12 and 24 h from discharge, whilst Fig. 2d shows the 95 percentile concentration of the patch plotted against time. The concentration is predicted to become safe after 20 h.
3. Uncertainty in model predictions 3.1. Mixing theory
Fig. 2. Predictions from the model for a discharge initially at position =; Ža. 6 h, Žb. 12 h and Žc. 24 h after discharge. Žd. Shows the time change in the predicted concentration Žthe dotted line represents the safe concentration..
The literature on modelling contaminant dispersion in estuaries and coastal areas shows a number of methods for parameterising the horizontal mixing. The application of these different methods can lead to a range of predicted concentrations and therefore to a range in the predicted time for concentrations to reach a safe level. The aim of this section is to quantify the differences resulting from using three different approaches to modelling the horizontal dispersion of a chemical spill in the Forth estuary, and to determine whether the formulations lead to a practical difference in the predictions. Three different theoretical formulations have been considered for the parameterisation of the horizontal mixing in the model; these are the following. Ø Fickian diffusion parameterised by a constant mixing coefficient as shown by Eq. Ž1. above; this is referred to as theory K h . Ø A mixing velocity formulation originally described by Joseph and Sender Ž1958. and used in an equivalent form by Elliott and Wallace Ž1989. and Rye et al., Ž1998.. This is parameterised as a constant diffusion velocity, Vm in Eq. Ž3. Žreferred to as theory V.. This method increases the horizontal mixing with time and simulates the effect of larger eddies contributing to the patch spread as the size of the patch increases. f x s f y s Vm 'td t
Ž 3.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Ø A method relating the mixing rate to the tidal current speed is based on the local application of the streamflow dispersion equation of Elder Ž1959. and is parameterised by a dispersion constant, k, in Eq. Ž4. Žsee Fischer et al., 1979, p. 125.. An example of the application of this formulation is given by Falconer Ž1986. in a two-dimensional depth averaged model for pollutant dispersion in Poole harbour England Žreferred to as theory U.: f x s f y s '2 kUHd t
Ž 4.
where U is the tidal velocity and H is the local water depth. This formulation specifies a spatially varying diffusion which can lead to the accumulation of particles in low diffusion areas; to compensate for this, a drift term dependent on the variation in the diffusion coefficient ŽHunter et al. 1993. has been incorporated into the model. The vertical mixing is parameterised by a constant mixing coefficient, K z , which simulates the vertical spread of material using Fickian diffusion ŽEq. Ž2.. as described by Lewis Ž1997.. Measurements of the vertical spread of dye patches in many areas around the UK ŽRiddle and Lewis, 2000. have shown applicability of this approach. In addition, the data for the Forth Žsee Section 3.2. indicate that mixing over depth will be complete in less than 4 h, which is a much shorter period than used for the model simulations. To calibrate each model results from a dye dispersion experiment, carried out in the Forth estuary on 13 July 1988 over a period of 6 h on the flood and ebb tides, were used. The data collected during the experiment consisted of horizontal patch size Žlength and width., patch depth and dye concentration. Beyond a time of 4 h the patch size data became unreliable due to difficulties in detecting the edges of the patch as it became dilute. The patch limits were determined using a submersible fluorometer with a detection limit of 0.001 mg ly1 . This has an influence of the fraction of the dye patch determined by the measurement, but in all cases the peak concentration was greater than a factor of 100 times the detection limit indicating that at least 99.75% of the patch had been measured Žassuming a Gaussian distribution of concentration in the horizontal direction..
273
Concentration data was only collected for the first 1.3 h of the experiment. The random walk model was used to simulate the dye dispersion in three different ways using the theoretical parameterisations for horizontal mixing represented by Eqs. Ž1., Ž3. and Ž4.. Eq. Ž2. was used to describe the vertical diffusion in each instance. All of the methods can adequately represent the measured dye patch spread and dilution. Fig. 3 shows the comparison of the methods against the measured data for dye patch size, expressed as standard deviation Ždefined as the square root of the patch length times the patch width. and the dye concentration. Each model gives a good prediction of the patch size over the 4-h timescale and the concentration over the shorter timescale of 1.3 h. For a model simulation that exceeds the length of the calibration period it might be expected that the three methods would give different concentration predictions. This has been investigated by running the calibrated models, using each of the mixing formulations, to simulate the rate of dilution of a chemical spill in the estuary. The conditions outlined in Section 2 have been used to predict concentrations over a 48-h period. Fig. 4 shows the three different predictions of concentration, which vary by up to a factor of 10 times Žat time 20 h after the spill.. By the end of four tides, the predictions have converged, because the effluent has mixed fully across the estuary cross-section area. The time required for the concentration to reach the safe concentration of 0.005 mg ly1 is 10, 22 and 36 h ŽFig. 4. for the diffusion velocity Ž V ., constant diffusion coefficient Ž K h . and tidal current ŽU . methods, respectively. Thus, the uncertainty in the duration of the incident is a factor of 3.6 times.
3.2. Mixing parameter Õalues In addition to uncertainty due to the theoretical modelling approach, the parameters used in simulation may be variable; for example, the tidal currents and mixing coefficients for the Forth estuary vary, as shown by the results of repeated studies in the estuary. This can lead to a significant variability in the predicted transport and dilution in the estuary
274
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
Fig. 3. Comparison of predicted and measured values for Ža. standard deviation of the horizontal patch size and Žb. peak dye concentration Žmg ly1 ..
from day to day. Data from seven different dye experiments ŽTable 1. show that the horizontal dis-
persion coefficient Ž K h . varies by a factor of 17 and vertical mixing coefficient Ž K z . by a factor of 16.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
275
Fig. 4. Predicted concentration Žmg ly1 . in the estuary from the three mixing models.
These data for horizontal and vertical mixing are shown in Fig. 5 as cumulative frequency plots against logŽ K h . and logŽ K z . ŽFig. 5a,b.. The quantity of data is small but the plots indicate a linear relationship between cumulative frequency and the log of the mixing coefficient. A positive lower bound would be expected for both parameters and minimum values have been estimated using linear regression fits to the data giving K hmin s 0.024 m2 sy1 and K zmin s 0.0002 m2 sy1 . These cumulative frequency plots correspond to an inverse relationship between the frequency of occurrence and the mixing coefficient. The random walk model has been used to simulate the effect of this variability on the spill dilution
and for this simulation the mixing has been simulated using the constant diffusion coefficient Ž K h . model. In the simulation, each model particle has been tagged with a horizontal and vertical mixing coefficient taken from the inverse distributions defined above. This is achieved as follows: LK h s LK hmin q a Ž LK hmax y LK hmin . LK z s LK zmin q b Ž LK zmax y LK zmin . where LK represents the log value of the mixing coefficient and a and b are random numbers from a uniform distribution between 0 and 1. For these
Table 1 Dispersion coefficients in the Forth estuary—from dye dispersion experiments Date
Time ŽBST.
Tide
Duration Žh.
Horizontal dispersion coefficient Ž K h m2 sy1 .
Vertical dispersion coefficient Ž K z m2 sy1 .
7r7r1986 7r7r1986 8r7r1986 9r7r1986 21r7r1986 13r7r1988 14r7r1988
0504 1127 0600 1212 1215 1120 0500
Ebb Flood Ebb Flood High water slack FloodrEbb Ebb
0.7 1.2 1.2 1.4 1.1 5.5 3.0
– – 0.18 0.05 0.16 0.75 0.85
0.0029 0.0008 0.0014 0.0009 – 0.0047 0.0003
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
276
Fig. 5. Cumulative frequency distributions Ž%. of the measured horizontal and vertical mixing values plotted against the logarithm of the mixing coefficient.
calculations, the regression equations have been used to compute values of 1.17 m2 sy1 for K hmax and 0.0047 m2 sy1 for K zmax . Repeated measurement of the tidal current at five positions in the estuary over 2 days on spring tides Ž25r26 June 1979. and 2 days on neap tides Ž2r3 July 1979. showed a variability of approximately 10% in the peak current amplitude in the model domain. This uncertainty in the tidal currents is modelled using the factor, C, to multiply the tidal currents: C s Cm q g Cs
Cm s 1
Cs s 0.05
where g is a random number from a standard normal distribution. Each particle has a unique values for the factor C as well as K h and K z throughout the simulation. At selected times, the information for
each particle was written to a file together with the particle position information. These data are then analysed using a postprocessing program to give concentration information for different ranges of the input parameters. Fig. 6 illustrates this process. The model output data, for a given time after the spill, can be represented in a three-dimensional space with axes K h , K z and U. This space can be divided up into elements; the particles within each element represent the chemical concentration and location for the given range Žsmall. of the parameter values. In Fig. 6, the right-hand element represents higher mixing conditions and gives rise to a larger but less concentrated patch of chemical than the left-hand element which represents lower mixing conditions. Taking the 95 percentile concentration from each element in the parameter space yields the overall concentration distribution from the uncertainty analysis. Fig. 7 shows the maximum and minimum predicted concentrations over the first 30 h of the simulation, together with the predicted concentration values obtained from a deterministic model run using the standard conservative values for the parameters. The largest range between the maximum and minimum predicted concentrations occurred 15 h after the spill. The predicted time to reach the safe concentration varied from 3 to 24 h. The worst prediction of time to reach the safe concentration is only 20% greater than the deterministic prediction, whereas the best prediction is a factor of seven times less, due to the fact that the initial choice of model parameters was made to estimate low mixing conditions which represent the worst case for the dispersion of a spill. 3.3. Discussion It is expected that different assumptions for modelling horizontal mixing in an estuary and variability in the mixing parameters would lead to differences in the predicted concentrations. The simulations reported have quantified these differences for a practical case in the Forth estuary. It is interesting to see that the parameter uncertainty could lead to extremely rapid dilution with the safe limit being achieved within 3 h due to the highly dispersive conditions that occur in the estuary
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
277
Fig. 6. Representation of the model output in three-dimensional parameter space, UrK h rK z .
at certain times. At the other extreme, the model representing mixing as a function of the local tidal current leads to the most conservative prediction of
36 h to reach the safe concentration. Variability in the mixing rates in the lower Forth estuary has been established from the dye experiments carried out
Fig. 7. Maximum and minimum predicted concentrations throughout the spill Žsolid lines., results of the deterministic predictions Ždashed line. and the safe concentration Ždotted line..
278
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
between 1986 and 1988. It would therefore be more productive to carry out a long-term dispersion experiment Ž24 h or greater. to determine the actual course of dilution over a longer timescale. This could then be used to evaluate the different approaches to modelling the horizontal mixing. For planning dispersion experiments in a new area, the results of these simulations indicate that: 1. A number of short-term experiments should be carried out to estimate the variability in the mixing parameters for the region. 2. The addition of a longer term tracer experiment would provide valuable extra information to determine exactly how the mixing should be modelled over a longer time period.
4. Conclusions A random walk tidal model of the Forth estuary in Scotland has been used to investigate the uncertainty in predicted estuarine concentrations resulting from a discharge of material to the estuary. Two aspects of uncertainty are investigated: that due to the theoretical formulation used to model the horizontal mixing and that due to variations with time in the mixing coefficients. Different methods are reported in the literature for modelling the horizontal mixing in estuaries and coastal waters. It might be expected that different modelling approaches would lead to differences in predicted concentrations resulting from a chemical spill into an estuary. In this paper, simulations using three different formulations for the modelling of horizontal mixing are reported in an attempt to quantify the likely effects in a real estuary. In each case, the model has been calibrated against data from a 6-h dye dispersion experiment and has then been used to model the spread of a chemical in the estuary over a period of 48 h. Different rates of dilution were predicted but by 48 h the concentrations had once again converged because the patch had become mixed across the whole cross-section of the estuary. A maximum difference in the predicted concentrations was found to be 10 times and occurred 20 h after the
discharge. This uncertainty in the predicted concentration leads to an uncertainty in the time before a safe concentration is reached which, in this example, ranged from 10 to 36 h. A series of seven dye dispersion experiments carried out in the estuary between 1986 and 1988 has provided information on the variability of the mixing, and data from fixed station surveys at five positions in the estuary over a period of 2 days on spring tides and 2 days on neap tides ŽJunerJuly 1979. have given information on the variability of the tidal current. The random walk model was used to simulate the uncertainty in the predicted concentrations resulting from the variability in the measured mixing rates and tidal currents, for the same spill conditions considered above. Post-processing of the model output showed the distribution of possible concentrations for a range of times from discharge. The predicted time needed for the spill to be diluted to below the safe concentration has a range from 3 to 24 h. This study shows that model formulation and parameter uncertainty are both important factors in estimating the uncertainty in model predictions. The uncertainty caused by the variations with time found in the measured mixing rates was found to mainly improve predicted dispersion due to the conservative choice of mixing coefficients for the operational model. Uncertainties associated with model formulation lead to the conclusion that it could take longer for a chemical spill to be sufficiently diluted to be safe. However, it is possible this could be reduced if a small number of longer timescale Že.g. 24 h. dispersion experiments were available. In addition, further data from short-term Ž; 3 h. dispersion experiments would give a better understanding of the distribution of mixing coefficients and how the mixing relates to other parameters such as tidal range and wind speed and direction.
Acknowledgements I would like to thank Avecia Ltd., Grangemouth, for permission to use the Forth estuary model and the two referees for their helpful comments.
A.M. Riddler Journal of Marine Systems 28 (2001) 269–279
References Bowden, K.F., 1983. Physical Oceanography of Coastal Waters. Ellis Horwood, Chichester, England. Elder, J.W., 1959. The dispersion of marked fluid in turbulent shear flow. J. Fluid Mech. 5 Žpart 4., 544–560. Elliott, A.J., Wallace, D.C., 1989. Dispersion of surface plumes in the southern North Sea. DT. Hydrogr. Z. 39, 114–137. Falconer, R.A., 1986. A two-dimensional mathematical model study of the nitrate levels in an inland natural basin. Water Quality Modelling in the Inland Natural Environment. Proceedings of the International Conference, Bournmouth, England, 10–13 June. Falconer, R.A., Chen, Y., 1991. An improved representation of flooding and drying and wind stress effects in a two-dimensional tidal numerical model. Proc. Inst. Civ. Eng. 91 ŽPart 2., 659–678, Dec. 1991. Fischer, H.B., List, E.J., Koh, R.C.Y., Imberger, J., Brooks, N.H., 1979. Mixing in Inland and Coastal Waters. Academic Press, New York. Hunter, J.R., Craig, P.D., Phillips, H.E., 1993. On the use of random walk models with spatially variable diffusivity. J. Comput. Phys. 106, 366–376. Joseph, J., Sender, H., 1958. Uber die horizontale Diffusion im Meere. Dt. Hydrogr. Z. 11, 49–77. Leendertse, J.J., 1967. Aspects of a computational model for long-period water-wave propagation. Memorandum RN-5294PR ŽMay 1967.. United States Air Force Project Rand, The Rand Corporation, Santa Monica, CA. Lewis, R.E., 1997. Dispersion in Estuaries and Coastal Waters. Wiley, Chichester, 312 pp. Ridderinkof, H., Zimmerman, J.T.F., 1992. Chaotic stirring in a tidal system. Science 258, 1107–1111.
279
Riddle, A.M., 1991. Spill model for the protection of the Forth estuary. In: Farmer, D.G., Rycroft, M.J. ŽEds.., Computer Modelling in the Environmental Sciences. Clarendon Press, Oxford, pp. 49–64. Riddle, A.M., 1998. Modelling tidal flows in the Forth estuary and over the intertidal mudflats at Grangemouth. In: McLusky, D.S. ŽEd.., Coastal Zone Topics: Process, Ecology and Management: 3. The Estuaries of Central Scotland. Joint Nature Conservation Committee and the Estuarine and Coastal Sciences Association, pp. 11–18. Riddle, A.M., Lewis, R.E., 2000. Dispersion experiments in U.K. coastal waters. Estuarine Coastal Shelf Sci. 51, 243–254. Rye, H., Reed, M., Ekrol, N., 1998. Sensitivity analysis and simulation of dispersed oil concentrations in the North Sea with the PROVANN model. Environ. Modell. Software 13, 421–429. Van Dam, G.C., Geurtz, M.J.J.M., 1994. Model particles as representatives of particles in nature. Neth. J. Aquat. Ecol. 28 Ž3–4., 317–328. Van Dam, G.C., Louwersheimer, R.A., 1992. A three-dimensional transport model for dissolved and suspended matter in estuaries and coastal seas. In: Prandle, D. ŽEd.., Dynamics and Exchanges in Estuaries and the Coastal Zone, vol. 40, Coastal and Estuarine Studies, Washington D.C., pp. 481–506. Van Veen, J., 1938. Water movements in the Straits of Dover. J. Cons., Cons. Perm. Int. Explor. Mer. 13, 7–38. Webb, A.J., 1982. A random walk model of the dispersion of caesium-137 in the Irish sea. Thesis for MSc, University of Wales, 1982. Webb, A.J., Metcalf, A.P., 1987. Physical aspects, water movements and modelling studies of the Forth estuary. The Natural Environment of the Estuary and Firth of Forth: Proceedings of the Royal Society of Edinburgh, vol. 93B, pp. 259–272.