Aquatic Ecology 38: 145–162, 2004. © 2004 Kluwer Academic Publishers. Printed in the Netherlands.

145

Lake circulation and sediment transport in Lake Myvatn Snorri P´all Kjaran, Sigurdur L´arus H´olm and Eric Matthew Myer Vatnaskil Consulting Engineers, Sudurlandsbraut 50, 108 Reykjavik, Iceland (e-mail: [email protected])

Key words: AQUASEA, Diatomite, Numerical modelling, Suspended sediment, Vatnaskil

Abstract Lake circulation and sediment transport in Lake Myvatn have been calculated using AQUASEA, a numerical model developed by Vatnaskil Consulting Engineers. The goal of the modelling was to calculate changes in sediment transport within the lake due to changes in lake bathymetry caused by diatomite mining. The model uses the Galerkin finite element method and consists of a hydrodynamic flow model and a transport-dispersion model. The flow model is based on the shallow water equations and the wave equation. The transport model is based on the conservation of mass for suspended sediment. The model was calibrated against measurements performed during the summer of 1992. These included measurements of water elevation, current velocity, wave height, and concentration of suspended sediment. After calibration, the model was run for different mining scenarios to determine their impact on the sediment transport in the lake.

Introduction Lake Myvatn, a shallow eutrophic lake in northern Iceland, is renowned for its abundant wildlife, especially waterfowl. The wildlife depends to a large extent on a rich fauna of deposit-feeding benthic invertebrates, of which chironomid midges and chydorid Cladocera are of fundamental importance (Jónasson 1979, Lindegaard and Jónasson 1979, Gardarsson 1979, Einarsson et al. 2002). Since 1967 dredging of the bottom deposits has been taking place in Ytrifloi (North Basin) for the manufacturing of diatomite, a fine powder of diatom frustules used for various industrial purposes. About 275,000 m3 of sediment are removed each year causing an increase in depth from about 1 m to between 2 to 6 m. The resources in Ytrifloi are now almost exhausted and plans have been made to expand the mining to the larger and somewhat deeper (up to 4 m) Sydrifloi (South Basin). Episodic wind events usually cause wave-induced resuspension of sediments in shallow lakes (Luettich et al. 1990; Douglas and Rippey 2000). There are two types of currents in Lake Myvatn, currents forced by inflow and outflow to and from the lake, which are relatively small in Lake Myvatn, and currents due to wind shear stress at the

surface of the lake, which can be significant during ice-free periods with high wind. The wind-induced currents can play a major role in the lateral transport of sediment (Sheng and Lick 1979). Mining pits in the flat, wave-exposed bottom of Lake Myvatn will inevitably trap fresh resuspended organic sediment derived from the surrounding bottom areas and possibly reduce the resources of the depositfeeding invertebrates. Preliminary studies of currents, resuspension, settling velocities and sediment characteristics (published in a series of internal reports by H. Jóhannesson et al.) indicated that such sediment focusing might reach the same order of magnitude as the annual net production of fresh organic sediment (detritus) in the lake basin. This might present a risk to the detritus-driven ecosystem of the lake. Hence, a more rigorous numerical modelling approach was called for to further investigate the lake circulation and sediment transport in the lake. Modelling of Lake Myvatn began in 1991 with the calculation of lake circulation for isolated strong wind events during the summer of 1990 (Vatnaskil 1991). In 1993 the lake circulation model was updated and sediment transport was first calculated (Vatnaskil 1993). This work was based on measurements proposed and funded by the Project Management Team for Research

146

Figure 1. Inflow and outflow to and from Lake Myvatn. From J´ohannesson (1991a).

on Lake Myvatn. The goal of the modelling was to calculate the effect of proposed diatomite mining by Kísilidjan Inc. at Bolir, the easternmost part of the Sydrifloi basin. The results of this modelling were published in an article later that year (Tómasson and Kjaran 1993). Further modelling of sediment transport was conducted for Kísilidjan Inc. to investigate proposed diatomite mining scenarios (Vatnaskil 1998 and 2000). The present paper describes the hydrodynamic modelling of Lake Myvatn and its application in predicting the impact of mining on the sedimentation dynamics of the lake.

Materials and methods Lake Myvatn is located in northeast Iceland, at 277 meters above sea level and latitude 65◦ 35 N. The lake is relatively shallow, with a maximum depth of 4.2 meters. The lake covers an area of 37 km2 and consists of two major basins, Ytrifloi and Sydrifloi. Inflow to the lake occurs mainly at groundwater springs along the eastern shore and outflow from the lake discharges into the River Laxá (Figure 1). The numerical model AQUASEA, which has been developed by Vatnaskil Consulting Engineers, was

