Abstract

Southern Iberian Peninsula coastal waters have a high ecologic and economic value. The objective of this work has consisted of describing numerical models which simulate the dispersion of particle-reactive tracers in these waters and reproducing measured concentration levels in water and bed sediment samples. Additionally, information on suspended matter concentration distributions and sedimentation rates is obtained. In the Alboran Sea, the model has been applied to simulate the transport of radionuclides introduced from atmospheric fallout. Transport pathways of heavy metals discharged from three rivers draining a large mining area have been obtained for the Gulf of Cadiz. It has been found that these rivers constitute a source of trace elements to the Mediterranean through the Strait of Gibraltar. Some characteristic times have been calculated, as well as fluxes of isotopes through the Strait of Gibraltar, as additional results of environmental interest.

1. Introduction

The Gulf of Cadiz (GoC) is responsible for 5–10% of fish and shell-fish catches of Spain and Portugal [1], holding important living resources of commercial and ecological interest. Similarly, the Alboran Sea (AS) is one of the most productive areas in the Mediterranean Sea [2]. This whole area (Figure 1) has a high ecological and economic value as it constitutes the only connection between the Mediterranean Sea and the Atlantic Ocean. The region is also exposed to potential pollutant discharges given the intense shipping activities through the Strait of Gibraltar, connecting the GoC and AS, and the usual adverse meteorological conditions here. We must also consider the presence of large industrial areas in the Spanish shores. Finally, the Guadalquivir, Guadiana, and Odiel-Tinto rivers, in the southern Iberian Peninsula (Figure 2(a)), present strongly enhanced heavy metal concentrations since they drain the Iberian Pyrite Belt [3], one of the most important mining areas in the south of Europe. Mineral resources have been extracted in the last 5000 years during two main periods: the Roman age and the last two centuries. During the last period, intensive exploitation has led to a relevant environmental impact, with vast surfaces covered with mining residues and subjected to erosion [3].

Consequently, it is relevant to study and understand geochemistry and dispersion patterns of contaminants in the system, since this will help assessing the potential influence of pollutants on ecosystem functioning. Thus the main objective of the present work has consisted of reviewing appropriate numerical models developed in the University of Seville, which simulate the dispersion of particle-reactive tracers in southern Iberian Peninsula waters and reproducing measured concentration levels in water and bed sediment samples [4, 5]. As an additional result, information on suspended matter distributions and sedimentation rates is obtained. Valuable environmental results concerning tracer fluxes through the Strait of Gibraltar and characteristic times in the Alboran Sea, not published before, are also obtained carrying out numerical experiments. This highlights the usefulness of this type of dispersion models for reactive pollutants.

Dispersion models are based upon appropriate hydrodynamic models which provide water circulation in each system (GoC and AS). Tracers which can be easily identified in each system are used for the study. Thus, heavy metals (Zn and Cu) have been simulated in the GoC. These metals are introduced by the three rivers mentioned above. Fallout radionuclides (Cs, Pu, and Pb) are considered for the AS.

Ideally, for a proper testing and validation of a dispersion model, we would need a well-known point source of pollutants. In the GoC, metal inputs from rivers may be used. In the AS we had to use fallout radionuclides due to the lack of such a well-defined source. Also ideally, we would need an spatial array of measurements to appreciate if the spatial distribution of pollutants is reproduced by the model. These requirements are satisfied in the GoC since there are metal measurements in water and in bed sediments, along the Spanish coast, where a clear spatial structure appears, but data is more limited in the AS. The situation is even more complicated since time variations will be apparent as well, more in the dissolved phase than in bed sediments. Only in limited cases there is enough information to allow a complete model validation. This is the case, for instance, of authorized discharges from nuclear fuel reprocessing plants to the sea, where wide monitoring programs are undertaken and the source term is well known [68]. But in practice, it is rather difficult that these requirements are satisfied and model results and validations have to be interpreted with care.

Published models describing trace metal dispersion in the GoC consider metals as conservative tracers, no interacting with sediments and without any other sources and sinks [1, 9]. Similarly, models describing trace elements dispersion in the AS, including water/sediment interactions, are not found in literature. These reactive tracers may provide information not only on water circulation, but on sediment processes as well.

Some background oceanography is provided in Section 2. After, models are described (Section 3). Next, results are presented and discussed for the AS and the GoC in Sections  4.1 and 4.2, respectively. Models have also been applied to evaluate radionuclide fluxes through the Strait of Gibraltar and characteristic times of the AS system. Some sensitivity analysis are finally presented.

2. Background Physical Oceanography

2.1. Alboran Sea

The water circulation in the Strait of Gibraltar is characterized by a surface inflow of Atlantic water and a deep outflow of more dense Mediterranean water. Exchanged flows are of the order of 1 Sv with a net inflow into the Mediterranean Sea of about 0.05 Sv [10]. This net inflow compensates the excess of evaporation over precipitation and river supply in the Mediterranean. The exchanged flows show some variability. In particular, Tsimplis and Bryden [10] gave an estimation of  Sv for the Atlantic inflow, where the large error is due to seasonal variability, and  Sv for the outflow. Nevertheless, these values are about 10% lower than the lowest previous estimates [11]. A review on transport estimations through the Strait of Gibraltar may be seen in such reference: from the very first carried out in 1877 and the widely cited value of 1.2 Sv for both outflow and inflow [12], to their own estimations of 0.72 and 0.68 Sv for inflow and outflow, respectively.

The water masses that may identified in the Alboran Sea are: the Surface Atlantic Water (SAW) that has its origin in the North Atlantic Central Water (NACW) that has been modified by air-sea interactions, thus strictly speaking is not a water mass [13]; the Surface Mediterranean Water (SMW) which occupies the surface layer in areas not reached by the SAW; the Levantine Intermediate Water (LIW) that flows from the Mediterranean to the Atlantic at depths between some 200 and 600 m and, finally, the Western Mediterranean Deep Water (WMDW), below the LIW and also flowing towards the Atlantic.

However, this description may be simplified to a two-layer system with water flowing in opposite directions [1416]. Thus, from an operative point of view, SAW and MSW constitute the upper layer, and LIW and WMDW constitute the bottom layer, flowing to the west. Indeed, a steady westward flow of LIW and WMDW has been recorded in several experiments [14, 17]. This approach has already been used to study the water exchanges between the Atlantic and the Mediterranean by means of numerical models [1821].

The surface circulation in the Alboran Sea is as follows. The inflowing surface Atlantic water penetrates the Strait of Gibraltar (see Figure 2(b) for geographic names) and reaches mean velocities of the order of 0.6 m/s. This water forms a jet that enters the Alboran Sea in a east-northeast direction. The jet flows along the Spanish coast and curves to the south at about longitude. Then part of it flows to the west, incorporating to an anticyclonic gyre, while the remainder flows to the African coast between Cape Tres Forcas and Alboran Isle [22]. The anticyclonic gyre is known as the Western Alboran Gyre (WAG). A detailed description may be seen for instance in [23, 24] and references cited in these papers, which include the classical works carried out since the 70s. The WAG is an almost permanent feature, although presents variations in intensity and even there are periods in which disappears [23]. The spatial scale of the WAG is about 100 km in diameter and some 200 m in depth.

Bryden and Stommel [25] have found that the outflow of water from the Mediterranean to the Atlantic occurs along the Moroccan continental slope, and that this circulation appears to be a nearly permanent feature of the southern side of the Alboran Sea.

With respect to the tides in the area, an important feature of the tidal flow in the Strait is that it can be considered, as a first approach, as barotropic. Indeed, 93% of the variance of current velocity in the semidiurnal band has a barotropic character in the Strait [26]. Tsimplis and Bryden [10] have pointed out that tidal currents are barotropic and larger than the mean inflow or outflow. The semidiurnal tide dominates ADCP records in the Strait, obscuring the expected two-layer character of the mean flow. The tidal signal is so strong that it reverses the currents near the bottom for a part of each tidal cycle. As a consequence, 2D depth-averaged models have already been applied to simulate surface tides in the Strait [27]. Tsimplis et al. [28] have even used a 2D barotropic model for simulating tides in the whole Mediterranean Sea.

