Hydrodynamic and sediment transport modelling in the canals of Venice (Italy)

Estuarine, Coastal and Shelf Science 75 (2007) 250e260 www.elsevier.com/locate/ecss

Hydrodynamic and sediment transport modelling in the canals of Venice (Italy) Elisa Coraci, Georg Umgiesser*, Roberto Zonta CNR-ISMAR, Consiglio Nazionale delle Ricerche, Istituto di Scienze Marine, San Polo, 1364, 30125 Venezia, Italy Received 15 October 2006; accepted 8 February 2007 Available online 9 July 2007

Abstract A framework of numerical models was applied to the canal network of Venice in order to simulate hydrodynamics and sediment transport. A two-dimensional finite element model (SHYFEM) was applied to the whole surrounding lagoon using water levels as a forcing. The computed levels were extracted along the contour of the city of Venice and then a link-node model was run in the canal network. The simulated variables were calibrated and compared with data from field measurements. Calculated water elevation displays a good agreement with the measured data, and current velocity is well reproduced. Simulations were initially carried out using a constant Manning parameter (MC) and successively a varying coefficient (VMC) was adopted, in order to account for the different depths of the canals. By applying the sediment transport module an evaluation of the sediment accumulation rate in the canal network was obtained, permitting an estimation of the yearly silting up of the canal bed. Different simulations were carried out considering material input from both the city and the lagoon. The model is able to reproduce the overall accretion rate of the canal bottom, and it is a useful tool for planning the dredging activities in the whole canal network. Ó 2007 Elsevier Ltd. All rights reserved. Keywords: Venice Lagoon; Venice canals; numerical modelling; sediment transport; silting

1. Introduction Venice Lagoon is the largest lagoon in the Mediterranean Sea. It covers a wide coastal area parallel to the sea in the NEeSW direction for about 50 km between the Piave River mouth and the delta of the Po River. The lagoon is separated from the sea by sandy barrier littorals, interrupted by three inlets (Lido, Malamocco and Chioggia) through which the tidal wave enters the lagoon (Fig. 1). Lagoons developing in sedimentary coastal zones are often unstable. If not destroyed by the sea, they tend to be filled up by sediment coming from the sea or supplied by rivers of the drainage basin. Geo-morphological variations have been naturally occurring in the Venice Lagoon since its formation. In the past the large sediment load delivered from the drainage basin threatened to transform the lagoon into a marshland. However, * Corresponding author. E-mail address: [email protected] (G. Umgiesser). 0272-7714/$ - see front matter Ó 2007 Elsevier Ltd. All rights reserved. doi:10.1016/j.ecss.2007.02.028

the diversion of the main freshwater tributaries out of the lagoon, as well as the increase in the lagoon depth, due to subsidence and eustatic water level rise (Gatto and Carbognin, 1981), has reversed the natural tendency of the lagoon to silt up, transforming it slowly into a more sea-like environment. The city of Venice is located in the central part of the lagoon and is built on more than a hundred small islands, bordered by an intricate network of almost 160 interlinked canals (known as ‘rii’) having width from a few to tens of meters, and a depth of 1e5 m. The Grand Canal is the main watercourse of this intricate network, shaping out a large ‘‘S’’: from there the majority of the other canals spread out in all directions until reaching the surrounding lagoon. The overall length of the canal system is about 40 km and its surface corresponds to, approximately, 10% of the total urban area (Zonta et al., 2005a). The hydrodynamic regime in the canals is driven by tidal forcing, and as a consequence of phase lags and level gradients occurring at the city boundary. During the flood tide the flow is

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

251

Fig. 1. The Venice Lagoon with the SHYFEM finite element grid.

predominantly from SE to NW, while the direction is reversed in the ebb. Current speed in the network is low, with average maximum values of up to 25 cm s1 throughout the tidal cycle. Besides constituting the principal way of transportation for human and public services, the canal network is the collector of all types of urban refuse from untreated domestic and commercial sewage effluents. Venice, in fact, does not have a real sewage system. Due to the sluggish circulation, the delivered organic and inorganic matter, as well as materials eroded from the urban surfaces, reach the network, allowing for a progressive silting up of the canal bottom. After a long period (1960e1995) where no dredging occurred, the Municipality of Venice started an overall series of dredging activities to permit navigation and maintain a healthy environment. Today, a general improvement of the sewage system has commenced including septic tanks and small wastewater treatment plants, as well as the restoration of building foundations and consolidation of canal walls (Gardin, 1999). As a consequence, a reduction of the loading of particulate, organic matter and contaminants to the canal network is expected to occur. Due to a growing interest towards the ecological problems of the canals, various research projects regarding behaviour and quality of the network were performed in recent years (Collavini et al., 2004; Zonta et al., 2005b, 2006; D’Este et al., 2006). Concerning hydrodynamics and modelling, De Marchi (1993) first applied a simple link-node model to the canal