147 used to model the currents and sediment transport in Lake Myvatn. The model consists of three parts. First, a hydrodynamic flow model is used to calculate the shallow water equations. It simulates water level variations and flows in response to wind forcing at the lake surface and inflow and outflow to and from the lake. The waterlevels and flows are approximated in a numerical finite element grid and calculated on the basis of information on the bathymetry, bed resistance coefficients, wind field and boundary conditions (Gray and Kinnmark 1986; Gray and Kinnmark 1987; Kjaran et al. 1988). Second, a simple wave model based on empirical and theoretical equations is used to predict the surface wave environment in the lake and the associated bottom shear stress. Finally, a transportdispersion model, which is based on the equation of conservation of mass, is used to calculate the transport and deposition of suspended sediments. It simulates the spreading of a substance in the environment under the influence of the fluid flow and the existing dispersion processes. Flow model The shallow water equations are based on the assumption that the depth of flow is small compared with the horizontal length scales involved, resulting in two-dimensional depth-averaged equations. Due to the shallowness of Lake Myvatn, this assumption is reasonable. The equation of continuity, which expresses conservation of mass, is given by ∂ ∂η ∂ (uH ) + (vH ) + =Q ∂x ∂y ∂t

(1)

where H =h+η

(2)

Explanation: h Mean water depth, m η Change in waterlevel, m H Total water depth, m u Velocity component in x-direction, m/s v Velocity component in y-direction, m/s t Time, s Q Injected water, m3 /s By assuming hydrostatic pressure variation across all verticals, the momentum equations which express con-

servation of momentum in the x and y directions are given by ∂u ∂u ∂u +u +v = ∂t ∂x ∂y ∂η g (3) −g (u2 + v 2 )1/2 u + fv − ∂x H C2 k Q + Wx |W | − u H H and ∂v ∂v ∂v +u +v = ∂t ∂x ∂y g ∂η (4) − fu − −g (u2 + v 2 )1/2 v ∂y H C2 k Q + Wy |W | − v H H The Coriolis parameter, f, is defined as f = 2ω sin ϕ

(5)

where ϕ is the latitude and ω is the Earth’s rate of rotation equal to 7.2722·10−5 s−1 . The wind shear stress parameter, k, is defined as ρa CD k= (6) ρ Explanation: g Acceleration of gravity, m/s2 ω The Earth’s rate of rotation, s−1 ϕ Latitude, deg C Chezy bottom friction coefficient, m1/2 /s ρa Density of air, kg/m3 CD Wind drag coefficient ρ Fluid density, kg/m3 Wx Wind velocity in x-direction, m/s Wy Wind velocity in y-direction, m/s |W| Wind speed, m/s The boundary conditions for the flow model are as follows. The shoreline of the lake is defined as a noflow boundary. The outflow from the lake to the River Laxá is assumed to be constant in time and equal to 33 m3 /s, which is the average measured discharge to the river. The distribution of the inflow to the lake, which is mostly from groundwater, is shown in Figure 1. The total inflow is 33 m3 /s, equal to the outflow to the River Laxá. Therefore, the model is forced by both the inflow/outflow to and from the lake and the external wind stress at the surface, which is based on wind measurements at Geiteyjarstrond. The bottom shear stress due to the currents is given by ρg (7) τ = 2 (u2 + v 2 )1/2 (u, v) C

148

Figure 2. Measurement locations.

The above equations are solved numerically using a Galerkin finite element method based on triangular elements (Gray and Kinnmark 1986; Gray and Kinnmark 1987; Kjaran et al. 1988). Wave model The wave generation model is based on equations developed by the US Army Coastal Engineering Re-

search Center (Shore Protection Manual 1984). We distinguish between four wind directions, southwest, northwest, southeast, and northeast. For each direction the lake is divided into 15–20 sub-areas. Each subarea is assigned a value for fetch and average water depth over the fetch, both of which are estimated from methods given in the Shore Protection Manual (1984).

149 The wave height and period of the significant wave are then given by gF 1/3

0.03790( 2 ) gT s UA = 7.54 tanh[γ2 ] tanh[ γ tanh[ UA 2] gH s

= 0.283 tanh[γ1 ] tanh[

UA2

γ1 = 0.530( γ2 = 0.833(

gd UA2 gd UA2

]

0.00565( gF2 )1/2 UA

tanh[γ1 ]

(8)

] (9)

)3/4 (10) )3/8