In the case of the main tidal constituent, , tide amplitude at the Atlantic entrance of the Strait is about 0.8 m. This amplitude decreases towards the east. Thus, at the Mediterranean entrance is only about 0.3 m. The amplitude of the tide is reduced more in the Alboran Sea, reaching 0.09 m at the east limit of our domain (see map in Figure 2(b)). The associated currents decrease in a similar way: from tidal currents of the order of 1 m/s in the Strait to currents of a few cm/s in the Alboran Sea basin. Similar behaviour is observed for the constituent.

Turning now to sedimentology, it has been reported [29] that essentially all the western Alboran basin and southern part of the eastern one (which indeed is the whole domain of the model) are covered with muds, with less than 7% of coarse material (>63 μm). The composition of these muds is controlled by river inputs, marine carbonate production in shallow areas, biogenic production in surface water, and water mass circulation [30]. Lithogenic particles are mostly introduced to the Alboran Sea by fluvial and eolian mechanisms and, to a small extent, by coastal erosion (less than 5% of the total; see [31]). Eolian supply is strongly controlled by the proximity of the Sahara Desert. Thus, some episodes of atmospheric dust deposition, mostly of Saharan origin, have been recorded [31]. The fluvial input of particles has been extensively investigated by Liquete et al. [32]. These authors have found that supply from the north coast is at least one order of magnitude larger than from the southern one. They provide detailed data on sediment supply by the main rivers discharging into the Alboran Sea.

2.2. Gulf of Cadiz

Recently, the GoC has been the subject of oceanographic studies dealing with surface and deep circulation, aimed at understanding the mechanisms of water exchange between the Atlantic Ocean and the Mediterranean Sea; as well as the behaviour of the dense plume of Mediterranean water [13, 3335]. Figure 3 shows a TS diagram obtained after a cruise in the GoC [13]. The different waters found in the samplings can be listed as follows. (1)North Atlantic Central Water (NACW): the linear behaviour in the TS diagram, characteristic of NACW, is found for the following and ranges: C; (although its most superficial region will be called SAW henceforth). Below a certain depth (from a certain isopycne), the TS diagram diverges from its former linear behaviour due to mixing with the underlying, salty Mediterranean Water (MW).(2)Surface Atlantic Water (SAW): this water has its origin in shallow NACW that has been modified by air sea interactions so, strictly speaking, SAW is not a water mass. In the Gulf of Cadiz, SAW can be characterized for practical purposes as a region of the TS diagram with temperatures above 16.0°C and salinity , which is found between the surface and a depth of approximately 100 m over the entire region except for the continental shelf.(3)Warm Shelf Waters (SW): at some of the sampled stations, mainly over the continental shelf, water warmer and fresher than SAW has been detected at the surface. This water comes from SAW that has been noticeably influenced by continental shelf processes, including heating and fresh water inputs from land. It corresponds to the points of the TS diagram between temperatures of 14.0°C and 18.0°C.(4)Mediterranean Water (MW): the analysis of the thermohaline properties of MW shows two main maxima in the TS diagram corresponding to two cores (upper and lower) with different densities, centered near 800 and 1200 m.(5)North Atlantic Deep Water (NADW): at a reduced number of deep stations (depth m) in the southwest part of the surveyed region, NADW, which is characterised by its depth-decreasing thermohaline properties [36], has been found.

The surface circulation of the Gulf of Cadiz, less studied than its deep circulation, is integrated into the general circulation of the Northeast Atlantic: the Azores Current, which transports some 15 Sv between latitudes 35.1°N and 40.1°N to feed the Canary Current, frequently forms meanders that separate themselves from the main flow [37]. Thus, the surface circulation of the Gulf of Cadiz could be understood as the last meander of the said Azores Current. Part of this meander enters the Mediterranean Sea through the Strait of Gibraltar to balance evaporation and buoyancy losses within this Sea.

The residual surface circulation in the northern GoC is characterized by a current directed to the SE [13, 35] along the Spanish coast. This circulation is a rather constant pattern during most of the year. Part of the flow enters the Strait of Gibraltar and part is deflected to the south. Below the surface, the Mediterranean waters flow into the Atlantic and mainly direct to the NW [38, 39].

The amplitude of the tide is about 1 m over all the GoC, decreasing near the Strait of Gibraltar. Associated currents are weak, with amplitudes below 0.10 m/s over most of the GoC. Indeed at the RAP (Red de Aguas Profundas, Spanish Institute of Oceanography) buoy position (see Figure 16), the barotropic tidal current is less than 0.03 m/s [33]. The current amplitude increases as approaching the Strait entrance, where currents about 0.8 m/s are produced. A similar behaviour is observed for the tide.

With respect to sediments, schematically, it may be said that sand and sandy mud dominate to a depth of 30 m. The middle shelf, to a depth of 100 m, is characterized by an extensive mud belt. Finally, outer sediments below 100 m are dominated by mixtures of sand and clay [40].

3. Model Descriptions

The tracer transport model is based upon an advection/diffusion equation with terms describing the exchanges of contaminants between the dissolved and the solid phases (suspended matter particles and bed sediments). Suspended matter concentrations and sedimentation rates are consequently required. They are obtained from a sediment transport model developed for this purpose. Water currents in the area are required to solve the advective parts of the equations and to calculate bed stress, which is needed to evaluate deposition and resuspension of particles from the seabed. Two computational domains were defined, Figures 2(a) and 2(b), for the GoC and AS, respectively.

3.1. Hydrodynamics

The water current at any position is obtained through the addition of the instantaneous barotropic tidal current plus a residual (mean or long-term) circulation. Both in the GoC and the AS, tides are computed using a 2D depth-averaged model, which is a reasonable approach [41, 42] already used successfully in the Strait of Gibraltar [43, 44]. Tidal equations are standard and may be seen, for instance, in [45]. The solution of these equations provides the water currents at each point in the model domain and for each time step. Currents are treated through standard tidal analysis [46] and tidal constants are stored in files that will be read by the transport models. The two main tidal constituents, and , are considered. Thus, the hydrodynamic equations are solved for each constituent and tidal analysis is also carried out for each constituent separately. Tidal constants allow a very fast calculation of the tidal current at any time and point in the domain.

From an operative point of view, circulation in the Strait of Gibraltar and the AS may be simplified to a two-layer system with water flowing in opposite directions, as commented above. Consequently, a 2-layer model has been adopted to calculate the residual circulation in the AS. Equations may be seen in Izquierdo et al. [18].

Complex mixing processes between several water masses occur in the GoC [13], thus a full 3D primitive-equation baroclinic hydrodynamic model is used. It is based upon the hydrostatic and Bousssinesq approximations on a plane. The model includes equations for salinity and temperature evolution and water density is calculated from them using a standard state equation. A one-equation turbulence model has been used to calculate the vertical eddy viscosity. Details on the 3D equations may be seen, for instance, in Kowalick and Murty book [50]. A summary of the main characteristics of hydrodynamic models may be seen in Table 1 and detailed descriptions and equations may also be seen in [44]. Residual circulations in the GoC and AS are again stored in files which are appropriately read by the dispersion codes.

It must be noted that the simplest approach (2-layer for AS, 3D for GoC) capable of reproducing the most relevant features of circulation in each area has been adopted.

3.2. Sediment Transport