network. Successively a new link-node one-dimensional model (DYNHYD), that is the hydrodynamic module of the water quality model WASP (Ambrose et al., 1993), was employed and calibrated with field measurements (Zampato and Umgiesser, 2000; Umgiesser and Zampato, 2001). This canal model was coupled with a finite element model of the lagoon (SHYFEM). The further step was to modify the link-node model from 1D to 2D structure in order to describe the sediment transport in the canal network (Umgiesser and Massalin, 2000). This paper presents the results of the recent modelling activity carried out in the canal network, in the framework of the project QUEST (QUality, Efficiency, Sedimentation and Transport in the Venice canal network). The aforementioned 2D canal model was improved creating a new link-node grid that has a greater accuracy in reproducing the canal shape. A new data set (water level and speed) was acquired to calibrate the model and to obtain a more complete investigation on the hydrodynamics of the canal system. Moreover, the sediment transport module was implemented in the erosion equation and tested, allowing for a comparison with field data derived from the bathymetric measurements. The sediment dynamic was finally investigated to quantify the material flux to and from the canal bottom. A first attempt at estimating the contribution of different sediment sources, both inside and outside the city, was made in order to evaluate erosion and deposition rates in the canals. A further objective of this study was to create a useful tool for planning the dredging activities in the whole canal network.

252

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

2. Methods 2.1. Field data For the purpose of calibrating the model of the canal network, both hydrodynamic and sedimentological data were needed. Time series of hydrodynamic data were acquired in the period from October to December 2004, in four zones of the city: San Polo, Santa Marta, Cannaregio, and Sant’Elena (Fig. 2). Three canals were chosen in every zone, and self-recording electromagnetic current meters (S4, InterOcean, USA) were deployed for a time period of 15 days in each of the canals to record water level and current speed. The sampling frequency of the data was 5 min. Data acquired were processed using a 1-h moving average to smooth peaks generated by the boat traffic. To get a measured value of the canal bottom, the bathymetric data acquired by Insula (the company which is in charge of the interventions on the canal network) were used. The depth variation was calculated from two field measurements carried out in a time period which varied from two to five years. The available data cover 50% of the canal network extent. The mean value of sediment density, 1540 kg m3, was derived from the available measured data regarding the sediment composition in the canal network of Venice (Venice Municipality database).

2.2. The hydrodynamic models The lagoon and the canal network are two environments with different morphological characteristics: the mean depth of the lagoon is about 1 m while the canal depth ranges about from 1 to 5 m; moreover, the canals have small width in comparison with their length. As a consequence, two different types of models were used.

2.2.1. The Venice Lagoon model For the lagoon, the 2D finite element shallow-water model SHYFEM (Umgiesser and Bergamasco, 1993; Umgiesser et al., 2004) was used. The finite element method permits to describe with a great detail the morphology and bathymetry of the area, and to represent with higher resolution the zones where hydrodynamic features are more interesting. The numerical computation was carried out on a spatial domain that represents the entire lagoon through a finite element grid (Fig. 1). The latter contains 4237 nodes and 7666 triangular elements. The model considers the three lagoon inlets as open boundaries, while the rest of the perimeter of the lagoon as a closed boundary. The model resolves, with a semi-implicit time stepping scheme, the vertically integrated shallow-water equations in their formulations with level and transports: vU vh þ gH þ RU þ X ¼ 0 vt vx

ð1Þ

vV vh þ gH þ RV þ Y ¼ 0 vt vy

ð2Þ

vh vU vV þ þ ¼ 0; vt vx vy

ð3Þ

where h [m] is the water level, g is the gravitational acceleration [m s2], H ¼ h þ h the total water depth, h is the undisturbed water depth, t is the time and R [s1] is the friction term that is expressed through a Manning formulation (see below). U and V [m2 s1] are the vertical integrated velocities (total or barotropic transports) defined as: Z n Z n U¼ u dz; V ¼ v dz; ð4Þ h

h

where u and v [m s1] are the velocities in x and y directions. The terms X, Y in Eqs. (1) and (2) contain all other terms like

Fig. 2. The city of Venice with the link-node grid. In black are shown the 12 canals where the hydrodynamic data were measured. The four zones are pointed out by letters: A e San Polo; B e Santa Marta; C e Cannaregio; and D e Sant’Elena.

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

