Modeling of 226Ra behavior in a Spanish estuary affected by the phosphate industry

Modeling of 226Ra behavior in a Spanish estuary affected by the phosphate industry R. Periáñez,1* A. Absi,2 M. Villa,3 H. P. Moreno,4 G. Manjón2 1 Departmento Física Aplicada I, EUITA, Universidad de Sevilla, Ctra Utrera km 1, 41013–Sevilla, Spain 2 Departmento Física Aplicada II, ETSA, Universidad de Sevilla, 3 Departmento Física Atómica y Nuclear, Universidad de Sevilla, 4 Departmento Física Aplicada II, EUAT, Universidad de Sevilla, The Odiel and Tinto rivers, southwest Spain, form a fully mixed estuary. An industrial area that includes a complex dedicated to the production of phosphate fertilizers is located by the Odiel River. This complex released phosphogypsum wastes directly to the Odiel River and also disposed them on open air piles located by the Tinto River. Due to new EU regulations, wastes are not directly released to the Odiel from 1998 on, although they are still disposed on the open air piles. The behavior of 226Ra in a system like this estuary is complex, since radionuclides are affected by tidal actions and interactions with sediments through adsorption/desorption reactions and erosion/deposition processes. A numerical 2D depth-averaged model of the estuary has been developed, including processes mentioned above. It has been applied to reproduce experimental data measured after a release from the industrial complex in the Odiel River and after an accidental release in the Tinto River from the gypsum piles. The model has also been applied to simulate the self-cleaning process observed in the estuary after the direct releases from the fertilizer complex were stopped. Introduction The Odiel and Tinto rivers, in the southwest of Spain, form a fully mixed tidal estuary which surrounds the town of Huelva (Fig. 1). Both rivers join at the Punta del Sebo. From this point, they flow together to the Atlantic Ocean. An industrial complex, containing a plant dedicated to the production of phosphoric acid and phosphate fertilizers, is located by the Odiel river. The phosphate rock used as a raw material by this industry contains significant amounts of natural radionuclides, mostly U, Th and Ra. The industrial processing of the phosphate rock leads to a redistribution of radioactivity. For instance, during the wet process for phosphoric acid production, 80% of the Ra follows the so-called phosphogypsum. Phosphogypsum was partially released directly into the Odiel river (20%), with the remaining 80% conducted with water through a pipeline to phosphogypsum piles located by the Tinto river (Fig. 1), where such material is stored in the open air. The gypsum piles cover some 12 km2 of the Tinto river margin. Since 1998, wastes are not released directly into the Odiel river due to new regulations from the EU, although phosphogypsum is still being disposed in the piles by the Tinto river. In our first paper,1 the time evolution of 226Ra activities in water and sediments over the years 1999– 2002 was studied. Results indicated that a self-cleaning process was taking place, as a consequence of the new waste policy, since a systematic and continuous decrease in activities was found in the water column and in the * E-mail: [email protected] bed sediments. Nevertheless, interpretation of experimental results is not an easy task due to the fact that several processes are affecting the radionuclide dynamics: advective transport of dissolved radionuclides along the estuary because of tides, turbulent diffusion and transfers of radionuclides between the dissolved and solid phases (suspended matter particles and bottom sediments). The solid phase may act as a sink or as a source of radionuclides to the water column, depending on the concentration of radionuclides in each phase and on the rates governing the adsorption and release reactions. Thus, a numerical model of the estuary, including all these processes, has been developed. It solves the hydrodynamics of the estuary and the dispersion of 226Ra, including the interactions between the dissolved and solid (bottom sediments) phases. The model has been initially tested through its application to simulate two dispersion scenarios. They consist of reproducing the measured 226Ra concentrations resulting as a consequence of releases carried out by the fertilizer complex in the Odiel river and a discharge into the Tinto river due to an accidental release from the phosphogypsum piles. Once tested, the model is applied to simulate the cleaning process of the estuary. Description of the model The 2D depth-averaged hydrodynamic equations are solved to obtain water currents and surface elevations required to solve the advection/diffusion dispersion equation for radionuclides. Equations are standard2 and solved by explicit finite differences. The two main tidal constituents M2 and S2 have been considered. The dispersion model is run off-line. Thus, tidal analysis is applied to calculate tidal constants that will be used by the dispersion code to obtain the current and water elevation at any time and position in the estuary. Adsorption and desorption reactions are described in terms of kinetic transfer coefficients. Thus, the adsorption process will be governed by a coefficient k1 and the inverse process by a coefficient k2. However, very recent experiments3 have shown that a two-step kinetic model consisting of two consecutive reactions is more appropriate than a one-step model, consisting of a single reversible reaction, to simulate both the sorption and release kinetics. This 2-step model has already been tested in other environments, and has been shown to describe the process of redissolution of radionuclides from contaminated marine sediments.4 The exchange model used3 to describe the experiments considers two successive reversible reactions. The first describes a reversible isotopic or ion-exchange process between dissolved radionuclides and some non-specific sites, S1, present on the particle surfaces. The second slower reaction represents a reversible sorption to more specific sites, S2. They can be represented as follows: XRSXSR k k m o 11 1 2(1) 1212 3 4SRSRSS k k m o (2) where k3 and k4 are the kinetic transfer coefficients, or sorption and release velocities, respectively, for the second reaction, R is the dissolved radionuclide, X is a competitive element that can be replaced by R on sites S1 and RSi is the radionuclide bound to sites Si of the solid particle. The adsorption process is a surface phenomenon that depends on the surface of particles per water volume unit into the grid cell. This quantity has been denoted as the exchange surface. As a first approach, assuming spherical sediment particles and a step function for the grain size distribution of particles, a mathematical formulation for the exchange surface can be obtained. Details are not given here, but may be seen elsewhere.2,4 Fig. 1. Map of the area of the estuary covered by the model. Numbered circles indicate the points where samples were collected. Lettered triangles indicate the points where currents measurements were available. Units on the axes give the grid cell number. The sea is approximately 1 km to the south of point 1 The equation that gives the temporal evolution of specific activity in the dissolved phase, Cd (Bq/m3), is: IU fLAkHCk y C HK yx C HK x y vHC x uHC t HC ssd d D d D ddd 21 )()()(   ¸ ¸ ¹ · ¨ ¨ © § w w w w  ¸ ¸ ¹ · ¨ ¨ © § w w w w w w  w w  w w (3) where H is the water depth, u and v are the two components of the water current, KD is the diffusion coefficient, L is a mixing depth in the sediment (the distance to which the dissolved phase penetrates the sediment), f gives the fraction of active sediment (particles with diameter <62.5 Pm, since adsorption by larger particles can be neglected), I is a correction factor that takes into account that not all the exchange surface of the particles is in contact with water because part of it may be hidden by other sediment particles, As (Bq/kg) is the activity concentration in the non-specific sites of the active sediment and U s is the sediment bulk density expressed in kg/m3. The external source of radionuclides should be added to this equation at the points where it exists. The equation for the temporal evolution of activity concentration in the non-specific sites of the active sediment fraction is: sss s ds AkAkAk fL HC k t A* 4321  w w I U (4) where *As is the activity concentration in the specific sites of the active sediment. The equation for the specific sites is: ss sAkAk t A* *43  w w(5) The total concentration of radionuclides in the sediment, Atot, is computed from: )*( sstot AAfA  (6) These equations are solved using explicit finite difference schemes with appropriate boundary conditions. In particular, second order accuracy schemes are applied for the advection and diffusion terms. Results and discussion The calibration of the hydrodynamic model consisted of selecting the optimum value for the bed friction coefficient. A reasonable agreement between observed and computed currents has been achieved5 for several points along the estuary (Fig. 1). The parameters required by the dispersion model, as kinetic rates or sediment density, were obtained from the experiments, direct measurements in the estuary and from the literature.6 The model was first applied to reproduce the 226Ra concentrations measured in water in three sampling campaigns carried out in 1990, 1991 (direct releases from the fertilizer complex were carried out in the estuary) and 1999 (no direct releases). Two source terms are considered: input from the factories in the Odiel river and from the phosphogypsum piles in the Tinto river. Unfortunately, the magnitude of the sources is not known. Thus, they had to be selected from trial and error exercises. For the year 1990, results are presented in Fig. 2a. The model gives a realistic distribution of 226Ra over the estuary, producing activity levels in good agreement, in general, with measurements. The dashed line corresponds to a simulation in which the source due to run-off from the phosphogypsum pile is not included. If it is not, activity levels in the region of sampling point 15 are slightly lower than the measured level. Model results are also in good agreement with observations for 1991 samples (Fig. 2b). Again, the dashed line corresponds to the simulation without the run-off source. Due to the lower input from the plant, the effect of runoff is clearly more relevant than in 1990 (effectively, it seems that run-off from the phosphogypsum pile is masked by the high direct releases from the complex in the 1990 campaign). The only source considered in the 1999 simulation was run-off, since direct discharges ceased in 1998. It can be seen that the model gives a realistic representation of activity levels along the estuary (Fig. 2c). The simulation has been repeated exactly in the same way but without considering interactions with the sediments (dashed line in Fig. 2c). It seems clear that the sediments are now acting as a source of 226Ra to the water column since computed activity levels are significant lower than the measured ones if the sediment is not included in the calculations. The difference between the simulations is smaller along the Tinto river, where the input from the phosphogypsum piles is the dominant source. It seems, on the other hand, that the most important source of 226Ra to the Odiel river is the redissolution from the contaminated sediment, as indicated by the difference in activity levels between the simulations. The computed 226Ra specific activities in bottom sediments resulting as the consequence of releases from the fertilizer complex are presented in Fig. 3 together with the measurements.7 Maximum concentrations are obtained in the Odiel river, in the area where discharges from the fertilizer complex occur, with decreasing activities going north and south from this area. An activity peak is also measured in the Tinto river, that is not reproduced by the model. It must be taken into account that sampling was carried out during a wet season and rains may produce run-off events from the phosphogypsum piles. These run-off episodes, not included in the simulation, are probably responsible for the activity peak detected in the Tinto river. It is also interesting to note that activity gradients in the transverse direction are apparent in the Odiel river. On the other hand, the Tinto river is relatively mixed in the transverse direction. The computed radium concentrations in the sediment are used as initial conditions to study the self cleaning processes in the estuary after the end of discharges from the fertilizer complex. The distribution of radium between the two sediment phases is also provided by the simulation: some 40% of radium in the sediment is in the nonspecific sites. These conditions are assumed to be representative of the moment at which releases ceased. A 4-year simulation is then carried out starting from these initial conditions and without external releases. This simulation provides the 226Ra distribution in sediments for the years 1999, 2001 and 2002 to be compared with measured distributions.1 The computed time evolution of 226Ra inventories in bed sediments for the whole estuary are presented in Fig. 4. It can be seen that the total radium content in the estuary bed decreases due to redissolution. The sediment half-life (time in which activity in the sediment decreases in by factor 2) is 510 days. This number agrees relatively well with the 630 days half-life estimated from the measurements.1 The computed distributions of 226Ra in sediments corresponding to years 1999, 2001 and 2002 are presented in Figs 5, 6 and 7, respectively, together with the measured concentrations.1 These three figures provide a general view of the sediment cleaning process over the estuary. There are some local activity maxima in the Tinto river not reproduced by the model which are probably due to episodic release events from the phosphogypsum piles, but, in general, the decrease in 226Ra concentrations is well reproduced by the model, especially in the area of the Odiel river where direct discharges took place. If these maps are compared with that presented in Fig. 3, it also seems that the cleaning process is more efficient in this region, which is the most contaminated. In the rest of the estuary, the specific activity in the sediment decreases more slowly. Fig. 2. Comparisons between measurements and model results for three sampling campaigns: (a) 1990, the dashed line corresponds to a simulation in which the run-off source is not considered; (b) 1991, the dashed line has the same meaning as above; (c) 1999, the dashed line corresponds to a simulation in which the sediments are not included in the model. Sampling points are indicated in Fig. 1 Fig. 3. Computed distribution of 226Ra in sediments (Bq/g) and measured concentrations Fig. 4. Time evolution of the 226Ra inventories for the whole estuary in bed sediments during the self-cleaning process Fig. 5. Computed distribution of 226Ra in sediments (Bq/g) and measured concentrations for the year 1999 Fig. 6. The same as Fig. 5 but for the year 2001 Fig. 7. The same as Fig. 5 but for the year 2002 Conclusions A model to simulate the dispersion of 226Ra in a tidal estuary is described. The model solves the hydrodynamic equations off-line and standard tidal analysis is carried out to determine the tidal constants required by the dispersion code to calculate the instantaneous currents over the estuary. The dispersion model includes advective/diffusive transport and exchanges of radionuclides between the dissolved phase and bed sediments. These exchanges are described by means of kinetic transfer coefficients. The model has been initially tested through the comparison of computed and measured activity levels resulting after external releases occurring into the Odiel and Tinto rivers from the fertilizer complex and phosphogypsum piles, respectively. Once tested, the model has been used to simulate the self-cleaning process of the estuary observed after the fertilizer complex ceased its direct releases into the Odiel river. In general, temporal trends of activity concentrations in bed sediments are reproduced by the model. The computed half-life of the sediment is 510 days, that is in relative good agreement with the 630 days value deduced from observations. As a general conclusion, although a reduction in activity levels through the years has been obtained, the process of cleaning of a given estuarine environment will probably be slower than might be expected after a reduction (or cessation) of waste disposal due to the redissolution of pollutants from the previously contaminated sediments. References 1. A. ABSI, M. VILLA, H. P. MORENO, G. MANJÓN, R. PERIÁÑEZ, Sci. Total Environ., 329 (2004) 183. 2. R. PERIÁNEZ, Modelling the Dispersion of Radionuclides in the Marine Environment, Springer-Verlag, Heidelberg, 2005. 3. P. CIFFROY, J. M. GARNIER, M. K. PHAM, J. Environ. Radioact., 55 (2001) 71. 4. R. PERIÁÑEZ, J. Environ. Radioact., 71 (2004) 243. 5. R. PERIÁÑEZ, Estuar. Coastal Shelf Sci., 54 (2002) 809. 6. R. PERIÁÑEZ, A. ABSI, M. VILLA, H. P. MORENO, G. MANJÓN, Sci. Total Environ., 339 (2005) 207. 7. A. MARTÍNEZ-AGUIRRE, M. GARCÍA-LEÓN, C. GASCÓ, A. TRAVESI, J. Radioanal. Nucl. Chem., 207 (1996) 357.