The transport of sediments is described by an advection-diffusion equation to which some terms are added. These are external sources of particles (river supply), terms describing particle deposition on the seabed and resuspension from the bed to the water column, and vertical settling. The formulation of these processes is based upon standard formulae. Thus, the erodability constant is used for the resuspension term. Particle settling and deposition are described using the settling velocity, which is obtained from Stoke’s law. Critical resuspension and deposition stresses are applied as usually. Details on the mathematical formulation may be seen elsewhere [5154]. Finally, it is also possible to calculate net sedimentation rates (SR) as the balance between the deposition and resuspension terms. All equations and terms are appropriately written both for the 2-layer and the full 3D approaches. They are presented in Appendix A, for illustrative purposes, in the full 3D form.

A number of approximations had to be assumed in both domains, mainly because the lack of field data to feed the model and to validate it. Thus, episodes of heavy rains, atmospheric deposition, and bedload transport are not included. Moreover, only the lithogenic particle fraction, that is about 80% of the total suspended load [30, 31], is simulated. Finally, atmospheric deposition events of particles coming from the Sahara Desert have not been considered in the model since they cannot be easily quantified [31].

3.3. Tracer Transport

Nonconservative pollutants are those which do not remain dissolved in the water column, but have a certain affinity to be fixed to particles (suspended matter and bed sediments) and/or those which may suffer any decay process (radioactive, biodegradation, and so forth). Thus, three phases are included in the model: dissolved, suspended particles and bed sediments. Only fine sediments are considered (particles with diameter < 62.5 μm) since metals and radionuclides are predominantly fixed to these [55]. The exchanges between the dissolved and solid phases may be described in terms of kinetic transfer coefficients. Thus, assuming that adsorption/release reactions are governed by a single reversible reaction, a coefficient governs the transfer from the liquid to the solid phase and a coefficient governs the inverse process. Also, the migration of tracers to the deep sediment is included. Thus, tracers deposited on the sediment surface will be buried by particle deposition and will migrate below the mixed sediment layer which directly interacts with the dissolved phase. This effect may be easily treated as a decay process with constant given by where is the sediment mixing depth (the distance to which the dissolved phase penetrates the sediment), is the sediment bulk density (dry mass divided by wet volume), and SR is the sedimentation rate provided by the sediment transport submodel. Note that the model implicitly assumes an instantaneous pollutant mixing over this sediment layer with thickness . A summary of all the processes involved may be seen in Figure 4, in the case of a 2-layered sea. There is advection and diffusion of suspended particles and contaminants in both water layers, as well as uptake/release reactions between dissolved pollutants and suspended sediments. There will also be some diffusion through the pycnocline, direct uptake/release reactions between the bottom water layer and bed sediments, and deposition of suspended particles on the seabed.

It is known that adsorption depends on the surface of particles per water volume unit at each point and time. This quantity has been denoted as the exchange surface [6, 7, 45, 56]. Thus, the kinetic coefficient is written as where and are the exchange surfaces for suspended matter and bottom sediments respectively (dimensions ) and is a parameter with the dimensions of a velocity denoted as the exchange velocity (see references cited above). The delta function is introduced to take into account that only the deepest water layer interacts with the bed sediment. Thus for the deepest layer and elsewhere.

Assuming spherical particles, the exchange surfaces are written as (again see references cited above) where and are particle radius and density, respectively, is the suspended matter concentration, is the fraction of fine particles in the sediment, is sediment porosity, and is a correction factor that takes into account that part of the sediment particle surface may be hidden by other sediment particles. Finally, is the thickness of the deepest water layer. This formulation has been successfully used in all modelling works cited above. Real particles are not spheres, but with this approach it is possible to obtain an analytical expression for the exchange surface [57]. The kinetic coefficient is considered to be constant.

Tracer dispersion equations are summarized in Appendix B, where they are presented in a general 3D form. They are appropriately modified in the case of the 2-layer AS model.

3.4. Numerical Solution

The topography of both domains was obtained from the NOAA (US National Ocean and Atmosphere Administration) GEODAS database, available on-line, with a resolution of two minutes in both longitude and latitude. Fifty vertical levels are considered in the GoC 3D model. Topography of both domains is presented in Figures 2(a) and 2(b). All the equations were solved using explicit finite difference schemes. Details on the solution of the hydrodynamic equations are given in [44] and are not repeated here. It is worth commenting that the simple addition of tidal and long-term currents implies that nonlinear interactions between both flow components are neglected. However, this is a common practice in pollutant transport modelling (see for instance [58]).

The second-order accuracy MSOU (Monotonic Second Order Upstream) scheme [59] is applied for all advective terms and a second-order scheme is used for diffusion as well [50]. Mass conservations has been careffully checked for all numerical schemes.

Time-step in the different submodels is selected with care to assure that stability conditions are satisfied. A time-splitting procedure was used to treat the water-sediment pollutant interactions. This is essentially the same as in hydrodynamic models, when time-splitting is used to solve the external and internal modes. The procedure is required because a very restrictive stability condition is introduced by such interaction terms [45]. Essentially, it consists of using a smaller time step to solve water-sediment tracer interactions. Thus, for each advection/diffusion time step, a number of smaller steps are needed to solve tracer partition.

Appropriate boundary conditions must be defined as well. They are described below. Finally, it must be commented that both tidal and residual currents are used in the sediment transport models, in order to describe appropriately deposition and resuspension events. However, only residual currents and mean SRs provided by the sediment models are used to simulate tracer dispersion, given the generally small tidal currents in the domains.

3.5. Parameters

The horizontal diffusion coefficient, 2.0 m2/s in both domains, has been selected according to the horizontal grid resolution [44]. The vertical diffusion coefficient in the 3D model is set equal to the vertical viscosity provided by the one-equation turbulence model, as usual in most models.

Suspended particles (diameter < 63 μm) are characterized by a mean size that is considered to be representative of suspended matter in the considered environment. The particle size controls, through Stoke’s law, the settling velocity of particles and, as a consequence, affects the sedimentation rate as well. It also affects the adsorption of tracers from the dissolved phase ((3) and (4)). The value used for particle diameter was 8 m for both the GoC and AS. Freitas and Abrantes [60] have found that particles with diameter <10 μm are dominant in all water masses in the Gulf of Cadiz. The GoC reflects the AS since surface waters in the first flow into the second and Mediterranean waters flow, in the deep layer, from the AS into the GoC. A sensitivity analysis (not shown) has also been carried out, trying different particle sizes. The best results were obtained with a diameter of 8 m.

The critical deposition stress typically ranges between 0.04 and 0.1 N/m2 [61], while the critical resuspension stress ranges 0.1–1.5 N/m2 [61]. In the present applications intermediate values of 0.06 and 1.0 N/m2 have been taken for the critical deposition and resuspension stresses, respectively. The erodability constant is fixed as  kg/m2 s. This parameter typically varies between and kg/m2 s [61]. The fraction of fine particles in bed sediments must be defined over the model domain. In the case of the GoC, this information has been compiled from [40]. As already mentioned, sand and sandy mud dominate to a depth of 30 m. The middle shelf, to a depth of 100 m, is characterized by an extensive mud belt. Finally, outer sediments below 100 m are dominated by mixtures of sand and clay. For the AS, since less than 7% of the bed sediments consists of coarse material [29], it has been taken an uniform value of .

The equilibrium distribution coefficient, , concept is used to estimate the exchange velocity as explained below. The describes the partition of a tracer between the dissolved and solid phases at equilibrium: where and are tracer concentrations in the solid and dissolved phases, respectively. This parameter is of course depending on the tracer geochemical behaviour and on environmental conditions, thus it is a site-specific parameter.

As described before [62], is very similar even for elements with a rather different geochemical behaviour, being (see (2)) the essential parameter describing the tracer geochemical behaviour. Thus, the same value is given to for all simulated tracers, which is  s−1. This value has been successfully used in earlier simulations in the Odiel-Tinto estuary, in the GoC [63]. It was measured for Cs by Nyffeler et al. [62].