Explanation: F Fetch, m D Average depth over the fetch, m Hs Significant wave height, m Ts Significant wave period, s UA Wind stress factor, m/s The wind stress factor relates to the wind speed measured at two meters height above the surface by UA = 0.71(1.25 U2 )1.23 = 0.93 U21.23

(11)

The above equations apply to fetch limited shallow water waves. For typical wave conditions at Lake Myvatn, the duration necessary for waves to become fetch limited is less than one hour. Thus, the assumption of fetch-limited waves is reasonable. As we are more interested in the average wave conditions rather than the conditions described by the significant wave, which relates to the highest one third of the waves, we compute the mean wave height and period from H = 0.626 Hs T = Ts based on information given in the Shore Protection Manual (1984). The bottom shear stress from the wave motions is estimated from linear wave theory. The wavelength L is computed from the linear dispersion relationship where ω2 = gk tanhkh ω = 2π/T k = 2π/L and h is the local water depth. The wave motions at the bottom are sinusoidal with excursion amplitude H /2 sinhkh and maximum velocity Ab =

Ubm

2π = Ab T

The maximum bottom shear stress is now given by ρ 2 τmax = fw Ubm (14) 2 with the corresponding average shear stress (averaged over one wave period) τmax τw = (15) 2 The wave friction factor fw is computed from the following relations (Delo and Ockenden, 1992): fw = max{fws , fwr } where




fws =  fwr =

2 Re−0.5 Re
5 × 105 −0.187 0.0521 Re Re > 5 × 105 0.3 Ab −0.19 ) )