the Coriolis terms, the wind stress, the non-linear convective terms and Reynolds stresses that need not to be treated implicitly in the time discretization. At open boundaries the water levels are prescribed. The SHYFEM model has been calibrated and reproduces the tidal oscillation well in most parts of the lagoon (Umgiesser et al., 2004). 2.2.2. The link-node model of the Venice canal network The model of the canal network has been developed starting from the one employed in a previous study about the circulation and sediment transport in the inner canals (Umgiesser and Massalin, 2000). It is a shallow-water 2D (i.e. along longitudinal canal axis and vertical directions) link-node model. The hydrodynamic module of the model computes water surface elevation and horizontal water velocity in the canal, whereas a sediment transport module calculates current bottom stress, erosionedeposition strengths, sediment concentration and bottom thickness. A new grid of the model was set up in order to reproduce accurately the canal shapes (Fig. 2). It is made out of 505 canal elements (links) and 431 nodes that connect the canals; each node is set in the intersection between two or more canals or in the points where there is an angle. Each canal is split into one or more links depending on its length and shape, therefore the number of links is higher than the number of canals. It has been supposed that each link is sufficiently represented by a constant width and consequently Coriolis acceleration is negligible, as well as other accelerations perpendicular to the flow direction. The model has a vertical structure implemented as s-coordinates that are based on a linear transformation function as described in Burchard and Petersen (1997). However, since the model has been used only in a one-layer approximation, the barotropic equations are reported here. The link-node model solves the layer-integrated equations where the horizontal diffusion terms have been neglected (Burchard and Baumert, 1998). All the equations are listed below. The continuity equation is: vh vðHuÞ þ ¼ 0; vt vx

ð5Þ