The exchange velocity for each tracer is obtained from the equation relating this parameter with and its corresponding equilibrium distribution coefficient, [45]: The mean values recommended by the International Atomic Energy Agency [64] for coastal waters have been used. This procedure has been widely used before [7, 51, 63, 65].

Distribution coefficients may vary in more than one order of magnitude depending on environmental conditions. Consequently, the model sensitivity to this parameter has been studied and some results are presented below (Section  4.2.4). For some tracers, the had to be obtained from calibration if there is not information about it or if results with the recommended value are not satisfactory. This may occur due to environmental variability.

Values for other parameters as , and have been obtained from previous modelling works. A summary of model parameters and their sources may be seen in Table 2.

4. Results and Discussion

Results are presented for the AS and for the GoC separately. In both cases results of the hydrodynamic model are initially discussed, later sediment transport results are presented. For the AS, dispersion patterns of Cs and Pu are presented in one Section  4.1.3 and, in a following one (4.1.4), Pb results are shown. Later, fluxes through the Strait of Gibraltar (Section  4.1.5) and characteristic times (Section  4.1.6) are analyzed. In the case of the GoC metal dispersion is studied in Section  4.2.3 and, in the next Section  4.2.4, sensitivity analysis results are discussed.

It should be highlighted that results from hydrodynamic models have been successfully implemented on rapid-response particle tracking models for radioactive, chemical and oil spills already developed for the AS and GoC [44, 6668].

4.1. Alboran Sea
4.1.1. Hydrodynamics

Computed values of the amplitude and phase of the water surface elevation due to the and tides have been compared with observed values [28, 69] at different locations in the Strait of Gibraltar and the Alboran Sea. Results are given in Table 3. Locations may be seen in the map in Figure 5. It can be observed that generally there is a good agreement between the measured and computed tidal constants for both constituents. A tidal chart for the tide is shown in Figure 6. The same color joins points with the same tide amplitude (corange chart). The amplitude decreases quickly in the Strait of Gibraltar, from some 0.75 m in the Atlantic entrance to 0.30 in the section of Gibraltar. There is a further, although slower, amplitude reduction in the Alboran Sea, reaching only some 0.10 m at the eastern boundary. It can be seen that corange lines are oriented in a south-north direction, while cotidal lines (not shown) are essentially in a northeast direction, in agreement with observations and the earlier computations in [28].

The spatial distribution of the current amplitude is presented in Figure 7 as an example, and a comparison of computed current ellipse parameters with those deduced from measurements [70] is presented in Table 4. In general, there is a good agreement between both set of data. Returning to Figure 7, the current amplitude in most of the Alboran Sea is some 0.03 m/s. Only in some specific areas larger currents are obtained. This is the case of the south-north section going from Cape Tres Forcas to Alboran Isle and the Spanish coast. Currents of the order of 0.06 to 0.09 m/s are computed close to the Cape, around Alboran Isle and close to the Spanish coast. This is in good agreement with the earlier computations presented in [28]. Current amplitude increases as approaching Gibraltar Strait, where, as expected, maximum values are obtained. Indeed, currents of the order of 0.8 m/s are computed in some locations. For clarity reasons, however, the scale maximum in Figure 7 is limited to 0.2 m/s. In the case of the tide, results show a similar decrease in the tide amplitude, from some 0.25 m at the Atlantic entrance of the Strait to some 0.04 m at the eastern open boundary of the domain. Current amplitudes decrease from maximum values of the order of 0.25 m/s in the Strait of Gibraltar to currents weaker than 1 cm/s in the Alboran Sea.

The calculated geostrophic flows for the upper and lower layers, for the mean water exchanges through the Strait of Gibraltar, are presented in Figure 8. In the upper layer, the jet of Atlantic water entering through the Strait of Gibraltar flows towards the east along the Spanish coast and partially curves to the south before reaching Alboran Island. Part of this flow continues to the east between Cape Tres Forcas and Alboran Island and the remaining rotates towards the west. A gyre of anticyclonic circulation is thus completed. Surface water velocity in the Strait reaches 0.53 m/s, in agreement with the current speeds of the order of 0.6 m/s reported from measurements [22] and models [71]. In the jet, along the Spanish coast, maximum surface velocities are of the order of 0.4 m/s. Measurements of Perkins et al. [22] in this area range between 0.1 and 0.53 m/s and the model of Werner et al. [72] produces maximum currents in the north of the gyre of 0.25 m/s. Vargas-Yáñez et al. [23] have measured eastward geostrophic velocities in the passage between Cape tres Forcas and Alboran Island (point F in Figure 5). They are , and cm/s at depths of 74, 117, and 178 m, respectively, (all depths correspond to the surface layer). Computed current at the same point for the surface layer is 3.8 cm/s.

A westward circulation is obtained for the bottom layer. In the western Strait of Gibraltar computed velocities are of the order of 0.20 m/s, in good agreement with the 0.15 m/s obtained by Béranguer et al. [73]. Maximum outflow current calculated in the Strait is 0.36 m/s; Sannino et al. [74] obtained a figure of 0.35 m/s. In the Alboran Sea, water velocities are reduced to some 0.1 m/s in the southern area and less in the northern part of the basin. This circulation pattern for the deep layer is in good agreement with the earlier calculations in [20], that also show that most of the deep water flows along the south coast with speed below 0.15 m/s, and with measurements in [14], which give velocities in the southern shelf, near the Strait, of 0.1 m/s. As commented before, Bryden and Stommel [25] have also found, from CTD transects, that the westward deep flow occurs along the Moroccan continental slope. These authors have measured a mean (using a 341-day record) outflow velocity of cm/s in a point located at about longitude and over the 500 m isobath in the Moroccan continental slope. They do not give the exact position of the current meter mooring, but the model predicts (for such longitude and a depth of 535 m) a westward current equal to 3.6 cm/s. In the northern part of the sea, Fabres et al. [31] have measured westward currents of 1.93 and and 0.95 cm/s in stations C and E (see Figure 5), respectively. The corresponding calculated currents are 1.0 and 0.60 cm/s.

As expected for geostrophic flow, the interface between both water layers in the Strait of Gibraltar is deeper in the south than in the north. The computed slope along a south-north section at the Atlantic entrance of the Strait is , in good agreement with [74] who found a value of .

Some nonlinear interaction between tides and the mean flow may be expected, as happens in the model of Nikiema et al. [75] for instance. However, despite such effects, it seems that, generally speaking, the present models give a representation of the Alboran Sea tides and circulation that is realistic enough to implement on it the sediment and pollutant transport models.

4.1.2. Suspended Particle Transport

The sediment transport model is integrated until a steady situation is achieved. The source of sediment is river supply: the 7 main rivers discharging in the Alboran Sea have been considered. The suspended particle concentration in the inflowing Atlantic waters is fixed as 0.5 g/m3 [76] and a zero gradient condition is used along the other open boundaries.

The computed steady distributions of suspended particles in the upper and lower layers may be seen in Figure 9. Particle concentrations in the upper layer are about 0.4 g/m3 in the Strait of Gibraltar and in the Atlantic jet in the Alboran Sea. Concentration decreases in the WAG and in the south part of the sea to some 0.2 g/m3. Slightly higher concentrations, about 0.4-0.5 g/m3, are obtained in most of the bottom layer. These values are in agreement with [76], who found suspended particle concentrations below 0.5 g/m3 in several stations sampled in the Strait of Gibraltar and Alboran Sea (exact values of concentrations are not reported in such paper). The plumes of the two main rivers, Guadiaro and Guadalhorce (see Figure 5), may be appreciated in both layers. It seems that the distribution of particles in the upper layer is controlled by water circulation, with larger concentrations of particles along the Spanish coast, in the Atlantic jet, and a decrease in the region occupied by the WAG. Nevertheless, it is necessary to point out that these results correspond to a steady situation with constant water fluxes through the Strait and no winds. In practice, these particle distributions will be altered by changing water exchanges through the Strait, winds, changes in river particle supplies, and eventual episodes of disappearance of the WAG.