0.00251 e(5.21( ks

Ab ks Ab ks




1.57 > 1.57

(16)  (17)

Here Ab Ubm (18) ν where Re is the wave Reynolds number and ks is a bottom roughness coefficient. The combined shear stress from the wave motions and the currents is estimated using the following relationship (Delo and Ockenden 1992): Re =

τb = τc + τw +ρ Bwc (

fc fw 0.5 ) Ubm |U | 2

(19)

where Bwc = 0.36 as suggested by Delo and Ockenden (1992) if no further measurements are available. Here |U| is the current speed and fc is the current friction factor defined as: 8g fc = 2 (20) C Typically the shear stress due to the waves is an order of magnitude larger than the shear stress due to the currents. However, the currents may be important due to the interaction term in Equation 19, which may be of the same order of magnitude as the wave shear stress term.

(12)

Sediment transport model

(13)

The sediment transport model is based on the conservation of mass for suspended sediments. Other means of sediment transport such as bedload are not

150

Figure 3. Measured wind speed at Geiteyjarstrond, May 30–June 25, 1992.

considered. Thus, the model applies only to the finegrained sediment (Jóhannesson and Birkisson 1989; Jóhannesson and Birkisson 1991; Jóhannesson 1991a) that is transported through suspension. As measurements at Lake Myvatn indicate a well-mixed water column with little concentration change over the water depth, a depth-averaged transport model is well justified. The resulting conservation equation is ∂ ∂ ∂c ∂c ∂ (H D x )+ (H D y )− (H cu) = ∂x ∂x ∂y ∂y ∂x (21) ∂ (H c) + S − Qco ∂t or alternatively, if we substitute in from the continuity equation, ∂ ∂ ∂c ∂c ∂c (H D x )+ (H D y ) − Hu = ∂x ∂x ∂y ∂y ∂x (22) ∂c H + S − Q(co −c) ∂t Explanation: c Concentration of suspended sediment u Velocity within each element taken from the solution of the flow problem, m/s Dx Longitudinal dispersion coefficient, m2 /s Dy Transversal dispersion coefficient, m2 /s H Total water depth, m Q Injected water, m3 /s co Concentration of suspended sediment in the injected water

S Flux of sediment between the lake bottom and the water column The above equations apply to a local coordinate system having the x-axis along the flow direction, so that v = 0. In the finite element approximation, such a local system can be introduced within each element at any given time. The flux of sediment between the lake bottom and water column may be separated into two parts. The deposition of sediment on the bed is expressed by     τb τb
τd  ws c 1 − (23) τd S1 =  τb > τd 0 where ws is the setting velocity, considered to be constant, τd is the critical shear stress for deposition and τb is the combined bottom shear stress from the waves and currents. The erosion of sediment from the bed is expressed by    0 n τb < τe  τb (24) S2 =  τb  τe E −1 τe where τe is the critical shear stress for erosion and E is the erosion rate coefficient. Both parameters are variable within the lake. The above expressions have been used by a number of authors including Delo and Ockenden (1992), Teisson (1991), Luettich et al. (1990), Mehta (1988), and Warren et al. (1992). The above equations are solved numerically using a Galerkin finite element method based on triangular elements.

151

Figure 4. Measured wind direction at Geiteyjarstrond, May 30–June 25, 1992.

Figure 5. Measured wind speed at Geiteyjarstrond, August 21–September 15, 1992.

Measurement data Wind measurements were taken at Geiteyjarstrond (Figure 1) at a height of two meters above ground level. Measurements of waterlevels, current velocity, and suspended sediment concentration were taken at a variety of locations throughout the lake (Figure 2). Although older measurements were used (Vatnaskil 1991), the main emphasis was put on calibrating against the measurements made in the summer of 1992. Within the summer, particular attention was

given to four strong wind ‘events’, June 3–6 and June 15–20 (Figures 3 and 4) with strong southwesterly winds and August 25–28 and September 5–12 (Figures 5 and 6) with strong northwesterly to northeasterly winds. Bathymetry measurements were used to enter water depth into the model (Figure 7). Measurements relating to the bottom sediments were taken from Jóhannesson and Birkisson (1989), Jóhannesson and Birkisson (1991), Jóhannesson (1991a).

152

Figure 6. Measured wind direction at Geiteyjarstrond, August 21–September 15, 1992.

Results Flow model The summertime circulation in the lake is mainly wind-driven, but is also affected by the inflow and outflow to and from the lake. Two parameters may be adjusted to calibrate the flow model, the wind shear stress parameter and the Chezy bottom friction coefficient. The wind shear stress parameter is relatively well known from measurements. It was calculated using Equation 6 by setting the wind drag coefficient, CD = 0.0013, the density of air, ρa = 1.2 kg/m3 , and the fluid density, ρ = 1000 kg/m3 . The resulting value of the wind shear stress parameter was k = 2.469 · 10−6 . Equation 6 is based on wind measured at a height of 10 meters, but here we have adjusted it to correspond to wind measurements taken at a height of 2 meters. The Chezy bottom friction coefficient controls the strength of the bottom friction and is the most important parameter in the calibration process. Although this parameter may be quite variable within the lake, there is no information available to quantify this variation. Therefore, we have assumed a constant value over the entire lake. The final calibration value was C = 60 m1/2 /s. Both of the above mentioned parameters agree with parameters used elsewhere (Shanahan et al. 1986). Figure 8 shows measured and calculated waterlevel fluctuations at Alftagerdi (vhm040). Although

the model does not simulate short-term fluctuations well, the overall trend shows good agreement with the measurements. Measured and calculated waterlevels at other locations in the lake are shown in Vatnaskil (1993, 1998, 2000). Figures 9 and 10 show comparison of measured and calculated current speed and direction at western Teigasund (V-T) during a period of southwesterly winds, which is the most common wind direction at Lake Myvatn during strong wind events. The current speed matches reasonably well where measurements are available and the direction agrees very well with the model reproducing frequent changes of direction. During northerly wind events, the agreement is somewhat poorer. Measured and calculated current speed and direction at other locations in the lake are shown in Vatnaskil (1993, 1998, 2000). The overall agreement between the model and the measurements can be considered quite good. The model reproduces quite well the main features of the wind-driven circulation in the lake, although some of the smaller scale features and localized effects are harder to reproduce. Still, some of the finer scale features, such as the rapid oscillations in the flow at western Teigasund during southwesterly winds are well reproduced by the model. The most important simplification in the model is the assumption that the wind measurements at Geiteyjastrond are representative for the entire lake. This assumption is likely to cause discrepancies between the model and the local measurements, as the wind conditions can be quite variable around the lake.

153

Figure 7. Water depth (m) in Lake Myvatn. Islands are shown in black.

Wave model No special attempt was made to calibrate the wave model, as no recent wave measurements are available for Lake Myvatn. However, some wave measurements were done in the summer of 1990. These measurements were compared with wave forecasts based on the same equations as the current wave model (Jóhannesson 1991a). The agreement was good, especially for the wave height, which is the most important parameter. The only parameter that has to be defined for the wave model is the bottom roughness coefficient, ks . From the value of the Chezy bottom friction coeffi-

cient in the flow model, we can estimate the bottom roughness coefficient for typical conditions at Lake Myvatn. For C = 60 m1/2 /s, which is the value used in the flow model, we get ks = 0.01 m. This value is used to compute the wave friction factor in the wave model. Sediment transport model In calibrating the sediment transport model, most emphasis is put on the measurements of sediment in suspension, with less emphasis on sediment trap results (see Jóhannesson and Birkisson 1991). The reason for this is the well-known fact that sediment traps

154

Figure 8. Measured and calculated waterlevel at Alftagerdi (vhm040).

Figure 9. Measured and calculated current speed at western Teigasund (V-T).

tend to overestimate the true deposition when there are large variations in the setting velocity, whereas the measurements of sediment in suspension are relatively simple and accurate (Kozerski 1994). The sediment measurements consist of two welldefined events, June 3–6 and June 15–20. A third event on August 25–28 is not as well defined as the other two. Thus, the first two events were used to calibrate the model, while the third event was used as a check to verify the calibration. Also, the wind direction in the first two events was southwesterly, whereas in the third event it was northwesterly. The southwesterly wind direction is by far the most important for strong wind events at Lake Myvatn (Jóhannesson 1991b).

The parameters available to calibrate the sediment transport model are the dispersion coefficients, Dx and Dy , the settling velocity, ws , the critical shear stresses for deposition and erosion, τd and τe , the erosion rate coefficient, E, and the power n in the erosion term. As very little information is available to quantify the areal distribution of these parameters (their variation within the lake) most of them are taken as constant over the entire lake. However, the parameters controlling erosion, τe and E, were varied within the lake as shown in Figures 11 and 12. This areal variation is based on information about the lake bottom as well as common knowledge, such as the fact that the erosion strength in coastal areas may be expected to be much larger

155

Figure 10. Measured and calculated current direction at western Teigasund (V-T). Table 1. Parameters for sediment transport model. Parameter

Value

Dx = Dy Ws τd τe – Sydri Floi τe – Bolir τe – Ytri Floi τe – Coastal areas E – Sydri Floi E – Bolir E – Ytri Floi E – Coastal areas n

10 m2 /s 1.0 · 10−4 m/s 0.8 N/m2 0.075 N/m2 0.050 N/m2 0.075 N/m2 0.4 N/m2 0.0055 g/m2 s 0.0115 g/m2 /s 0.0075 g/m2 /s 0.0035 g/m2 /s 0.67

than in deeper areas due to the higher erosion stresses there. However, no quantitative data is available on the relative bottom strength between these areas. Therefore, the values used for these parameters are purely a product of the calibration process. The most important calibration parameters for the sediment transport calculations are the settling velocity and the parameters controlling erosion. The final calibration parameters are given in Table 1. From measurements, Jóhannesson (1991c) suggested the value ws = 1.16 · 10−4 m/s for computations at Lake Myvatn, which is similar to the value used in the model. The measured and calculated suspended sediment concentration at Station B1 are shown in Figure 13. Overall, the fit is reasonable although the model slightly overestimates the higher concen-

trations. Measured and calculated suspended sediment concentration at other locations in the lake are shown in Vatnaskil (1993, 1998, 2000). It should be emphasized here that the data available to calibrate the sediment transport model is limited, therefore affecting the reliability of the calculated results. Still, given the reasonably good agreement between measured and calculated concentration, it can be assumed that the model will simulate the concentration levels in the lake reasonably well throughout the entire summer. The values for the parameters used here all fall within the range of values used by other modelers (see Delo and Ockenden 1992; Luettich et al. 1990; Mehta 1988). The value of critical shear stress for deposition is much higher than that given by Delo and Ockenden (1992) and Mehta (1988), but others have omitted this parameter altogether (setting τd = ∞). Thus, our value lies somewhere between the two extremes. Most modelers have taken n, the power in the erosion term, equal to one. However, during the calibration process it was found necessary to use a smaller value, n = 0.67. This value is within the range used by Luettich et al. (1990). The reason for this lower value is due to a simplification in the modeling of the lake bottom sediment. The model treats the bottom sediment as non-layered with uniform erosion strength. In reality, the bottom sediment is layered with variable erosion strength between the layers. In general, layers with lower erosion strength (less compressed) overlay layers with higher erosion strength (more compressed). A value of n = 0.67 compared with the more commonly

156

Figure 11. Calibrated values of critical shear stress.

used value of n = 1.0 may reflect this simplification. The model uses a lower value of erosion strength that corresponds to the upper, weaker layers of sediment.

Discussion The main purpose of the model was to simulate changes in sediment transport within the lake due to changes in the lake bathymetry caused by proposed diatomite mining of two areas in Bolir, the easternmost

part of the Sydrifloi basin (Figure 14). The sediment deposition into mining Areas 1 and 2 was calculated for the years 1988, 1989, 1990 and 1992 for the three different deepening scenarios (case A, B, and C) shown in Table 2. These four years were chosen because these were the only years for which wind measurements at Geiteyjarstrond existed. The calculated sediment deposition is given in Figure 15 for natural and mining conditions. Sediment deposition increases with increasing mining depth in Areas 1 and

157

Figure 12. Calibrated values of erosion rate coefficient.

2. The maximum increase is in 1992 (about 6,000 tons from the natural deposition) for case C. In addition to Geiteyjarstrond, wind measurements were also available at a nearby station, Reykjahlid (Figure 1). Geiteyjarstrond is an automated monitoring station, and measurements were taken on a 10-minute basis for the years 1988, 1989, 1990, and 1992. Measurements at Reykjahlid exist from 1961 until the present and are taken manually on an 8hour basis. The Geiteyjarstrond measurements were used during the calibration of the model because of their higher frequency. However, since only four years

could be run with these measurements, the question arose as to whether the wind measured during these four years was representative of an extended period of time (several decades). Therefore, the model was run for the time period 1961 to 1997 with the Reykjahlid wind measurements and the total sediment transport in the lake was calculated. The total sediment transport is computed by summarizing all sediment erosion in the lake throughout the summer. The net sedimentation is defined as: ds =D−E dt

158

Figure 13. Measured and calculated suspended sediment concentration at Station B1. Table 2. Deepening of mining scenarios. Deeping

Area 1 m

Area 2 m

Case A Case B Case C

3.2 4 4.5

2 4 6.5

where s = net sedimentation t = time D = deposition E = erosion The total sediment transport, TS is defined as: TS =

T 

E(t) t

lake t =0

where T = summer period This was compared with the total sediment transport calculated for 1988, 1989, 1990, and 1992 using Geiteyjarstrond wind measurements (Figure 16). The comparison shows that the years in question (1988, 1989, 1990 and 1992) have a slightly higher average sediment transport than the entire period. Figure 16 also shows that the Reykjahlid wind produces more extreme values of total sediment transport than the Geiteyjarstrond wind. This is due to the lower frequency of wind measurements at Reykjahlid. The

model assumes a constant wind direction during the time between measurements. Therefore, the Reykjahlid wind direction is constant for 8-h periods while the Geiteyjarstrond wind direction is constant for only 10-min periods. This causes the calculated current velocity for the Reykjahlid wind to remain constant for longer periods which in turn causes the sediment transport calculations to have more variance, giving both higher and lower estimates for the sediment transport than the Geiteyjarstrond wind. The increase in sediment deposition into areas 1 and 2 from the natural conditions is about 2% of the total sediment transport in the lake for the years shown in Figure 15 (the deposition rate into areas 1 and 2 combined is about 2% of the total sediment transport in the lake). It is perhaps more realistic to compare this increase in sediment deposition with the annual accumulation of fresh sediment, however, this value is difficult to calculate. According to an earlier report (‘Committee of Experts’ 1991) this value is approximately 28,000 tons. This means that the increase in sediment deposition due to the deepening of the two mining areas is roughly 20% of the annual accumulation of fresh sediment (dry weight). The flow model may be assumed to be quite accurate, at least with respect to the overall wind-driven circulation in the lake, although smaller scale and localized effects may be harder to simulate. The main cause of discrepancies between the measurements and the calculations is likely to be the assumption of homogeneous wind conditions over the entire lake. A three-dimensional version of the AQUASEA flow

159

Figure 14. Location of proposed mining areas within the lake.

model was developed in order to verify the results of the two-dimensional model, as well as calculate the three-dimensional flow patterns in the lake. The results of the three-dimensional model compared well with the results of the two-dimensional model. Flow patterns were examined for periods of southwesterly winds, the most common direction during strong wind events. In the shallow areas of the lake (near the east and west shores) the current direction is northerly over the entire water column. In the deeper areas of the lake (in the middle) the current direction is not constant

over the entire water column. The current direction is northerly in a thin layer at the top of the water column (upper 30 cm). Below this upper zone, the current direction is southerly. Southwesterly winds cause a build up of waterlevels at the northern edge of Sydrifloi and a lowering of waterlevels at the southern edge, causing gradient induced flow to the south in the deeper areas of the lake. The three-dimensional model was not used to calculate sediment transport, however, it would be expected to give similar results as the two-dimensional model. First, the sediment erosion is controlled mostly

160

Figure 15. Calculated combined sediment deposition (dry weight) in areas 1 and 2.

by wave action, which is similar in both models. Therefore, the amount of sediment erosion should be similar between the models. Second, the reversal of the flow direction in the middle of the lake from south to north is confined to a very thin layer at the top of the water column. The majority of the sediment transport occurs below this upper layer where the flow is to the south (as in the two-dimensional model). Therefore, the sediment transport and deposition should be similar between the models. On the other hand, the dependability and accuracy of the sediment transport model has to be considered much less than that of the flow model. The reasons for this are numerous; some of the more important ones are as follows: 1. In general, the theory behind the sediment transport processes is much less developed than that for the flow processes. Thus, our knowledge of the basic processes such as deposition and erosion is quite incomplete, and mainly based on empirical results. 2. The sediment transport model is based on certain simplifications that may cause considerable errors. An example of this is the assumption of a constant settling velocity. It is well known from measurements that the settling velocity may vary over several orders of magnitude. 3. The data available to verify the sediment transport model is much more limited than that for the flow model. For example, very little is known about the characteristics of the lake bottom, such as the erosion strength and layering.

Therefore, the results of the sediment transport calculations should be viewed with some care. This is particularly true for the quantitative results regarding the sediment fluxes within the lake. However, results regarding the relative change in sediment transport from one situation to another, such as the natural conditions and the situation after mining at Bolir, are expected to be considerably more accurate. The model is least dependable in predicting where the erosion takes place. This is due to the limited knowledge of the characteristics of the lake bottom and their areal variability within the lake. However, if the concentration of sediment in suspension is calculated accurately, as the results of the calibration indicate, the sediment deposition will be well predicted, as it depends mainly on the concentration levels and the settling velocity. Thus, the amount of sediment deposited in the mining areas at Bolir is reasonably well predicted, although the results regarding origin of the sediment and the sediment transport between regions within the lake may be less accurate. There are ways to improve the sediment transport model. First, the current wave model is rather primitive and could be improved by defining the fetch on a finer scale. Second, the model currently assumes a homogeneous, non-layered lake bottom. The model could be improved by adding a multi-layered, heterogeneous lake bottom with varying erosion strength. However, in order to calibrate such a model, more extensive measurements of the lake bottom would be necessary (see Teeter et al. 2001; Madsen et al. 2001).

161

Figure 16. Calculated total sediment transport (dry weight) under natural conditions.

Conclusions Erosion from the lake bottom is a very dynamic process that occurs mainly in shallow areas with relatively loose sediment. As erosion changes the lake bathymetry, erosion areas are shifted from one place to another. Therefore, major changes in sediment erosion at any particular place in the lake are unlikely over an extended period of time, as the lake will seek a new equilibrium. The increase in sediment deposition due to the deepening of the two mining areas is roughly 20% of the annual accumulation of fresh sediment (dry weight). Sediment deposition increases with increasing mining depth in areas 1 and 2. The maximum increase is in 1992 (about 6,000 tons from the natural deposition) for case C. The increase in sediment deposition into areas 1 and 2 from the natural conditions is about 2% of the total sediment transport in the lake for the years shown in Figure 15.

Acknowledgement We would like to thank the Icelandic Research Council for their support through Grant No. 996230000.

References Committee of Experts for Lake Myvatn Research 1991. Effects of the operations of Kísilidjan Inc. on the Lake Myvatn biota. Ministry for the Environment, mimeogr. report. Reykjavik 83 pp.

Delo E.A. and Ockenden M.C. 1992. Estuarine Muds Manual. HR Wallingford Report No. SR 309. Douglas R.W. and Rippey B. 2000. The random redistribution of sediment by wind in a lake. Limnol. Oceanogr. 45: 686–694. Einarsson A., Gardarsson A., Gislason G.M. and Ives A.R. 2002. Consumer-resource interactions and cyclic population dynamics of Tanytarsus gracilentus (Diptera: Chironomidea). J. Anim. Ecol. 71: 832–845. Gardarsson A. 1979. Waterfowl populations of Lake Myvatn and recent changes in numbers and food habits. Oikos 32: 250–270. Gray W.G. and Kinnmark I.P.E. 1986. Evolution of TwoDimensional FE Wave Equation Models, Proceedings of VI International Conference on Finite Elements in Water Resources. Springer Verlag, Berlin, pp. 29–47. Gray W.G. and Kinnmark I.P.E. 1987. Finite Element Surfaces Flow Simulation using Wave Equation. Cunge J.A. and Ackers P. (eds), Proceedings of XIIth Congress IAHR, Technical Session B. EPFL, Lausanne, pp. 113–119. Jóhannesson H. 1991a. Setflutningar í Mývatni, Nidurstödur. Vegagerd Ríkisins (Sediment transport in Lake Myvatn. Results) Icel Road Administration, mimeogr. report, in Icelandic. Jóhannesson H. 1991b. Setflutningar í Mývatni, Mælingar á vindhrada og vindátt. Vegagerd Ríkisins (Sediment transport in Lake Myvatn. Measurements of windspeed and direction). Icel Road Administration, mimeogr. report, in Icelandic. Jóhannesson H. 1991c. Setflutningar í Mývatni, Mælingar á sethrada og kornastærd. Vegagerd Ríkisins (Sediment transport in Lake Myvatn. Measurements of settling velocity and grain size). Icel Road Administration, mimeogr. report, in Icelandic. Jóhannesson H. and Birkisson S.F. 1989. Setflutningar í Mývatni, Yfirbordssetid í Mývatni og botnrof í Ytriflóa. Vegagerd Ríkisins (Sediment transport in Lake Myvatn. Surface layer of sediment in Lake Myvatn and bottom erosion in Ytrifloi). Icel Road Administration, mimeogr. report, in Icelandic. Jóhannesson H. and Birkisson S.F. 1991. Setflutningar í Mývatni, Svifaursmælingar og setgildrur. Vegagerd Ríkisins (Sediment transport in Lake Myvatn. Suspended sediment measurements and sediment traps). Icel Road Administration, mimeogr. report, in Icelandic.

162 Jónasson P.M. (ed) 1979. Ecology of eutrophic subarctic Lake Myvatn. Oikos 32, 308 pp. Kjaran S.P., Holm S.L. and Sigurdsson S. 1988. Testing of Finite Element Schemes for Linear Shallow Water Equations, Proceedings of VII Internaztional Conference on Computational Methods in Water Resources. Kozerski H.P. 1994. Possibilities and limitations of sediment traps to measure sedimentation and resuspension. Hydrobiologia 284: 93–100. Lindegaard C. and Jónasson P.M. 1979. Abundance, population dynamics and production of zoobenthos in Lake Myvatn, Iceland. Oikos 32: 202–227. Luettich Jr. R.A., Harleman D.R.F. and Somlyódy L. 1990. Dynamic behavior of suspended sediment concentrations in a shallow lake perturbed by episodic wind events. Limnol. Oceanogr. 35: 1050–1067. Madsen J.D., Chambers P.A., James W.F., Koch E.W. and Westlake D.F. 2001. The interaction between water movement, sediment dynamics and submersed macrophytes. Hydrobiologia 444: 71– 84. Mehta A.J. 1988. Laboratory studies on cohesive sediment deposition and erosion. In: Donkers J. and van Leussen W. (eds), Physical Processes in Estuaries. Springer-Verlag, Berlin, pp. 427–445. Shanahan P., Harleman D.R.F. and Somlyódy L. 1986. Windinduced water motion. In: Modelling and Managing Shallow Lake Eutrophication: With application to Lake Balaton. Springer-Verlag, Berlin. Sheng Y.P. and Lick W. 1979. The transport and resuspension of sediments in a shallow lake. J. Geophys. Res. 84 (NC4): 1809– 1826.

Teeter A.M., Johnson B.H., Berger C., Stelling G., Scheffner N.W., Garcia M.H. and Parchure T.M. 2001. Hydrodynamic and sediment transport modeling with emphasis on shallowwater, vegetated areas (lakes, reservoirs, estuaries and lagoons). Hydrobiologia 444: 1–23. Teisson C. 1991. Cohesive suspended sediment transport: Feasibility and limitations of numerical modeling. J. Hydraulic Res. 29: 755–769. Tómasson G. and Kjaran S.P. 1993. Setflutningar í Mývatni. (Sediment transport in Lake Myvatn), Árbók VFÍ (VFI Annual Report), Reykjavík: 294–305 US Army Coastal Engineering Research Center 1984. Shore Protection Manual, Vol. 1. Department of the Army, Washington, DC. Vatnaskil 1991. Straumar í Mývatni (Currents in Lake Myvatn). Vatnaskil Consulting Engineers, mimeogr. report, in Icelandic. Vatnaskil 1993 Mývatn – Lake circulation and sediment transport Vatnaskil Consulting Engineers, mimeogr. report, in Icelandic. Vatnaskil 1998 Mývatn – Lake circulation and sediment transport Vatnaskil Consulting Engineers, mimeogr. report, in Icelandic. Vatnaskil 2000 Mývatn – Lake circulation and sediment transport Vatnaskil Consulting Engineers, mimeogr. report, in Icelandic Warren R., Jensen O.K., Johnsen J. and Johnson H.K. 1992. Modelling of the morphological processes in Venice lagoon. Proceedings of the Second International Conference on Hydraulic and Environmental Modelling of Coastal, Estuarine and River Waters. Volume 2. Ashgate.