Ocean model simulations of sediment plume behaviour

(1)

Volume Two May 2014

Appendix 23

Ocean model simulations of sediment plume behaviour (Hadfield 2013)

(2)

Ocean model simulations of sediment plume behaviour

Prepared for Chatham Rock Phosphate Limited

March 2011 (Updated April 2013)

(3)

© All rights reserved. This publication may not be reproduced or copied in any form without the permission of the copyright owner(s). Such permission is only to be given in accordance with the terms of the client’s contract with NIWA. This copyright extends to all forms of copying and any storage of material in any kind of information retrieval system.

Whilst NIWA has used all reasonable endeavours to ensure that the information contained in this document is accurate, NIWA does not give any express or implied warranty as to the completeness of the information contained herein, or that it will be suitable for any purpose(s) other than those

specifically contemplated during the Project or agreed by NIWA and the Client.

Authors/Contributors:

Mark Hadfield

For any information regarding this report please contact:

[email protected] Contracts Manager +64-4-386 0369 Neville Ching

National Institute of Water & Atmospheric Research Ltd 301 Evans Bay Parade, Greta Point

Wellington 6021

Private Bag 14901, Kilbirnie Wellington 6241

New Zealand

Phone +64-4-386 0300 Fax +64-4-386 0574

NIWA Client Report No: WLG2010-71

Report date: March 2011(Updated April 2013)

NIWA Project: WSE11301

(4)

Executive summary

A limited-area, three-dimensional ocean Chatham Rock Phosphate Limited model has been set up for a 500 × 440 km domain over central Chatham Rise, encompassing the

Widespread Energy mineral permit area. The model was forced by ocean eddies at its lateral boundaries and by winds and (for the sediment transport calculations) by tides. By

comparison with measurements it is found to give an approximate reproduction of the true currents over Chatham Rise.

Idealised sediment plumes were introduced into the model, at three locations at the centre and the eastern and western edges of the mineral permit area. The sediment tracer was implemented in the form of virtual particles. Particles were released at several different heights at each location and at seven different periods during the simulation, in order to sample the model’s temporal variability. Particles sank relative to the modelled fluid at four different sinking-velocities: zero and 10−3, 10−4 and 10−5 m s−1 (or 86.4, 8.64 and 0.864 m d−1).

The sediment-plume simulations give semi-quantitative and illustrative information about aspects of sediment dispersion in the area. Important aspects of the problem are found to be:

 the role of the mean flow and eddies in advecting and stirring the particle plumes;

 vertical mixing within ~100 m of the surface and ~50 m of the bottom by wind- driven shear and tidal bottom stress, respectively;

 the importance of fall speed, with a fall speed of 10−3 m s−1 (86.4 m d−1) having marked effects on plume behaviour and a fall speed of 10−4 m s−1 (8.64 m d−1) having modest effects.

(7)

1 Introduction

Building on previous NIWA work in the area (Hadfield et al. 2007) a three-dimensional, hydrodynamic model has been set up for a 500 × 440 km domain at 2 km horizontal resolution over central Chatham Rise (Figure 1). The model is forced by:

 currents, temperature and salinity at its lateral boundaries from a larger-scale ocean model of the New Zealand region;

 surface stresses and heat fluxes from a global atmospheric model, the NCEP Reanalysis.

A base simulation was run for two model years, producing time-varying temperature, salinity and velocity fields that were saved regularly. There was also a series of sediment-plume simulations, described below. These were run for shorter periods, initialised from the base simulation and exposed to the same forcing, plus

 tidal fluctuations in velocity and sea surface height from a NZ region tidal model, also applied at the lateral boundaries.

For simplicity, the tidal forcing consisted of only one constituent, the lunar, semi-diurnal tide (M2).

The model was ROMS (Haidvogel et al. 2008), a widely used ocean/coastal model with the option of an embedded model of suspended sediment and sediment bedload processes (Warner et al. 2008).

Figure 1: Bathymetry of Chatham Rise with an outline of the ocean model grid.

(8)

2 The base model