The vertical particle fluxes calculated by the model have been compared with those obtained from sediment traps [31] located at stations C and E (lower layer). Results are shown in Table 5. It may be seen that differences are, as maximum, a factor 1.5. Both measurements and the model give a higher flux in station E, which is more distant from the coast than C. This distribution is unusual for continental margins, where fluxes at a given depth generally decrease with distance from the coast. It has been attributed [31] to a funneling of particles from the periphery of the WAG towards its centre. However, it seems from Figure 9 that such process is not occurring in the upper water layer since suspended particle concentrations in the WAG are lower than in the Atlantic jet. It can be seen in Figure 9 that, in the lower layer, particle concentrations are slightly lower (and consequently also the fluxes) in the area of point C. Nevertheless, this effect is less apparent in the model results than in measurements and a clear explanation cannot be given.

A map of computed sedimentation rates (SR) is presented in Figure 10. Except in the river plumes, they are of the order of  g/cm2 year and are obviously correlated to the suspended matter concentrations in the deep layer (see Figure 9). The smallest SR is obtained in some areas of the Strait of Gibraltar and the Alboran Sea, where larger currents are present in the lower layer (Figure 8) and/or stronger tidal currents (Figure 7) exist. A comparison between measured and computed SR may be seen in Table 5. It seems, looking at Table 5, that essentially uniform SR are given by the model. However a clear structure is produced, with SR varying over one order of magnitude from the south to the north of the sea (Figure 10). The order of magnitude of the SR is correctly given by the model, although there is a general underestimation of SR over the Alboran Sea. This is not surprising given the approximations made in the model (Section  3.2).

4.1.3. Cs and Pu

In the case of Cs, the model has been run for the period 1985 (one year before Chernobyl accident) to 1997, when computed Cs concentrations in the bed sediment are compared with measurements in Masqué et al. [30]. Initial Cs concentrations for both water layers are obtained from the long-term box model in [77] and are assumed to be homogeneous over the entire domain. They are 3.0 and 2.0 Bq/m3 for the upper and lower layer, respectively. Concentrations in suspended matter and bottom sediments are assumed to be initially in equilibrium with the dissolved phase. Atmospheric and Chernobyl inputs are obtained from [77] as well. A simulation over 20 years has been carried out for Pu. The atmospheric deposition has been compiled from [76, 78], and zero activity concentrations in all phases are considered as initial conditions for the simulation. For both radioisotopes, concentrations in inflowing waters are specified [77] and a zero gradient boundary condition is applied where outflow occurs.

The computed distribution of Cs in the surface bed sediments at the end of the simulation may be seen in Figure 11 and a comparison with measurements [30, 79] is presented in Table 6. Computed concentrations are in the range 3–6 Bq/kg (the color scale in Figure 11 is limited to 5 Bq/kg for more clarity), in generally good agreement with measurements (Table 6). There are higher concentrations along both coasts of the Alboran Sea, particularly along the Spanish coast, where SR are higher due to the riverine input of particles. Thus, it could be deduced that Cs incorporated to the surface water from the atmosphere is attached to suspended particles that are later deposited on the bed. Indeed, lower concentrations are in general obtained in the Strait of Gibraltar and along the main path followed by the deep outflow current, where SR are smaller. Nevertheless, direct adsorption of Cs from the dissolved phase to the bed sediment is occurring simultaneously and an absolutely clear correlation between SR and Cs concentration in sediments cannot be established. Indeed, direct adsorption is higher in shallower areas, as can be seen from (4), and consequently the map in Figure 11 is affected by topography with a general trend of increasing Cs concentrations with decreasing water depth. Moreover, direct adsorption of radionuclides (and trace metals in general) is affected by the water-sediment contact time [80]. Thus, direct adsorption is smaller in regions of stronger currents in the bottom water layer, as it is the case with SR. It may be concluded that a clear correlation between contamination of the bed sediment and SR cannot be established, although a correct estimation of SR is of course required for an appropriate simulation of particle-reactive pollutant transport.

The mean inventory of Cs in bed sediments over the Alboran Sea has been calculated as well. This inventory includes radionuclides in the surface sediments plus radionuclides that have migrated down the sediment core according to (1). The obtained value, 158 Bq/m2 (Table 6) is in good agreement with the estimations carried out in [79] from measurements:  Bq/m2. The total inventory in bed sediments is 6.70 TBq. The total inventory in the water column (for both layers) has been calculated as well, resulting 75.2 TBq. However, there are not experimental estimations to compare with in this case.

Gascó et al. [81] have measured Cs concentrations in unfiltered water samples from the Strait of Gibraltar (surface and deep water) in 1997. Model results have also been compared with these measurements. Results are presented in Table 7, where it may be seen that the model gives realistic results for the dissolved phase as well. Essentially the same Cs concentrations are computed in surface and deep waters, in agreement with measurements [81].

The computed radionuclide concentrations corresponding to unfiltered water (in dissolution plus suspended particulate matter) are presented for both water layers in Figures 12 and 13, for Cs and Pu, respectively. The distributions in the upper layer are similar for both radionuclides, showing higher concentrations in the area occupied by the WAG. This is due to the fact that some water is trapped in the gyre (this fact will be discussed in detail in Section  4.1.6) and thus its radionuclide concentration increases as a result of the atmospheric input. This situation is in contrast with the distribution of suspended particles in the upper layer (Figure 9): lower concentrations in the WAG area may be seen in this case. The reason of the difference is that in the case of suspended matter the source is not atmospheric. Thus, waters in the WAG, which remain relatively isolated from their surroundings, have lower particle concentrations because of settling. Nevertheless, the same comment as when the distribution of suspended particles was discussed should be made: simulations are made under steady flow conditions, thus these patterns will be destroyed by episodes of migration and disappearance of the WAG and winds. However, it is worth commenting that the pattern of higher radionuclide concentrations (for Pb, whose main source is atmospheric deposition as well) in the centre of the oceanic gyres (in both the Pacific and North Atlantic) has already been detected [82, page 384].

A highly reactive radionuclide as Pu is quickly fixed to suspended matter and then, by settling, introduced into the deep water layer. Finally radionuclides are deposited on the seabed. Although it cannot be clearly seen in Figure 13 because of the scale in the colorbar, it is evident from Table 7 that computed concentrations in the deep layer are about a factor 2 higher than in the surface layer. This result is in agreement with the observations presented in Table 7, although the model is slightly underestimating concentrations. Thus, in the lower layer, distributions are clearly different for both radionuclides due to the different efficiency of the removal process of dissolved radionuclides by settling suspended particles. For sediments, only one measurement has been found to compare model results. Measured Pu concentration in the surface sediment at point J [76] is 0.550 Bq/kg, and the computed one is 0.379 Bq/kg.

The computed average inventory of Pu over the model domain is 45.8 Bq/m2. Gascó et al. [79] have obtained a value of  Bq/m2, although only one sediment core was analyzed. García-Orellana et al. [83] report values between 11 and 15 Bq/m2. Consequently, it can just be said that the order of magnitude seems correct. The total inventory in the sediment results 1.9 TBq. The computed inventory of Pu in the water column is 0.68 TBq, which agrees with the estimation of 0.64 TBq in [76].

4.1.4. Pb

The case of Pb is extremely complex. It is a member of the U radioactive decay chain and is supplied to the sea water from the atmosphere, rivers and in situ decay in the water column of its parent radionuclide Ra. In the bed sediment, Pb is produced by decay of Ra. Because of their desintegration periods and since both remain fixed to the sediment, these radionuclides are in secular equilibrium (their activities are the same) in the seabed. This amount of Pb, which is in equilibrium with Ra, is known as the supported fraction. Nevertheless, Pb is also incorporated to the bed sediment due to deposition of Pb containing particles from the water column. This is known as the excess Pb fraction. In practice, the difference between total Pb and Ra activities in the sediment gives the Pb in excess. Excess Pb has been measured in sediments of the Alboran Sea [30] and will be used to test the model results, since all the Pb is incorporated to the bed sediment from the water column in the model.