where H [m] is the total height of the water column and u [m s1] is the horizontal velocity (x direction). The momentum equation (horizontal advection and pressure gradient) reads: vðHuÞ vðuHuÞ þ  t[z ðuÞ þ tYz ðuÞ ¼ Px ; vt vx

ð6Þ

where tz is the vertical diffusive fluxes described below. The arrows [ and Y indicate that the quantity is calculated at the surface and the bottom of the water column, respectively. The depth-integrated pressure force term Px has the form:

Px ¼ gH

vh gH 2 vr þ : vx 2r0 vx

253

ð7Þ

The first term on the right-hand side represents the external part (barotropic term) of the pressure gradient; the second term is the internal part (baroclinic term) that results from integrating the horizontal varying density over the total water depth. Here r [kg m3] is the density and r0 [kg m3] is a reference density. The term tz ðf Þ describes the vertical diffusive flux of the quantity f: on f depends the meaning of nt (vertical turbulent diffusion coefficient) and nm (kinematics viscosity if f hu or molecular diffusivity if f is concentration or salinity). The salinity S [psu] transport equation is: vðHSÞ vðuHSÞ þ ¼ 0: ð8Þ vt vx The simplified equation of state that gives the density in function of the salinity for a constant temperature value of 10  C is:   ð9Þ r ¼ r0 1 þ aS þ bS2 ; where a and b are numerical coefficients (a ¼ 7:776  104 and b ¼ 3:351  109 ). The suspended particulate matter (SPM) transport equation reads: vðHCÞ vðuHCÞ þ þ Fe  Fd ¼ 0; vt vx

ð10Þ

where C [kg m3] is the SPM concentration. The erosion and deposition equation is: vzb 1 ¼ ðFd  Fe Þ; vt rs where (c

ð11Þ

e

ðjtb j  tce Þ if tb > tce rs ( 0 elsewhere ws Cb ðtcd  jtb jÞ if tb < tcd and Cb > 0 Fd ¼ tcd 0 elsewhere Fe ¼

ð12Þ

are the erosion and deposition strengths [kg m2 s1], respectively, and ws [m s1] is the constant settling velocity. In these equations zb [m] is the sediment bed thickness, rs [kg m3] is the sediment density, Cb [kg m3] is the near-bottom SPM concentration, tb [kg m1 s2] is the bottom shear stress, tcd and tce [kg m1 s2] are the critical stress threshold for deposition and erosion, respectively, and ce [kg s m4] is the proportionality factor. A new equation, proposed by Li and Amos (2001), of the critical stress for erosion was employed: therefore the erosion threshold value is not constant but increases as sediment is eroded away. Its variation with the down core depth z is calculated by: tce ðzÞ ¼ tce ð0Þ þ Aðrs  rÞgz tan fi ;

ð13Þ

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

where tce ð0Þ is the critical erosion stress at the sediment surface, A is an empirical coefficient for down core sediment resistance and is equal to 0.01, and fi [degrees] is the internal friction angle of cohesive sediment. The further sediments are eroded, the more difficult it becomes to continue erosion. The new deposits are considered as superficial sediment without consolidation, so their erosion threshold is the critical erosion at the sediment surface tce ð0Þ. The bottom shear stress is calculated using the following parameterization: tb ¼ rCD juju;

ð14Þ

where CD is the drag coefficient given by (barotropic assumption): CD ¼

gMC2 ; H 4=3

ð15Þ

where MC [m1/3 s] is the Manning coefficient (Zampato and Umgiesser, 2000). 2.2.3. Implementation of the models The two hydrodynamic models are coupled in a way that the lagoon model provides the boundary conditions for the canal model. The coupling occurs on the boundary of the city in the points where the canals flow into the lagoon: in such points the surface elevation values calculated by the lagoon model are imposed for the whole duration of the simulation, providing realistic boundary conditions. The latter include the physics of tidal propagation in the lagoon such as phase shift, tidal amplification and the effect of wind stress. SPM concentration values have to be imposed at the lagoonecanal interface when the external load has to be simulated. In the SHYFEM simulation the Manning coefficient was used in the friction term of the momentum equations: the coefficient values varied from 0.023 to 0.031 for the lagoon, depending on the characteristics of the area, e.g. canals, mud flat and vegetation. The drag coefficient for the momentum transfer of wind was set to 2:5  103 . All the simulations were run with a time step of 5 min. The lagoon model reproduces accurate water levels, and there was no need for further calibration. The time step of integration for the canal network model was set to 3 s. The wind forcing was neglected in the city canals and it is only taken into account through the simulated water levels of the lagoon model since many buildings of Venice hinder the wind action on the canal flow. All the simulations of the canal model were run with only one vertical layer, since the SHYFEM is a barotropic model. The calibration of the link-node model was performed varying the MC, as described in the next section.

measurements were carried out. In the second part, results concerning the erosion and sedimentation rates are described. 3.1. Hydrodynamic calibration The water level in the canals was calculated and compared with field data. For instance, Fig. 3 shows the case of the Rio delle Becarie (A zone, San Polo). Water elevation is reproduced with high accuracy by the model and this result is representative of the other 11 canals investigated. This good agreement with the experimental data confirms that the canal model does not require further calibration with respect to the water level. As far as the current velocities are concerned, the canal model was calibrated by tuning the MC (Eq. (15)) around the 0.035 value, which was found to be the best fit in a previous study (Zampato and Umgiesser, 2000). Five coefficients were selected ranging from 0.025 to 0.045, with an increment of 0.005, in order to carry out the simulations with the same values for all canals; the last simulation, instead, was done with a Manning coefficient (VMC) varying with the water depth H, as proposed by De Marchi (1997): values ranged from 0.015 (H > 3 m) to 0.065 (H < 0.9 m). In order to determine the model goodness-of-fit on the basis of mean, variance and correlation, a comparison was performed by using a linear regression analysis of the prediction error, i.e. the difference between matched pairs of observed and predicted values (Smith and Rose, 1995). The comparison was based on the decomposing of the sum of squared prediction errors (SSDs) into three components: 2

SSD ¼ Nd 2 þ Nðsobs  smod Þ þ 2Nð1  rÞsobs smod ;

ð16Þ

corresponding to accuracy (bias), precision (variance differences), and degree of association (correlation difference from 1.0), respectively. N is the number of matched pairs, d is the mean of the di values (i.e. the error in prediction), defined as di ¼ yobs  ymod , i ¼ 1.N; sobs and smod are the stani i dard deviations of values yobs and ymod, and r is the correlation coefficient between the two data sets. 1.2 data model

1.0 0.8

Liv [m]

254

0.6 0.4 0.2 0.0 -0.2

3. Results and discussion The results of the simulations in the canal network are presented. The first part of the discussion involves water elevation and current velocity in the 12 canals where the field

-0.4

0

24

48

72

96

120 144 168 192 216 240 264 288

Time [h] Fig. 3. Simulation from 7 to19 October 2004 in the Rio delle Becarie: water elevation computed by the model and measured.

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

The statistical results are represented in a 2D plot on which the correlation coefficient and the root-mean-square difference between the two data sets (observed and calculated) are indicated by a single point (Taylor, 2001). This type of diagram summarizes the degree of correspondence between the field data and the model. In fact, Taylor proposed a similar method to quantify the differences in two data sets using the rootmean-square difference E, defined as: sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi PN i 2 i i¼1 ðyobs  ymod Þ : E¼ N

ð17Þ

The quantity E can be resolved into two components to isolate differences in phase and amplitude from differences in the means of the two fields: 2

E2 ¼ d 2 þ E0 ;

ð18Þ

where d is the bias, and E0 is computed by standard deviations and the correlation coefficient, which are the second and third terms in Eq. (16) Therefore Eqs. (16) and (18) are equivalent with the exception of a factor N. Taylor diagrams of the four investigated zones of the canal network are reported in Fig. 4, showing the statistical results of the simulations. Six points (simulations) are shown for each zone, corresponding to five different MC values (from 0.025 to 0.045, step 0.005), and the case with the variable coefficient VMC. The diagrams refer to current speed only, since the value of MC has a negligible effect on the water level. The standard deviation is reported along the two axes, so the radial distance from the origin is proportional to the standard deviation of a single simulation. The point representing the observed data is plotted on the abscissa and the position of other points (i.e. simulations) depends on their standard deviation and correlation coefficient with respect to the field data. The correlation is given by the azimuthal position. The distance between the simulated and field data is the error term E0 and is due to the differences in amplitude and phase. Referring to Eq. (16) the diagram gives information about the second and third terms but nothing about the bias; the best simulations display points are as close as possible to the data point on the abscissa. In order to compare the different simulations, standard deviations have been normalized by the standard deviation of the field data. In the first zone (A San Polo), the current velocity of all the three canals is well reproduced using the MC value of 0.040. Here simulated and observed standard deviations are equal. Moreover, the correlation coefficient is quite high, about 0.9 for the three canals; the time series of current velocity for Rio delle Becarie is shown in Fig. 5. In the simulation with MC equal to 0.040 the velocity is reproduced with greater accuracy, particularly in spring tide conditions, even if the model overestimates the amplitude of the data. In the second zone (B Santa Marta) the simulations carried out with MC equal to 0.045 show the best agreement with the data collected from two canals (Rio dell’Arzere and Rio di

255

Santa Margherita) in which the standard deviations of data and model are equal. In the third canal (Rio del Tentor), however, the model overestimates the current values. Here with the VMC the standard deviation is 1.5, which is closest to that of the measured data. The lowest correlation, around 0.35, is found in the two canals (Rio dell’Arzere and Rio del Tentor) characterized by slow currents (up to 10 cm/s). The three canals in the third zone (C Cannaregio) have some peculiar aspects and are quite different from each other (Fig. 6). The Rio di San Giobbe shows low unidirectional current: it does not invert the direction during the ebb and flood tides (Carrera, 1999). The model was not able to reproduce this atypical behaviour, which occurs also in other canals of the network. The canal Rio San Felice is characterized by higher current velocities than the other two canals, with maxima of up to 50 cm s1. In the Taylor diagram for the two canals with bi-directional flow, the current velocities are reproduced with smaller error than in the Rio di San Giobbe, especially in the Rio di San Felice where the correlation is particularly high (0.9). The fourth zone investigated (D Sant’Elena) is located in the south-eastern part of the city, near the Lido inlet, where the tidal wave has a small lag between the southern and northern boundaries, and consequently the elevation gradient is low. This zone was chosen in order to test the model under hydrodynamic conditions which are difficult to reproduce. The model was able to reproduce the current amplitude in Rio della Tana and Rio di Sant’Isepo using MC values equal to 0.25 and 0.40, respectively. In the third canal (Rio dei Giardini) the correlation is low due to the phase shift between the two time series. Based on this analysis, one characteristic observed in the statistical parameters is that the most important contribution in the error estimation is due to E0 . In fact the component E is always smaller because the average values of measured and predicted current velocities frequently have the same value. For example, the values of E0 are between 0.04 and 0.05 m s1 in the San Polo and Santa Marta zones but they increase up to 0.12 m s1 in the other two zones investigated where there are canals with unidirectional flow (Cannaregio) and the circulation induced by the tidal force is difficult to reproduce (Sant’Elena), due to the geographical position. Selecting one single Manning coefficient suitable for all these canals is not easy. On the other hand, notable results were not obtained using a VMC to take into account the different depth values in the canals. However, data show that the minimum error is obtained with a MC of 0.035 (Fig. 7), therefore this value was chosen for the further simulations of the sediment transport in the canal network. 3.2. Sediment transport The link-node model calculates the erosion and deposition rates generated by the stress exerted by currents only. Factors such as boat traffic or other anthropogenic influences cannot be described. All the simulations were run using the settling velocity ws equal to 102 cm s1, value in the range proposed

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

256

Fig. 4. Taylor diagram of current velocity in the four investigated zones: A e San Polo; B e Santa Marta; C e Cannaregio; and D e Sant’Elena. The Manning coefficients are indicated by the following letters: a ¼ 0.025; b ¼ 0.030; c ¼ 0.035; d ¼ 0.040; e ¼ 0.045; and v ¼ variable.

in Amos et al. (2004), and the sediment density rs equal to 1540 kg m3. The critical erosion threshold is calculated as in Eq. (13) with tce ð0Þ ¼ 0:30 [kg m1 s2] and with the internal friction angle fi ¼ 20 ; the critical deposition threshold has the value of 0.25. The two thresholds have been chosen based on a study carried out in the Venice Lagoon (Amos et al., 2004). All simulations are extrapolated to all of 2004; input water level values around the city were obtained from a SHYFEM

simulation. This simulation utilized measured tidal and wind data taken from the Lido inlet. The lagoon plays an important role in the sedimentation process occurring in the canal network since the material is brought into the city by the tidal flux and settles there, especially during slack tide. At the same time, organic and inorganic matter directly enters the canals from the urban effluents. Another input of matter is represented by the masonry

0.4

0.2 data model

0.2

Velocity [m/s]

Vel [m/s]

0.1

0.0

-0.1

-0.2

0.0 -0.2 San Giobbe Misericordia San Felice

-0.4 -0.6 0

24

48

72

96

120 144 168 192 216 240 264 288

Time [h] Fig. 5. Simulation from 7 to 19 October 2004 in the Rio delle Becarie: current velocity computed by the model and measured.

0

12

24

36

48

60

72

84

96

Time [h] Fig. 6. Current velocity measured from 9 to 12 November 2004 in the three canals of the Cannaregio zone.

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

257

Sediment thickness [mm]

35 1A 1B 2A 2B

30 25 20 15 10 5 0 -5 0

30

60

90

120 150 180 210 240 270 300 330 360

Time [d]

Fig. 7. Taylor diagram of current velocity in all the 12 canals investigated.

Fig. 8. Average sediment thickness calculated in the simulations for all the year 2004 with different inputs: 1A ¼ with lagoon (23 mg/l) and urban load (2217 þ 1120 m3/y); 1B ¼ only with lagoon (23 mg/l); 2A ¼ with lagoon (13 mg/l) and urban load (3350 þ 1700 m3/y); 2B ¼ only with lagoon (13 mg/l).

debris originating from erosion of building surface due to weathering and capillary movement of salt water through the structure. Due to the absence of measured values related to SPM concentration in the water entering the city, yearly sewage load, and yearly load due to the erosion of the urban surfaces, some reasonable values must be chosen to test the model. The estimated value for the sewage input is 2217 m3 y1 (Carrera et al., 1999). This load was calculated using water consumption of the Venice inhabitants during the years 1997e1998 and estimating a density value of 300 g m3 for solids in the wastewater. Moreover, the contribution from the building erosion was calculated to be about 1120 m3 y1, considering the average value of different building deterioration rates. Two different simulations were carried out. In the first simulation the value of 23 mg l1 (Umgiesser and Massalin, 2000) for the suspended particles concentration at the boundary nodes of the canal network was imposed in the flood tide to simulate the lagoon input. This SPM value was kept constant during all the simulation due to the unavailability of time series of concentration at the boundary nodes. This assumption is supported by the presence of strong currents along the city boundary that carry away the suspended matter outgoing from the city. Moreover, the SPM could settle almost completely in about 3 h (i.e. half an ebb tide) with the used settling velocity (0.01 cm s-1), thus the possibility of the SPM re-entering might be low. In the second simulation the urban load (i.e. sewage and surfaces erosion) was added in all nodes at every time step to represent the real situation. Fig. 8 shows the average sediment thickness in the two scenarios, obtained from the ratio between the total volume of the sediment at the bottom and the total canal network area. The major contribution in the deposition is due to the input coming from the lagoon: with only this load the average sediment

thickness in the city is 2.6 cm y1. In the simulation with the three sediment sources the mean thickness is 3.1 cm y1. The difference between the two lines in the plot represents the contribution of the urban load which is about 0.5 cm y1. So the results of these simulations show that the more important input for the sedimentation in the city comes from the lagoon, which is responsible for 85% of the overall silting up. However, results of a sedimentation study in a test canal (Zonta et al., 2005b) and data calculated from depth variations measured in the 250 links furnish a sedimentation rate (ca. 2 cm y1) which is lower than the obtained 3.1 cm y1. The model was therefore applied to estimate the sediment bottom thickness using different load values, which were considered more realistic with respect to the previous ones. For the suspended particles at the boundary nodes a concentration of 13 mg l1 was used, which is in agreement with average values measured in some monitoring activities (Alberotanza and Zucchetta, 1992; Scattolin et al., 2003; Zonta et al., 2006). Concerning the urban inputs, load values proposed by Carrera et al. (1999) were multiplied by a factor of 1.5. The values 3350 m3 y1 and 1700 m3 y1 were used for the sewage and erosion loads, respectively, in order to take into account both the large daily presence of tourists and commuters that almost doubles the city population (Zonta et al., 2005a), and other matter sources such as food refuse from domestic sinks and runoff from roofs and paved areas, which were not considered in the previous load estimate. The results of simulations performed with the new load values are also shown in Fig. 8. The obtained sedimentation rates are 1.2 cm y1 only when the lagoon contribution is considered, and 1.8 cm y1 when including the urban inputs. The urban loads produce a sedimentation rate of 0.6 cm y1, i.e. 50% with respect to the lagoon input. This latter value is reduced by 40% with respect to the first simulation, while the lagoon input still has the main role in the silting up of the canal network.

258

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

The new estimated values of the sedimentation rates are in a better agreement with data resulting from past studies and field measurements. Fig. 9 shows a map of the erosion and deposition rates predicted in the whole canal network, which are referred to the simulation carried out with the second set of load values. Generally the canals are characterized by the sedimentation rate between 1 and 3 cm y1, whereas only in 20% of the canal segments erosion is observed. Erosion occurs in canals characterized by high current speed, such as in the Cannaregio zone, as well as to the east of the Grand Canal, where strong NeS level gradients occur at the boundaries nodes. The canals with the highest sedimentation rate (about 5 cm y1) are located mainly in the centralewestern zone of the city (San Polo and Santa Marta zones). A sedimentation rate of about 3 cm y1 finally characterizes the canals in the Sant’Elena zone. Since the results depend on the imposed SPM concentration value, another set of simulations was carried out to calculate the range of the sedimentation rate. The available measured data of the suspended matter were analyzed to find their variability: the average value is 13 mg l1 (used in the last simulations) and the standard deviation is 5 mg l1. Thus the new simulations were performed using the minimum (8 mg l1) and maximum (18 mg l1) values for the lagoon input only. The other simulations were carried out taking into account the internal input also: in three simulations the internal loads were kept constant at the mean value used before (i.e. 3350 m3 y1 for sewage and 1700 m3 y1 for erosion) while the minimum, mean and maximum values of the lagoon input were imposed. Afterwards it was assumed that the internal loads also have a similar variability as the lagoon input. In the two last simulations therefore two input values were calculated from the mean loads of sewage and erosion contributions and they are 2080 and 4600 m3 y1 for the sewage effluents, 1050 and 2340 m3 y1 for erosion and runoff from roofs and

Table 1 Sedimentation rates [cm y1] calculated by the model with different inputs, depending on the imposed concentration coming from the lagoon. The columns show the results obtained without urban input, with the mean urban input and with the variable urban input

Min. Mean Max.

No urban input

Mean urban input

Variable urban input

0.4 1.2 1.9

1.1 1.8 2.5

0.8 1.8 2.7

paved areas. The average sedimentation rates calculated in the last set of simulations are presented in Table 1. The results obtained without the urban input show a sedimentation range between 0.4 and 1.9 cm y1 depending on the variable lagoon suspended concentration. The same variability is obtained when also the mean urban input was added. Finally considering all the particulate inputs that vary between the minimum and maximum values, the estimated rates are between 0.8 and 2.7 cm y1, respectively. The results show that with minimum load lagoon and internal loading account for approximately 50% each (0.4 cm y1); on the other, the lagoon is responsible for 70% of the overall silting up when the input is incremented up to the maximum value. The lagoon plays the main role in the sedimentation process since the variation in the lagoon input produces a variation in the deposition rate of about 0.7 cm y1, while the difference due to the urban load variability is smaller, about 0.2 cm y1.

4. Conclusions An existing shallow-water finite element model of the Venice Lagoon (SHYFEM) was used to reproduce the tidal elevation around the city of Venice. The simulated data were used

Fig. 9. Map of the erosion and deposition rates calculated by the model in the 505 links of the grid.

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

as the boundary condition in a link-node model of the Venice canal network. A good reproduction of the tidal wave propagation in the Venice Lagoon by the finite element model is very important in order to simulate the hydrodynamic features in the canal network. Calculated water elevations show a good agreement with the measured data collected in the years 2003e2004 in four zones of Venice. As far as the current velocities are concerned, a calibration was performed using the Manning coefficient as the tuning parameter. Results show that the coefficient providing the best agreement with the data has different values in the investigated canals. The value of 0.035 was taken as the best value for all the canals, but for instance a lower Manning coefficient produced the best results for canals with current velocity higher than 20 cm s1. Major differences between measured and simulated data were observed in canals located in the south-eastern zone of the city, which are a consequence of the low level gradient generally occurring between the boundaries. In this case the resulting water circulation pattern is difficult to be reproduced with the model. The canal network collects particulate matter entering from the lagoon with the tidal currents. Another main contribution of particulate matter is directly delivered from the domestic and commercial sewage effluents and from the erosion of building surfaces and paved areas. These material fluxes increase the natural trend of the canals to accumulate sediment transported in the network by tidal currents. The link-node model was run to estimate the average sedimentation rate in the whole network. For this purpose a set of simulations were carried out to calculate the single contribution of the three aforementioned main sediment sources. Results show that the major contribution in the sedimentation comes from the lagoon input. The model gives a good description of the sediment dynamics in the Venice canals and it is able to estimate with a reasonable approximation the average sedimentation rate. It is therefore a useful tool in order to distinguish the sediment bottom thickness variation in the different zones of the city, providing a tool in the planning dredging activities. However, results strictly depend on the imposed value of the matter inputs. It is therefore necessary to collect new data, in order to have a better estimation of these inputs. For instance, some measuring activities are presently carried out with the aim of acquiring better estimates of the sewage inputs and sedimentation fluxes, and temporary varying SPM concentration values to impose at the city boundary nodes. Acknowledgements The study was realized in the ambit of the project QUEST (QUality, Efficiency, Sedimentation and Transport in the Venice canal network), financed by Insula, Venice. Authors thank Antonella Zaupa for the bathymetric data and the technicians of CNR-ISMAR (M. Meneghin, R. Ruggeri, F. Simionato and G. Zamperoni) for the precious assistance in the field data acquisition.

259

References Alberotanza, L., Zucchetta, G., 1992. Caratteristiche delle acque della laguna di Venezia, (Bacino centrale). In: Dati relativi al monitoraggio effettuato nel corso dell’anno 1990. CCID, CNR, 50 pp. (in Italian). Ambrose, R.B., Wool, T.A., Martin, J.L., 1993. The Water Quality Analysis Simulation Program, WASP5. Part A. Environmental Research Laboratory, Athens, Georgia. Amos, C.L., Bergamasco, A., Umgiesser, G., Cappucci, S., Cloutier, D., Denat, L., Flindt, M., Bonardi, M., Cristante, S., 2004. The stability of tidal flats in Venice Lagoon e the results of in-situ measurements using two benthic, annular flumes. Journal of Marine Systems 51, 211e241. Burchard, H., Petersen, O., 1997. Hybridization between s- and z- co-ordinates for improving the internal pressure gradient calculation in marine models with steep bottom slopes. International Journal for Numerical Methods in Fluids 25, 1003e1023. Burchard, H., Baumert, H., 1998. The formation of estuarine turbidity maxima due to density effects in the salt wedge. A Hydrodynamic Process Study. Journal of Physical Oceanography 28, 309e321. Carrera, F., 1999. Il comportamento idrodinamico dei canali interni di Venezia. In: Venezia la citta` dei rii. UNESCO/INSULA, pp. 197e207 (in Italian). Carrera, F., Borrelli, A., Horstick, J., 1999. Quantification of Sediment Sources in the City of Venice, Italy. Project Number: FC-VEI2. Faculty of Worcester Polytechnic Institute, Massachusetts, USA. Collavini, F., Coraci, E., Zaggia, L., Zaupa, A., Zonta, R., Zuliani, A., 2004. Environmental characterisation of the Venice canals: first results of the ICARO project. In: Proceedings of the 37th Congress C.I.E.S.M., Barcelona, Spain, 7e11 June, 182 p. De Marchi, C., 1993. A Hydrodynamic Model of the Network of Canals of the City of Venice. Master Thesis, University of California, Davis, CA, USA. De Marchi, C., 1997. Un modello matematico del regime idraulico della rete dei canali interni della citta` di Venezia. UNESCO project ‘‘I canali interni di Venezia’’. Draft Report 2, 72 pp. (in Italian). D’Este, C., Zonta, R., Collavini, F., 2006. Water quality of the Venice canals (Italy): results of a first monitoring activity. In: Proceedings of the ECSA 41st International Conference, Venice, 15e20 October, 93 p. Gardin, P., 1999. La ripresa dello scavo dei rii. In: Venezia la citta` dei rii. UNESCO/INSULA, pp. 91e105 (in Italian). Gatto, P., Carbognin, L., 1981. The Lagoon of Venice: natural environmental trend and man-induced modification. Hydrological Sciences Bulletin 26 (4), 379e391. Li, M.Z., Amos, C.L., 2001. SEDTRANS96: the upgraded and better calibrated sediment transport model for continental shelves. Computers and Geosciences 27, 619e645. Scattolin, M., Zonta, R., Botter, M., Tagliapietra, D., Rismondo, A., 2003. Indagine sulla qualita` delle aree marginali di Venezia. Rapporto Tecnico, Comune di Venezia Eds., 80 pp. (in Italian). Smith, E.P., Rose, K.A., 1995. Model goodness-of-fit analysis using regression and related techniques. Ecological Modelling 77, 49e64. Taylor, K.E., 2001. Summarizing multiple aspects of model performance in a single diagram. Journal of Geophysical Research 106 (D7), 7183e 7192. Umgiesser, G., Bergamasco, A., 1993. A Staggered Grid Finite Element Model of the Venice Lagoon. In: Morgan, K., Ofiate, E., Periaux, J., Zienkiewicz, O.C. (Eds.), Finite Elements in Fluids. Pineridge Press, pp. 659e668. Umgiesser, G., Massalin, A., 2000. Modeling the Sediment Transport in the Channels of Venice. Final Report of the Project ‘‘Development of a Sediment Transport Model of the Inner Canals of Venice’’. UNESCO Venice Office, Venezia, Italy, 54 pp. Umgiesser, G., Zampato, L., 2001. Hydrodynamic and salinity modeling of the Venice channel network with coupled 1De2D mathematical models. Ecological Modelling 138 (1e3), 75e85. Umgiesser, G., Melaku Canu, D., Cucco, A., Solidoro, C., 2004. A finite element model for the Venice Lagoon. Development, set up, calibration and validation. Journal of Marine Systems 51 (1e4), 123e145.

260

E. Coraci et al. / Estuarine, Coastal and Shelf Science 75 (2007) 250e260

Zampato, L., Umgiesser, G., 2000. Hydrodynamic modeling in the channel network of Venice. Il Nuovo Cimento C 23 (3), 263e284. Zonta, R., Zaggia, L., Collavini, F., Costa, F., Scattolin, M., 2005a. Sediment contamination assessment of the Venice canal network (Italy). In: Fletcher, C.A., Spencer, T. (Eds.), Flooding and Environmental Challenges for Venice and Its Lagoon: State of Knowledge. Cambridge University Press, pp. 603e615.

Zonta, R., Cochran, J.K., Collavini, F., Costa, F., Scattolin, M., Zaggia, L., 2005b. Sediment and heavy metal accumulation in a small canal, Venice, Italy. Aquatic Ecosystem Health & Management 8 (1), 63e71. Zonta, R., Zuliani, A., Coraci, E., D’Este, C., Pesce, A., 2006. Fluxes of particulate, metals and nutrients in a test canal of Venice (Italy). In: Proceedings of the ECSA 41st International Conference, Venice, 15e20 October, 106 p.