For an indication of how well the model reproduces the true currents over the Chatham Rise, model output statistics are compared below with satellite altimeter data. (The slope of the sea surface, which can be estimated by satellite altimeters, is related to surface currents by the geostrophic relationship.) Figure 2 shows the mean surface height from 700 days of the base model run, compared with an estimate of the same quantity synthesising satellite altimeter and other data. Strictly speaking, the quantity in the lower panel is the mean dynamic topography (MDT) of the ocean surface and the specific product is the best one currently available, MDT_CNES-CLS091 version 1.1, produced in March 2010. Both sea surface height fields are relative to an arbitrary datum, so the absolute values are not expected to agree, but ideally the differences will agree. In both the model output and the observation-based product, there is a difference of almost 0.5 m from south to north of Chatham Rise, implying a broadly eastward surface flow. On the southern side of the Rise, the sea surface height gradient implies a flow to the ESE, and here the model matches the observations quite well. On the northern side of the Rise, the observations show two local maxima near the northern boundary of the figure, at 178.5° E and 177° W, extending ridges into the model domain. Somewhat similar ridges appear in the model output, but the eastern one is too far west, implying that the pattern of the mean currents on the northern side of Chatham Rise is not right in the model.

A similar comparison is presented in Figure 3, but here the quantity is the geostrophic eddy kinetic energy, an estimate of the intensity of the variability in the surface currents, but not excluding tides or short-term wind-driven fluctuations2. Both the model and the observations have a minimum along the Rise and larger values north and south. However the model has a somewhat different pattern, particularly north of Chatham Rise, where there are large values in the northwestern corner of the model domain, with a tongue of elevated values extending southeastward from there.

High-quality, long-duration, in situ measurements of the currents in the Chatham Rise area are rare, but some data are available from a set of four current meters installed in 1996–1997 (Chiswell 2001). There were two sites, on the northern (178.63° E, 42.70° S) and southern (178.61° E, 44.61° S) flanks of the Rise at the 1500 m isobath and meters were installed at 250 m and 1000 m depths. The record length was between 140 and 365 days. Scatter plots of the detided currents are compared with colocated model data in Figure 4. The red line in each panel represents the variance ellipse: the centre of the ellipse shows the mean current vector and the size and shape of the ellipse represents the variability about that mean.

Agreement in the mean is reasonably good except for the upper current meter at the northern site where the observed mean flow vector is directed eastward at 0.05 m s−1 whereas the modelled one is directed southward at 0.05 m s−1. Modelled variability (the size of the ellipse) is usually, but not always, a little low.