The atmospheric deposition of Pb in the western Mediterranean is 81 Bq/m2 year [84]. A more detailed evaluation may be seen in [85]; however the use of a mean value (which also is in good agreement with the average calculated by these last authors) is enough for our modelling purposes. In normal river water, Pb concentrations are low, about 0.2 Bq/m3 [86], and most of it is present in suspended particles. We have used this value for the Pb supply from river waters, assuming that activity in suspended particles discharged by the rivers is in equilibrium with water (thus deduced from (5)). It is not easy to calculate in situ production of Pb. The reason is that water containing Ra from the Odiel-Tinto estuary may reach the Strait of Gibraltar [81], thus a constant and uniform water Ra concentration cannot be assumed. The approach used has consisted of specifying Pb concentrations at the open boundaries, where inflow occurs, from estimations of Pb mean concentrations in the Atlantic inflow and Mediterranean outflow [81], and assuming that these concentrations contain the in situ fraction of Pb that is subsequently transported into the model domain. The limitations are obvious, since some of this in situ Pb will be lost by radioactive decay in the domain and its content in the water column will be underestimated. Moreover, in the case of the artificial Cs and Pu the model could be started from clean sediments, but this is obviously not the situation with the natural Pb. One possibility consists of using an initial concentration of excess Pb in both the mixed and deep sediment. However, this may look tricky since such initial concentrations could be fitted to obtain results that compare well with measurements. Thus, it was decided to start the model from a clean sediment, keeping in mind these limitations when analyzing the model results. Consequently results have to be interpreted with care and will be useful to have a general overview of the dynamics of Pb in the system, but a good comparison with experimental data cannot be expected.

The dispersion model is run for 35 years, time in which a steady distribution of Pb is achieved. Moreover, the effects of wrong initial conditions will have been reduced, given that the radioactive desintegration period of Pb is 22 years. The computed distribution of excess Pb inventories is much like the sediment Pu distribution (not shown). This is not surprising since their geochemical behaviour is similar (s of and  m3/kg for Pu and Pb resp.). It also seems that the contribution of rivers is negligible in comparison with atmospheric deposition. Finally, the Pb distribution in the water column presents the same enhancement in concentrations in the area of the WAG as we have seen before and is not shown.

A comparison between measured and computed inventories in sediments may be seen in Table 8. The inventories measured in [30] are included in the table. However, these authors have found that they are larger than those expected from a balance model for the water column (details may be seen in the reference). Thus the expected inventories, deduced from the balance model by such authors, are also given in Table 8. Masqué et al. [30] have concluded that the discrepancy may be due, at least partially, to the presence of a bottom nepheloid layer and to the occurrence of turbiditic flows in some areas. This highlights the difficulties of the study of Pb dispersion. Of course, these problems will be present when studying dispersion of all reactive tracers, although in some cases may not be as relevant as in the case of Pb and better agreement with observations would be obtained.

Nevertheless, it may be seen in Table 8 that, at least, the order of magnitude of the calculated inventories seems reasonable. Moreover, Masqué et al. [30] state that there is a general trend of increasing excess Pb inventories with water depth in the area of Málaga, and this trend is given by the model too.

The computed mean inventory of excess Pb in the sediments of the Alboran Sea is 2.98 kBq/m2 and the total inventory is 126 TBq.

4.1.5. Fluxes through the Strait of Gibraltar

The fluxes of the three studied radionuclides through the Strait of Gibraltar (through a north-south section at the longitude of Tarifa) have been evaluated from the model results. Indeed, several estimations of fluxes for some radionuclides and metals have already been done from measurements [9, 76, 81, 87, 88].

The water transport used in the estimations ranged from 0.69/0.65 Sv (inflow/outflow) (used by [89]) to 1.77/1.73 Sv (used by [90]). As a consequence of the variability of the water exchanges in the literature, radionuclide and metal fluxes may present differences up to a factor 3 for the same concentration of the substance [88]. A comparison of the model calculations with previous estimations may be seen in Table 9, where the water exchanges used in each calculation are also given. In spite of the variability of the water exchanges all results are in reasonable agreement and differences reach a factor 2 as maximum. In general, it seems that the model calculation tends to overestimate the inflow of Cs and both flows of Pb. It should be taken into account, however, that estimations are sometimes based on a single measurement of radionuclide concentration in the Strait.

4.1.6. Characteristic Times

The turn-over-time is defined as the time in which the tracer inventory in the water column decreases in a factor [91] and the sediment half-life is defined as the time in which the tracer inventory in the sediment decreases in a factor 2 [65]. These parameters are relevant for the water quality of a system. It is important to know the time scale for a pollutant discharged into a water body to be transported out of the system [92]. Thus, flushing times (that are equivalent to turn-over-times in the case of constant flow) have been recently determined for a number of water bodies using numerical experiments [92, 93]. On the other hand, it is also known that a contaminated sediment may act as a long-term delayed source of previously released contaminants [94]). Consequently, it is also relevant to have estimations of the sediment half-life [65].

Two numerical experiments have been carried out to determine flushing times and sediment half-lives in the cases of Cs and Pu. The model is started, for each radionuclide, from the computed radionuclide distributions in water, suspended matter and bed sediments presented above and simulations over 15 years are carried out. Atmospheric radionuclide inputs are cancelled and clean water, with zero radionuclide concentrations, is assumed to enter the domain through the open boundaries. The radionuclide inventories within the Alboran Sea in the region comprised between and longitude are evaluated each time step and results are written to a file. Thus, the system-wide (as defined by [93]) flushing time has actually been determined. These authors found that this parameter is useful for determining the long-term water quality of a system. Note that flushing times give an indication on the time scale required to flush-off pollutants, determined by the physical characteristics of the studied system and contaminant (water circulation, geochemistry of the pollutant, typical suspended matter concentrations and sedimentation rates). To evaluate the flushing time we need to assess at which rate a contaminant initially inside a given region decreases its concentration, without any additional source which would modify the result. Thus, all external radionuclide sources are removed.

As an example, the temporal evolution of the sediment inventories of Cs and Pu are shown in Figure 14. From numerical fitting to exponential decay curves, sediment half-lives are 136 and 778 days for Cs and Pu, respectively. Half-life for plutonium is about a factor 6 larger than that of caesium, revealing the higher affinity of plutonium to remain fixed to the solid phase.

Flushing times have been evaluated for each water layer separately and for the dissolved phase and suspended matter. The temporal evolution of the mean radionuclide concentration in each water layer is presented in Figure 15 for both radionuclides (dissolved phase). Choi and Lee [93] have found that the system-wide flushing time is better obtained from a double exponential decay curve rather than a single exponential decay. The temporal evolution of the mass within the system, , is thus given by the following equation: where is the initial mass within the system. If the three parameters , , and are determined from numerical fitting to the decay curves in Figure 15, the system-wide flushing time is given by [93]

In the case of dissolved Cs in the surface layer, for instance, a flushing time of 176 days is obtained if the corresponding curve in Figure 15 is fitted to a single exponential decay and the regression coefficient is . However, if the double exponential is used (7), the regression coefficient increases to and the flushing time is 292 days.

Flushing times have been determined for Cs and Pu, for both water layers. Results are summarized in Table 10. It can be seen that Pu flushing time is about one order of magnitude larger than for Cs, specially for the surface layer. This is again due to the higher affinity of Pu to be fixed to the solid phase, which makes it less mobile than Cs. Essentially the same values are obtained for the dissolved phase and suspended matter.

In the case of constant flow through a region, as is the case in our simulations, the flushing time can be obtained as [91], where is the region volume and is the flow. A flushing time of some 80 days is obtained if this equation is applied to the surface water layer with Sv. This is shorter than the flushing times presented in Table 10 because two reasons: firstly the fact that the WAG traps water (which stays longer in the AS) and, secondly, radionuclide interactions with the solid phase make them less mobile, which is of course more evident for Pu than for Cs.

4.2. Gulf of Cádiz
4.2.1. Hydrodynamics

A map showing the different locations and sampling points mentioned in the paper is presented in Figure 16. Computed tidal constants have been compared with established values for both and tides. Results are shown in Table 11, and it may be seen that, generally, there is a good agreement between both set of data. As an example, maps showing the tide amplitude and the amplitude of the tidal current are presented in Figure 17. The amplitude of the tide is about 1 m over all the GoC, decreasing near the Strait of Gibraltar. Associated currents are weak, with amplitudes below 0.10 m/s over most of the GoC. Indeed at the RAP (Red de Aguas Profundas, Spanish Institute of Oceanography) buoy position (see Figure 16), the barotropic tidal current is less than 0.03 m/s [33]. The computed current at this position is 0.034 m/s. The current amplitude increases as approaching the Strait entrance, where currents about 0.8 m/s are produced (the color scale in Figure 17 is limited to 0.5 m/s for more clarity). A similar behaviour is observed for the tide.

The residual surface circulation in the northern GoC is characterized by a current directed to the SE [13, 35] along the Spanish coast. This circulation is a rather constant pattern during most of the year. Part of the flow enters the Strait of Gibraltar and part is deflected to the south. The residual circulation computed with the baroclinic model at the sea surface is presented in Figure 18. The current is effectively directed to the SE over the continental shelf and part of this flow enters the Strait of Gibraltar. Maximum currents are of the order of 0.3 m/s, in agreement with [33]. The anticyclonic eddy at the east of Faro (see Figure 16) has been described by [34]. Also, the cyclonic eddy in front of the Strait of Gibraltar appears clearly in models [1, 96]. These last authors have attributed it to the strong convergence occurring in this area.

Below the surface, the Mediterranean waters flow into the Atlantic and mainly direct to the NW [38, 39, 97]. As an example, the computed circulation 590 m below the surface is presented in Figure 19. Only the northern part of the GoC is shown to appreciate details more clearly. These currents are in agreement with the geostrophic velocities below 400 m, referenced to 300 m, provided in [13], and with the calculations in [96] for summer. Water velocity is higher close to the Strait and then slows to about 0.1 m/s, in agreement with [97].

Temperature and salinity profiles in the water column calculated by the model have been compared with those obtained from observations carried out under the TOROS project [9] in summer 1997. Examples for 4 points shown in Figure 16 may be seen in Figure 20. There is an acceptable agreement between model results and experimentally obtained vertical profiles of and .

Finally, two vertical south-north sections of water salinity are presented in Figure 21, at longitudes of and . In the first case, the core of more dense Mediterranean water is clearly seen, with the salinity maximum below 600 m. This flow is aligned with the Spanish continental slope, in agreement with [96, 98]. Close to the Strait, the situation may be characterized as a two-layer exchange flow, with the interface tilted down southward [98].

Although seasonally averaged values have been used for open boundary conditions and wind stress, the main features of water circulation in the GoC are, generally speaking, reproduced by the model.

4.2.2. Suspended Particle Transport

The sediment transport model is run until steady particle distribution and SR are obtained. Open boundary conditions are the same as in the AS domain.

As expected, the river-discharged suspended matter is transported to the southeast along the Spanish coast, which is the residual current direction in this area of the GoC. As widely discussed [40], there is a dominant eastward transport throughout the entire northern GoC. Some authors [99] have postulated the existence of some westward transport of sediments released by the Guadalquivir River in the inner shelf, although little indication of it has been found in other works [40].

Maps of suspended matter concentrations at the surface and in the deepest water layer are presented in Figure 22. Logarithmic scale is used to appreciate differences. The suspended particle plumes produced by the three rivers can be clearly seen in the maps. Indeed, concentrations of the order of 10 g/m3 are obtained near the river mouths. Concentrations of the order of g/m3 are obtained in part of the northern GoC and much smaller values are apparent to the south. Computed surface particle concentrations are in agreement, by order of magnitude, with measured surface concentrations. Indeed, concentrations of the order of g/m3 have been measured at the surface in the northern GoC [100]. Also, concentrations in the range 20–45 g/m3 have been detected in the Guadiana River plume [101] over a narrow area (<10 km from the coast). In the deepest water layer concentrations are slightly higher than at the surface over most of the northern GoC, in agreement with [100]. In some areas as the Strait of Gibraltar and the Spanish continental slope, concentrations are more significantly enhanced in the deep water. This is probably due to resuspension produced by stronger tidal currents (in the area of the Strait) and also produced by the Mediterranean Water current.

Computed net sedimentation rates (not shown) are larger, as expected, in the vicinity of the river mouths, where they are of the order of  g/cm2 year. Much smaller values are obtained far from the river discharge influence. In particular, there is a region where no sedimentation occurs at the west of the Strait of Gibraltar. This is due to the strong currents in the Mediterranean outflow water, which keep particles in suspension. Nevertheless, all these results have to be interpreted with care due to the approximations which have been made in the model (the same as in the AS), described in Section  3.2. Also, although constant suspended matter concentrations have been defined in the three estuaries, there will be seasonal variations depending on pluviometry, for instance. It must be taken into account that suspended particles will affect tracer transport only close to the coast, given the extremely low suspended matter concentrations which are measured offshore [100], and river supply is the main sediment source in this coastal area.

4.2.3. Metal Transport

Metal concentrations in the dissolved phase have been defined at the three main sources: Guadalquivir, Guadiana, and Odiel-Tinto estuaries. Values are given in [1] and may be seen in Table 12. The metal transport model is started from background dissolved concentrations, which correspond to the open Atlantic Ocean metal concentrations reported in [9]. The corresponding background in the solid phase is obtained from the metal through (5). Equations are integrated until steady metal distributions are obtained in all phases. Open boundary conditions are again the same as in the AS domain.

The concentrations of several metals have been measured in the fine sediment fraction (<63 μm) along the Spanish coast from the Guadiana to the Guadalquivir mouths [3, 102] and also at some points closer to the Strait of Gibraltar [103]. Samples were collected at an approximate distance of 500 m from the shoreline. A comparison of measured and computed metal concentrations along the northern coast of the model domain can be seen in Figure 23. In general, there is a good agreement between the measured and calculated concentrations. Metal concentrations are very low westward from the Guadiana River. There is an increase in concentrations here since, as has already been mentioned, the three rivers drain the Iberian Pyrite Belt. Maximum concentrations exist in the mouth of the Odiel-Tinto rivers. Although river flows are much smaller than those of the Guadiana and Guadalquivir, the Odiel-Tinto rivers are considerably more contaminated [40] and, indeed, they have been recognized as the main source of metals along the coast [102].

A map of computed dissolved Zn concentrations at the surface for summer 1997 is presented in Figure 24 together with a contour plot obtained from empirical data [104], which shows essentially the same banded structure. It may be clearly seen that the impacts from river outflow are restricted to a narrow band along the shore. Samples were not collected in the coastal area from Cadiz to the Strait, thus the impact is apparently restricted (in the case of observations) in Figure 24 to the zone located to the north of Cadiz. The highest concentrations are obtained in the mouth of the Odiel-Tinto rivers, obviously as in the case of bed sediments. The plume of dissolved metals, however, reaches the Strait of Gibraltar. Indeed, it has already been found [9] that coastal waters transport dissolved metals from the Odiel-Tinto rivers to a distance of more than 200 km. Moreover, it has been found that these rivers constitute a source of natural radionuclides into the Mediterranean through the Strait of Gibraltar [81].