1 MDT_CNES-CLS09 was produced by CLS Space Oceanography Division and distributed by Aviso, with support from Cnes (http://www.aviso.oceanobs.com/)

2The altimeter data from which the geostrophic eddy kinetic energy was calculated were produced by Ssalto/Duacs and distributed by Aviso, with support from Cnes

(http://www.aviso.oceanobs.com/duacs/).

(9)

8 Ocean model simulations of sediment plume behaviour It is important to note that because the larger-scale model, which provides the lateral

boundary conditions, does not assimilate data3 the model cannot reproduce the fluctuations in currents in the area at the time they actually occurred.

There are various ways in which the model could be developed to improve its match with reality. The first option to investigate would be to generate lateral boundary conditions from the Bluelink Reanalysis, a data-assimilating ocean model and analysis system maintained by a consortium of Australian public agencies.

Given the limitations of the data against which we can validate it, we can say that this model gives an approximate reproduction of true currents over Chatham Rise. It has eastward jets on the northern and southern flanks of the Rise, and weaker, generally eastward currents on top of the Rise. In the mineral permit area the mean surface current (Figure 5 top) follows the bathymetry and is directed towards the east or southeast, turning to the northeast at the eastern edge of the area. The mean bottom current (Figure 5 bottom) is weaker and directed to the southeast.

3 In the open ocean, much of the variability is associated with eddies that are generated by instability processes and to reproduce these features at the correct place and time a model must assimilate data (esp. from satellite altimeters). In the absence of data assimilation, the best a model can do is to simulate eddies with the correct statistical properties.

(10)

Figure 2: Time-mean sea surface height from the model (upper) compared with a mean dynamic topography (lower) estimated from satellite altimeter and other data—see text for details

(11)

10 Ocean model simulations of sediment plume behaviour Figure 3: The logarithm (base 10) of the geostrophic eddy kinetic energy from the model (upper) and from satellite altimeter data (lower).

(12)

(13)

12 Ocean model simulations of sediment plume behaviour Figure 4: Scatter plots of detided velocities from current meters (left) and model (right) with variance ellipses.

Figure 5: Time-mean velocity vectors from the base-model at the uppermost (top) and lowest (bottom) model levels.

(14)

3 Sediment-plume

To simulate the transport and deposition of sediment produced by hypothetical dredging operations, we introduce a tracer into the model. Broadly speaking, there are two

approaches to doing this. One approach is to treat each of the classes of sediment as a continuous quantity, similar to the ones the model already deals with (temperature and salinity). This is the approach that was taken for the fine-scale sediment modelling with Gerris, described elsewhere. The other approach is to treat the sediment as a collection of virtual particles that are moved around by the modelled velocity fields. Both approaches can be used with ROMS; for this exercise the virtual particle technique was used as it easier to set up and faster to run for exploratory simulations with multiple sources.

4 Virtual-particle simulations

Three locations were considered, at the centre and the eastern and western edges of the mineral permit area (Figure 6). At each location there were six sources, at heights of 4, 12.5, 25 and 50 m above the bottom, and 200 and 10 m below the surface, each releasing a stream of particles. There was a series of seven of these simulations, starting at t_{= 100,} 200, …, 700 days relative to 1 January on the first model year. Each was forced by tides at its boundary, in addition to the other forcings, and had its hydrodynamic state initialised from a snapshot of the two year (non-tidal) simulation. Release of particles started 5 days after the beginning of the simulation (giving time for the tides to develop) and continued for 10 days, with one particle released from each source every 30 minutes. The particle-release

simulation then continued until 10 days after the final particle was released, for a total duration of 25 days.

The virtual particles in these simulations followed the three-dimensional, time-varying currents in the model but with three extra pieces of “behaviour”:

 Each particle was subjected to a vertical random walk (i.e., at each time step it was displaced vertically by a random amount) to simulate the effect of turbulent diffusion (Hunter et al., 1993). The turbulent diffusion rates were calculated by the hydrodynamic model using the ROMS adaptation of the Large et al. (1994) turbulence scheme, which allows for both near-surface and near-bottom turbulent boundary layers, and also allows for turbulence to develop in the interior of the ocean (though this does not usually happen).

 Particles sank relative to the modelled fluid at a rate that was fixed for each simulation. Four different sinking-velocities were used: zero and 10−3, 10−4 and 10−5 m s−1 (or 86.4, 8.64 and 0.864 m d−1). Compared with settling velocities of the sediments from the DISCOL experimental area (Thiel et al., 2001) the three non-zero sinking velocities correspond to grain sizes of approximately 500, 40 and 20 μm (Jankowski et al. 1996).

 Particles that contacted the bottom in the model were considered to have been deposited, i.e., they remained immobile thereafter. Resuspension was not considered.

(15)

14 Ocean model simulations of sediment plume behaviour With seven starting times and four settling velocities, there were 28 particle simulations, each involving 8658 particles.

Figure 6 shows particle trajectories from one of the particle simulations, one with zero sinking velocity beginning at 300 days. Particles released near the surface (upper panel) are swept in different directions from each source by ocean eddies. The cycloidal, or wavy, pattern that can be seen in some of the trajectories (but is actually present in all of them) arises because of the tidal oscillations, which move the particles in an ellipse of approximately 2.5 km in size.

Particles released near the bottom (lower panel) follow broadly similar paths but do not move so far and are not so widely dispersed. The question of horizontal dispersion will be

examined in more detail below with reference to scatter plots of particle positions after 10 days.

(16)

Figure 6: Trajectories up to 10 days after release of non-sinking particles from the three source locations in a single simulation (305–325 d) at 10 m below the surface (upper) and 50 m above the bottom (lower).

(17)

16 Ocean model simulations of sediment plume behaviour Vertical dispersion of the particles is illustrated in a series of figures (Figure 7 to Figure 10) each showing particle height versus time after release for the six different release heights (the six panels in each figure) and for sets of simulations with different sinking velocities (the different figures). Each panel in these figures overlays the time series from all the particles with the given release height and sinking velocity, from three different source locations and seven different simulations and at different release times, giving 10101 time series in each panel. With zero sinking velocity (Figure 7), particles released 10 m below the surface are mixed rapidly through the surface wind-mixed layer, which is up to 140 m thick. (It is often thinner than this, but the figure emphasises the outer envelope of particle heights.) Particles released within 25 m of the bottom are influenced by a similarly turbulent near-bottom layer, typically 30 m thick, this one driven by shear from tidal currents. In between the near-surface and near-bottom turbulent layers, vertical mixing is much weaker (in the model, and probably in the real world). Particles released 200 m below the surface remain within a relatively narrow vertical layer; the slight acceleration in the growth of this layer that occurs at 7 days is not a result of turbulent mixing, but of vertical divergence of particle trajectories as the flow they are in encounters varying bathymetry. Particles released 50 m above the bottom are not immediately affected by the near-bottom turbulence, but some are eventually mixed down after 8 or so days. Again, much of the vertical growth of the layer of particles is a result of the varying bathymetry.

The effect of vertical sinking at the highest velocity considered, 10−5 m s−1, are dramatic (Figure 8). Of the particles released near the surface, most sink to the bottom after 4–6 days;

the particles that take the longest are the ones that have been detained in the near-surface turbulent layer the longest. Particles released at 200 m depth all take 2.5 days to reach the bottom Particles released at 50 m or less above the bottom sink to the bottom within one day, some being briefly held up by turbulent mixing.

When the sinking rate is reduced by a factor of 10, to 10−4 m s−1 (Figure 9) the effects of sinking are qualitatively similar but much less marked. The plume of particles released at 200 m depth sinks by ~ 80 m over 10 days, as expected. Of the particles released 10 m below the surface, some sink well below the surface turbulent layer, but many remain close to the surface. For particles released within 50 m of the bottom, the plume settles down to a quasi-steady state, mixed up to 20 m above the bottom, but progressively depleted as particles encounter the surface and are deposited.

When the sinking rate is further reduced to 10−5 m s−1 (Figure 10) the effect of the sinking is barely perceptible in these graphs (cf. Figure 6).

(18)

Figure 7: Particle height vs time for particles with zero sinking velocity released at 6 different heights. Release heights are 10 and 200 m below the surface and 50, 25, 12.5 and 4 m above the bottom.

(19)

14 Ocean model simulations of sediment plume behaviour Figure 8: As Figure 7 but for a sinking velocity of 10−3 m s−1.

(20)

Figure 9: As Figure 7 but for a sinking velocity of 10−4 m s−1.

(21)

16 Ocean model simulations of sediment plume behaviour Figure 10: As Figure 7 but for a sinking velocity of 10−5 m s−1.

(22)

Figure 11 to Figure 14 are organised in the same way as Figure 7 to Figure 10, but show the particle horizontal locations, relative to the source, at 10 days after release. For the release 10 m below the surface of particles with zero sinking velocity (Figure 11) the particle cloud is dispersed quite widely, up to 170 km from the source, in a generally eastward direction. The generally streaky and lumpy nature of the cloud of particles would doubtless smooth out with more releases in different flow conditions, but it indicates the patchy nature of dispersion by eddies. For the release at 200 m depth the particle cloud is more compact, with a maximum distance of 90 km. For the near-bottom releases, the particle clouds are more compact again and much the same size and shape for all these release heights; as the release height decreases the fraction of deposited particles (blue symbols) increases,

A sinking velocity of 10−3 m s−1 (Figure 12) has a marked effect on horizontal dispersion. By 10 days after the release, all particles at all release heights have been deposited (as was apparent from Figure 8) and the horizontal spread is consequently much reduced. For the near-bottom releases, the elliptical shape of the cloud of deposited particles is a result of tidal oscillations, as the particle travel times are less than one tidal period.

A sinking velocity of 10−4 m s−1 (Figure 13) has much more modest effects on horizontal dispersion. For the releases at 10 and 200 m depth, the particle cloud is of a similar size and shape to that seen in the zero-sinking case (Figure 11). For the near-bottom releases, most of the particles are deposited before 10 days, which limits horizontal spread.

With the lowest sinking velocity of 10−5 m s−1 (Figure 14), horizontal dispersion is very similar to that seen in the zero-sinking case (Figure 11) but does increase the number of particles deposited for releases near the bottom

(23)

18 Ocean model simulations of sediment plume behaviour Figure 11: Scatter plot of horizontal locations at 10 days for particles with zero sinking

velocity released at 6 different heights. Release heights are 10 and 200 m below the surface and 50, 25, 12.5 and 4 m above the bottom. Deposited particles are represented by blue symbols and non-deposited particles by red symbols.

(24)

Figure 12: As Figure 11 but for a sinking velocity of 10−3 m s−1.

(25)

20 Ocean model simulations of sediment plume behaviour Figure 13: As Figure 11 but for a sinking velocity of 10−4 m s−1.

(26)

Figure 14: As Figure 11 but for a sinking velocity of 10−5 m s−1.

(27)

5 Conclusions

The hydrodynamic simulation with virtual particle releases has given semi-quantitative and illustrative information about aspects of sediment dispersion at the Widespread Energy mineral permit area. Important aspects of the problem include:

 the role of the mean flow and eddies in advecting and stirring the particle plumes;

 vertical mixing within ~ 100 m of the surface and ~ 50 m of the bottom by wind- driven shear and tidal bottom stress, respectively;

 the importance of fall speed, with a fall speed of 10−3 m s−1 (86.4 m d−1) having marked effects on plume behaviour and a fall speed of 10−4 m s−1 (8.64 m d−1) having modest effects.

Behaviour of real sediment plumes released by undersea mining is more complicated than the virtual particles with specified sinking velocities considered here. For example, a more complete model would consider (Jankowski and Zielke 2001):

 the effect of the suspended particle load on the density of the fluid and thereby on the flow;

 particle flocculation.

However with the addition of quantitative information about source strengths and sediment parameters the current model would give preliminary indications of sediment concentrations and deposition rates.

(28)

6 References

Chiswell, S.M. (2001). Eddy energetics in the Subtropical Front over the Chatham Rise. New Zealand Journal of Marine and Freshwater Research 35: 1-16.

Hadfield, M.G.; Rickard, G.J.; Uddstrom, M.J. (2007). A hydrodynamic model for Chatham Rise, New Zealand. New Zealand Journal of Marine and Freshwater Research 41: 239-264.

Hunter, J.R.; Craig, P.D.; Phillips, H.E. (1993). On the use of random walk models with spatially variable diffusivity. Journal of Computational Physics 106(2): 366- 376.

Haidvogel, D.B.; Arango, H.; Budgell, W.P.; Cornuelle, B.D.; Curchitser, E.; Di Lorenzo, E.; Fennel, K.; Geyer, W.R.; Hermann, A.J.; Lanerolle, L.; Levin, J.;

McWilliams, J.C.; Miller, A.J.; Moore, A.M.; Powell, T.M.; Shchepetkin, A.F.;

Sherwood, C.R.; Signell, R.P.; Warner, J.C.; Wilkin, J. (2008). Ocean forecasting in terrain-following coordinates: Formulation and skill assessment of the Regional Ocean Modeling System. Journal of Computational Physics 227(7): 3595–3624.

Jankowski, J.A.; Malcherek, A.; Zielke, W. (1996). Numerical modeling of

suspended sediment due to deep-sea mining. Journal of Geophysical Research 101(C2): 3545-3560.

Jankowski, J.A.; Zielke, W. (2001). The mesoscale sediment transport due to technical activities in the deep sea. Deep Sea Research Part II: Topical Studies in Oceanography 48(17-18): 3487-3521.

Large, W.G.; McWilliams, J.C.; Doney, S.C. (1994). Oceanic vertical mixing: a review and a model with a nonlocal boundary layer parameterization. Reviews of Geophysics 32(4): 363-403.

Thiel, H.; Schriever, G.; Ahnert, A.; Bluhm, H.; Borowski, C.; Vopel, K. (2001). The large-scale environmental impact experiment DISCOL--reflection and foresight.

Deep Sea Research Part II: Topical Studies in Oceanography 48(17-18): 3869- 3882.

Warner, J.C.; Sherwood, C.R.; Signell, R.P.; Harris, C.K.; Arango, H.G. (2008).

Development of a three-dimensional, regional, coupled wave, current, and sediment-transport model. Computers & Geosciences 34(10): 1284-1306.

Ocean model simulations of sediment plume behaviour

Appendix 23

Ocean model simulations of sediment plume behaviour

Prepared for Chatham Rock Phosphate Limited

Contents

Executive summary

1 Introduction

2 The base model

3 Sediment-plume

4 Virtual-particle simulations

5 Conclusions

6 References