Once it seems that the model is giving a realistic representation of reactive transport in the GoC, a numerical experiment has been carried out to assess the effect of the outflow of Mediterranean Water on the Pu contamination of bed sediments in the GoC. The model is started from clean water, suspended matter, and bed sediments. Pu concentration in water is specified, as a boundary condition, in the Mediterranean Water entering the GoC through the Strait of Gibraltar. Such Pu concentration is obtained from the output of the AS model: 0.012 Bq/m3. The model is integrated until a steady situation is achieved. The computed Pu distribution in bed sediments is presented in Figure 25. These results could not be compared with measurements. However, the path followed by the Mediterranean Water is clearly marked on the bed sediment as a tongue of enhanced Pu concentrations directed to the noth-west. Its shape is in agreement with the computations in [38], Figure 4(b) in this paper, for salinity.

4.2.4. Sensitivity Analysis

The main parameter which describes the geochemical behaviour of tracers is the exchange velocity , obtained from the equilibrium distribution coefficient through (6), as already described. The model sensitivity to changes in this last parameter has been investigated since the is poorly defined and presents a high natural variability depending on environmental conditions (mainly pH, temperature and salinity, although other parameters as light intensity may be relevant as well).

Obviously, a higher means that the substance has a higher affinity to be fixed to the solid phase and vice versa. Thus, surface sediment metal concentrations are linearly correlated with the , as should be expected from the model formulation.

The distribution of metals in coastal sediments obtained with several values have been computed and are presented in Figure 26 in the case of Zn as an example. The value has been measured in coastal waters of the USA [105]. A smaller value, by one order of magnitude, has also been tested (). Effectively, too high metal concentrations, when compared with empirical data, are obtained with . On the other hand, concentrations are too low if a of is used. The best agreement between measured and computed metal concentrations are obtained with the recommended by IAEA [64], which is indicated in Figure 26 as the nominal simulation.

The sensitivity of the formulation of the water-sediment interaction processes to other parameters as particle size, and density, sediment mixing depth, fraction of small particles and the correction factor has already been studied in detail [56, 106] and will not be repeated here.

5. Conclusions and Future Work

Tracer dynamics in the Gulf of Cadiz and Alboran Sea have been studied through numerical modelling. Water circulation, required to solve dispersion, is obtained from appropriate hydrodynamic models for each domain. These models have been tested by means of comparisons of computed tides and currents with measurements in the GoC and the AS. Sediment transport models are also needed.

An interesting effect has been observed. It consists of an enhancement of radionuclide concentrations in the WAG. This effect has also been observed in the centre of oceanic gyres in the Atlantic and Pacific. It is due to the fact that water remains trapped in the gyre region, remaining relatively isolated from the surroundings, and thus accumulating radionuclides introduced from the atmosphere.

In the case of the AS, the model has been applied to three radionuclides: Cs, Pu, and Pb. They are mainly introduced from atmospheric deposition. Generally speaking, the model produces results in agreement with measurements. A correct estimation of Cs levels in bed sediments and in the water column is obtained. Also, the model reproduces the observed fact that essentially the same Cs concentrations exist in both water layers. Although the model slightly underestimates Pu concentrations in the water column, concentrations in the deep water layer are essentially a factor 2 higher than in the surface layer because of the efficient removal of dissolved plutonium by settling suspended particles. Excess Pb inventories in the sediment have been calculated. They are underestimated with respect to measurements, although this is due to the difficulties in modelling Pb: this is a natural radionuclide which is also produced by in situ radioactive decay of its parent radionuclide Ra. Also, initial conditions for the simulation cannot be appropriately defined.

Radionuclide fluxes through the Strait of Gibraltar have been evaluated. They are in relative good agreement with previous estimations based upon radionuclide concentration measurements in the Strait. Finally, system time-scales that are useful from a management perspective have been calculated. They are the sediment half-life and the system-wide flushing time. Both time scales are significantly larger for Pu than for Cs, which is due to the fact that reactive elements are less mobile in the environment than those remaining dissolved. Thus, flushing time for Pu in the surface water layer is one order of magnitude larger than that of Cs, for instance. Also, it has been found that the presence of the WAG increases the flushing time in the surface layer.

In the case of the GoC, the contamination of bed sediments by metals discharged by the main rivers draining the southern Iberian Peninsula has been adequately reproduced. The impact from river outflow is restricted to a narrow band along the shore, the area of the Odiel-Tinto river mouth being the more contaminated. The plume of dissolved metals reaches the Strait of Gibraltar, confirming that coastal waters transport dissolved metals from the Odiel-Tinto rivers to a distance of more than 200 km and that these rivers may constitute a source of metals and radionuclides into the Mediterranean through the Strait of Gibraltar.

The dispersion of Pu entering the Gulf of Cadiz with the outflow of Mediterranean Water has been simulated with the model. The path followed by the Mediterranean Water is clearly marked on the bed sediment as a tongue of enhanced Pu concentrations. Its shape is in agreement with earlier salinity computations.

Although both models give a generally good overview of pollutant dynamics in GoC and AS, a clear weak point is that they run using steady residual circulations. Models should be adapted to run with evolving circulation. This would be interesting to evaluate, for instance, the role of Mediterranean mesoscale eddies in pollutant and suspended sediment dynamics. Of course, modelling the behavior of the biotic components of the ecosystem would also be of great interest, given the high productivity of the region.

Appendices

A. Sediment Transport Model

In an three-dimensional form, the equation which provides the suspended matter concentration, , is where is the settling velocity of suspended particles and and are, respectively, the horizontal and vertical diffusivities. , and are water velocities along the , - and -axis. The deposition and resuspension terms are incorporated into the sea bed boundary condition of the equation. The deposition rate is written as where is particle concentration evaluated at the sea bottom, is bottom stress due to tides and the residual current, and is a critical deposition stress above which no deposition occurs since particles are maintained in suspension by water turbulence.

The settling velocity is determined from Stokes’s law: where and are suspended particle density and diameter, respectively, is kinematic viscosity of water and its density. The resuspension rate is written in terms of the erodability constant: where is the erodability constant, gives the fraction of fine particles in the bed sediment, and is a critical resuspension stress below which this process does not occurs. The model can also calculate sedimentation rates (SR) as the balance between the deposition and resuspension terms.

B. Tracer Dispersion Model

The equation that gives the time evolution of tracer concentration in the dissolved phase, , is where and are, respectively, the concentrations of tracers in suspended matter and the fine fraction of bottom sediments (particles with diameter <62.5 μm, as mentioned in Section  3.3). All remaining parameters have been defined before (Section  3.3). The last three terms of the equation correspond to the tracer transfer from water to suspended matter and bed sediment, from suspended matter to the dissolved phase and from the bed sediment to the deepest water layer dissolved phase.

The equation that gives the time evolution of tracer concentration in suspended matter is where is the particle settling velocity. The last three terms are the tracer transfer from the dissolved phase to suspended matter, from suspended matter to the dissolved phase and the deposition of metals (fixed to suspended particles) from the deepest water layer to the sediment, evaluated according to

Note that means that the corresponding magnitude is evaluated at the deepest water layer.

The equation for the temporal evolution of tracer concentration in the bottom sediment mixed layer is where deposition is now calculated as

The total tracer content, , in the sediment below the mixed depth is given by the following equation:

Radioactive decay, , should be added to the equations to simulate a radioactive isotope.

Acknowledgments

The author is indebted to the Agencia Española de Cooperación Internacional for partially funding this work and to the Spanish Ministerio de Educación y Ciencia for a fellowship to stay during three months at the University of Wales, Bangor, where part of this work was carried out. This work is partially supported through the project “Determination of scavenging rates and sedimentation velocities using reactive-particle radionuclides in coastal waters; application to pollutants dispersion modeling,” Ref. CTM2009-14321-C02-01-02, funded by the Spanish Ministerio de Ciencia e Innovación.