The Mathematical Modelling Of Radionuclide Transport By Surface Water Flow From The Vicinity Of The Chornobyl Nuclear Power Plant

The hierarchy set of the models { box model, one-dimensional model, two-dimensional model and three-dimensional model are used to simulate local phenomena in the Prypiat river at the Chornobyl Nuclear Power Plant and radionuclide transport into the whole Dnieper river reservoirs system. Each model at its speciic level of resolution simulates the ow dynamics, suspended sediment transport, radionu-clide transport in dilute and on suspended sediments, radionuclides in bottom depositions. The modelling results were used to support the post-accident measures.


Introduction
The Chornobyl accident heavily contaminated the Dnieper river watershed (the territory of Byelorussia, the Russian Federation and Ukraine) including the watershed of the Prypiat river { a right tributary of the Dnieper.The Chornobyl Nuclear Power Plant (NPP) is located on the bank of the Prypiat river, approximately 30 km from the Prypiat river out ow to the Kiev reservoir of the Dnieper river.The oodplain territory near the Chornobyl NPP and the surrounding watersheds are heavily contaminated by 137 Cs and 90 Sr.The spots of 137 Cs are in the upper Dnieper watershed on the Russian and Byelorussian territory and in the whole Prypiat watershed.This surface contamination leads to the permanent in ux of 137 Cs and 90 Sr into the Kiev reservoir and consequent radionuclide transport through the whole cascade of six Dnieper reservoirs 900 km downstream to the Black Sea.The magnitude of this distributed source increases during each spring ood generated by snow melting and during high rainfall oods in the Prypiat river.As a result, the Ukrainian population consuming the Dnieper water downstream the from the Prypiat mouth to the Black Sea is under the impact of Chornobyl radionuclides.Also shery and irrigation of agricultural products are linked with the water usage from the Dnieper reservoirs.Therefore, even 10 years after the accident the population is very sensitive to this problem.
Since May 1986 in IPMMS, Kiev has been developing a computerized system for processing the Dnieper basin radiological monitoring data and modelling radionuclide dispersion in rivers and reservoirs.The objectives of the computerized system are formulated as follows: 38 M.J.Zheleznyak seasonal and long-term prediction of the surface water radioactive contamination; decision support for the aquatic post-accident measures directed to diminishing radionuclide uxes from the Chornobyl area through the Prypiat river and the Dnieper reservoirs; decision support for the measures directed to changes in the water consumption.Due to these objectives the computerized system should simulate radionuclide transport within a wide range of temporal and spatial scales { from ten meters to thousand kilometers and from seconds to decades.The multiple scale analysis of the problems based on the hierarchy set of hydraulics and pollutant transport models is provided 1-4].Radionuclide transport models have been formulated, of which equations were derived from the equations of the basic three-dimensional model by their averaging over di erent scales.The rst operational version of the computerized system has been running since September 1986.It has been used for forecasting radionuclide dynamics in the Dnieper reservoirs during the period of autumn heavy rainstorms, and also for the operative evaluation of the e ciency of aquatic measures.In 1987 the system was for the rst time used for the prediction of the Dnieper reservoirs contamination during the spring ood.Within a decade after the accident the system was further developed and used for the prediction of radionuclide concentration in water bodies and also for the water protection decision support.A brief overview of the methods and the results of the analysis based on this model is presented in the paper.

Radionuclide transfer processes in surface water
The development of the model set of di erent scales was done on the basis of the consideration of radionuclide transport processes.A most comprehensive analysis of the problem was made before the accident by Y. Onishi 5].A unique set of the post-accident data lead to better understanding the radionuclide fate in aquatic systems 6].
In water bodies radionuclide is transferred as a dissolved contaminant, a particular contaminant (pollution adsorbed by sediment) and also is transported by biota.The latter phenomenon is important mainly for lakes.
The main ways of radionuclide transfer in water and sediments are as follows: Dispersion of dissolved radionuclides by the water ow.The process is driven by the ow hydraulics.It includes advection and turbulent di usion.
Dispersion of radionuclides inherent to suspended sediments.The process is driven by the suspended sediment transport in the river/reservoir ow.A rough quantitative estimation of the ratio of the ux of radionuclides carried by the suspended sediment to the horizontal ux of radionuclides in dilute can be done by the formula where n is the number of typical grain size fractions of the suspended sediments by which the suspended sediment fractional distribution (i.e., clay, silt, sand) could be represented, S i Kg/m 3 ] { the concentration of i-fraction of suspended sediment in river water; K di m 3 /Kg] { the equilibrium distribution coe cients for i-fraction of suspended sediments.
As an example, for plain rivers a typical range of the total suspended sediment concentration is 10 1 ?10 2 Kg/m 3 .It is clear from formula (1) that radionuclide suspended sediment transport is very important for 137 Cs (typical K d value range is 1 ?10 m 3 /Kg), and it has minor signi cance for 90 Sr (typical K d order of magnitude is 0:1 m 3 /Kg).The isotopes of plutonium whose K d order of magnitude is 10 3 5] are carried out practically only on suspended sediments.Deposition/resuspension.This process can accumulate or deplete radionuclides in the bed.It is controlled by hydraulic factors (e.g., river ow, sediment transport), and signi cantly depends on the typical sediment size fractions (e.g., clay, silt, sand) in the ow.
Di usion at water-bottom interface.Radionuclide di usion between insitual water in the bed and overlying water can create a long term ux of pollutant directed to less contaminated compartment of this system.
Physico-chemical exchange processes in the system \water-suspended sediments".A pollutant transfer between river water and suspended sediment is driven by the adsorption-desorption processes.
Physico-chemical exchange processes in the system \water-bottom sediments".The main chemical mechanism is radionuclide adsorption to and desorption from the surface bed sediment going simultaneously with the above presented di usion to insitual water.For radionuclide transfer between the bottom and water column most important exchanges occur within a thin top layer of the river bed.The process is controlled by geochemical reactions of radionuclides with the river water and sediment, and is not always completely reversible.Chornobyl studies show the importance of the transformation of chemical species of radionuclides in sediments, i.e. transfer of non-exchangeable forms into exchangeable ones 6].
Uptake and subsequent decay and excretion of radionuclides by aquatic biota.The process transfers radionuclides from water to bed sediment.Sometimes water-bottom interface processes (di usion, adsorption-desorption, biota uxes) can be de ned as direct uptake processes.
Transfer processes between the upper bottom layer and deep buried sediments.Through the thin top bed layer radionuclides can be accumulated in or depleted from the deeper bed sediment.These radionuclides are further mixed within the deeper river bed layers by di usion, bioturbation and movements of bed sediment formation.

General Approach
The modelling of radionuclide dispersion in water bodies is part of a more general problem of water quality modelling, including the modelling of hydrodynamic (hydraulic) processes, and the simulation of sediment and pollutant transport driven by such hydrodynamics (see, e.g., 7]).
The radionuclide transport models, independently of a spatial resolution, include two main types of submodels { hydraulic ones, describing water, suspended sediment and bottom dynamics, and submodels concerning the fate of radionuclides in di erent phases driven by these hydraulic processes.
Modelling the fate of radionuclides in three di erent phases { in solution, in suspended and in sediment deposition is particularly important.
To describe the adsorption/desorption and di usion contamination transfer in the systems \solution { suspended sediments" and \solution { bottom deposition" the K d approach has been used.The exchange rates between the solution and the particles are taken into account for simulation kinetics processes.The adsorption and desorption rates are assumed not to be equal.The contamination exchange between the bottom deposition and suspended sediment is described taking into account sedimentation and resuspension processes.The developed approach is similar to the one developed in 5,8,9].The two step kinetic model describing the transfer of non-exchangeable forms into exchangeable ones was also used in some cases.The transport equations of the 3-D model THREETOX and the 2-D model COASTOX are presented below.The models of a lower dimension were derived from them by consequent averaging.

THREETOX: 3-D model
It is reasonable to provide a three-dimensional modelling of radionuclide transport in water bodies under the conditions of a complex con guration of a coast line and bottom topography, and large vertical and horizontal gradients of hydrophysical and pollution elds.Such conditions could be present in the vicinity of a heavily contaminated bottom or in the area of release, in strati ed water bodies, such as deep reservoirs, estuaries, cooling ponds.
The THREETOX code was developed to simulate 3-D hydrodynamic elds, suspended sediment and radionuclide transport in water bodies 9].The hydrodynamics is simulated on the basis of the three-dimensional, timedependent free surface, primitive equation, Princeton Ocean Model 10].
The model equations are written in the Cartesian coordinates (x; y; z).The prognostic variables of the hydrodynamic code are three components of the velocity eld (U; V; W), temperature, salinity and surface elevation h.A water body is assumed to be hydrostatic and incompressible.The concept of eddy viscosity/di usiveness and the Prandtl hypothesis with the variable turbulence length scale are used to de ne the turbulence stresses.At free surface all uxes (momentum, heat, etc.) are prescribed.At the bottom and the land boundaries the conditions of non-di usive uxes of any property is used.The open lateral boundary conditions are modi ed by radiation ones.
The advection-di usion equation governing the transport of suspended sediment is @S @t + U @S @x + V @S @y + (W ?W o ) @S @z = @ @z ( @S @z ) + A( S); ( where S is the sediment concentration, W o { settling velocity of the particles, and A are vertical and horizontal di usiveness, respectively.At the free surface no vertical ux of the sediments is assumed, i.e. (W ?W o )S = @S @z ; z = : (3) The vertical ux of suspended sediments at the bottom is assumed to be equal to the di erence of the resuspension and sedimentation rates.It yields the following formulation of a boundary condition: @S @z + W o S = q s ?q b ; z = ?h? a (4) The mathematical modelling of radionuclide transport 41 where a stands for the level of the boundary between suspended sediment and bottom sediment transport, q s ; q b are sedimentation and resuspension rates, respectively.It is supposed that in the case of ne noncohesive sediments these rates may be estimated as q s = ( W o (S ?S o ); S o > S ; 0; q b = ( 0; S o > S ; E r W o (S ?S o ); S o < S : (6) Here E r is the erosion coe cient, S o is the actual sediment concentration at the bottom level z = ?h?a, S { the near-bottom equilibrium sediment concentration that corresponds to the sediment capacity of the steady and uniform ow with the same local parameters.The semi-empirical formulas of Bijker or Van Rijn are used for the calculation of the equilibrium concentration.
The submodel of radionuclide transport describes the radionuclide concentration in solute C, that is in the suspended sediments as well as in the bottom deposition.The transport equation of dissolved radionuclides is @C @t + U @C @x + V @C @y + W @C @z = @ @z ( @C @z ) + A( C) ?C ? a 1:2 S(K d C ? C s ); (7) where K d is the distribution coe cient in a solid particle-water system in equilibrium conditions, K d = C s =C when t ! 1 under hydraulic steady conditions; a 1:2 is the rate of a water-suspended sediment exchange; { the decay constant for a radionuclide under consideration.
The boundary condition at the free surface is @C @z = WC; z = : The di usion ux into the bottom is taken as @C @z = s (1 ?")Z a 1:3 (K db C ? C b ); z = ?h?a; (9) where " is the bottom porosity, Z is the e ective thickness of the contaminated upper bottom layer, a 1:3 is the rate of a water-bottom exchange, s { the sediment density, K db { the distribution coe cient for the bottom sediments.
Radionuclide transport by suspended sediments is described by the equation @SC s @t + U @SC s @x + V @SC s @y + (W ?W o ) @SC s @z = @ @z @SC s @z !+ A( SC s ) ?SC s + a 1:2 S(K d C ? C s ); (10) Boundary conditions imply that there is no ux of C s through the water surface, and also that its ux is equal to the di erence of amount of eroded and deposited pollutants at the bottom boundary: (W ?W o )SC s ?@SC s @z = 0; z = ; W o SC s + @SC s @z = C s q s ?C b q b ; z = ?h?a: (12) The thickness of the upper layer of the contaminated bottom deposition is governed by the equation of the bottom deformation s (1 ?") @Z @t = q s ?q b : (13) Radionuclide concentration in the upper bottom layer is described by the equation @(Z C b ) @t = a 1:3 Z (K db C ? C b ) ? 1 s (1 ?") (C s q b ?C b q s ): The governing equations together with the boundary conditions were solved by the nite di erence method (FDM).The splitting of the barotropic and baroclinic modes imposed on the code going to sigma co-ordinates was provided to avoid di culties in the numerical solution for realistic bottom topography.A horizontally and vertically staggered mesh used for the computations.

COASTOX: 2-D lateral-longitudinal model
Integration of equations ( 2), ( 7), (10) over the depth with the presented above boundary conditions derives the 2-D equations which, together with equations ( 13),( 14), create the model COASTOX 1,2].The hydrodynamic part of the model is 2-D Saint-Venant equations (shallow water equations).The 2-D model is used to analyze radionuclide dispersion in shallow water bodies with signi cant spatial variations of the concentrations (vicinity of the releases, transport above the inhomogeneously contaminated bottom, river oodplains, and shallow reservoirs).Two-dimensional hydrodynamics equations are as follows: @U j @t + U k @U j @x k + g @ @x k = ?b U j j ?!U j + w W j j ?!W j; (15) @ @t + @(hU k ) @x k = R; (16) where b is the bottom friction coe cient; w is the wind friction coe cient on free water surface, k = 1; 2; j = 1; 2; U j is the vertically averaged horizontal ow velocities; W is the horizontal wind velocity above the water surface; and R are the sinks and sources distributed on the free ow surface (precipitation and evaporation).The depth-averaged equation of the suspended sediment transport is derived from equations ( 2)-( 6) as follows: @hS @t + @ @x k (hSU k ) = @ @x k hE ik @S @x i !+ w o (S ?S); (17 where E ik is a component of the dispersion coe cient; is the parameter that includes E R from equation (6).The quantity describes a ratio between the near bottom and depth-averaged concentration.
Taking into account boundary conditions at the bottom and free water surface, the depth-averaged equation of dissolved radionuclide transport can be derived from equation ( 7): @hC @t + @ @x k (hCU k ) = @ @x k hE ik @C @x i !?hC ? ha 1:2 S(K ds C ? C s ) ? s (1 ?")Z a 1:3 (K db C ? C b ); where C and C s are the depth-averaged radionuclide concentrations in solution and suspended sediments, respectively; C b is the radionuclide concentration in the bottom deposits averaged over the contaminated layer thickness Z .
Integrating equations ( 10)-( 12) over the water ow depth we can write the transport equation for C s as @hSC s @t + @ @x k (hSC s U k ) = @ @x k hE ik @SC s @x i !?hSC s + ha 1:2 S(K ds C ? C s ) + C b q b ?C s q s : ( Equations ( 15){( 19) and ( 13){( 14) are the basis of the computer code COASTOX.The nite-di erence method is used to solve the model's equations { the Godunov scheme for a shallow water equation, the second order scheme or the 4-th order Holly-Preissmann scheme for transport equations.

VERTOX: 2-D vertical-longitudinal model
The model is developed to simulate radionuclide transport in river canals in the vicinity of the abrupt depth changes { such as bottom traps for contaminated sediments and other engineering constructions 1,2,4].The model equations have been derived by averaging equations of the basic 3-D model ( 2)-( 14) over the river channel width.The numerical solution of the model equations is obtained by FDM.

RIVTOX: 1-D river channel model
The model is used to simulate radionuclide transport in networks of river canals generated by the direct pollutant release into the rivers or by the washing out of radionuclides from the catchments 3].The hydraulic part of the model are Saint-Venant equations or their \di usive wave" approximation.The equations of the RIVTOX model are derived by averaging equations of the 3-D model ( 2)-( 14) over the river canal cross-sections.
The FDM developed for the computer code CHARIMA 7,13] was used for the numerical solution of the Saint-Venant equation for river network.The fourth order of accuracy schemes were used to solve the advectiondi usion equations.

WATOX: box model
The model was derived to simulate radionuclide transport in the whole set of the Dnieper reservoirs on seasonal and long-term temporal scales.Equations ( 12){( 14) are averaged over compartments by representing the whole reservoir or its large section 1,2].The resulting set of ordinary di erential equations was solved numerically together with the water balance equation.The results of the prediction of the radionuclide concentration within the ood period depend on the operation mode of the reservoirs, i.e. from the changes in the water levels at the Hydropower Plants dams.Optimization methods are used to provide choices of the reservoir system operation mode under the water quality criteria.

Modelling analyses of water contamination
The post-accident analyses of water contamination require to combine the consideration of local e ects of constructions designed to diminish radionuclide uxes from the area surrounding the Chornobyl NPP into the Prypiat River with the forecasting of radionuclide dynamics in the Kyiv Reservoir (capacity 3.7 km 3 ) and in the whole Dnieper cascade of six reservoirs.

Floodplain scale
A danger of signi cant Prypiat River radioactive pollution came from a heavily contaminated territory on the left-bank oodplain of the Prypiat River upstream the Chornobyl NPP opposite the town of Prypiat.The results of COASTOX simulations predicted 2] that, in case of this territory ooding, the 90 Sr concentration in the river can increase up to 10000 Bq/m 3 ( gure 1) which was con rmed by the eld measurement in January, 1991, when the oodplain was over own due to an ice jam 3].These values exceed the accepted Maximum Permissible Level { 3700 Bq/m 3 .
The search for e ective measures for such a situation by the method of mathematical modelling revealed that the best one is the creation of a special 10-km dam on the left-side oodplain around the contaminated area 2,3].The dam was constructed before the 1992 spring ood.
As the continuation of this water -protection measure, the project of the construction of a short right-side dam to prevent the ooding of contaminated spots upstream the Yanov Bridge was considered in the recent years.Simulation of the oodplain ow demonstrated that the most dangerous situation causing a large increase in radionuclide concentrations arises during a spring ood with a maximum discharge at 2000 m 3 /sec (25 % probability of exceeding).The comparison of the model results in the case of such a ood for the existing situation (one dam) and in the case of the second dam construction ( gure 1, lower chart ) demonstrates that this supplementary measure could decrease 90 Sr cross-sectionally averaged concentration at the outlet from this territory from 2150 Bq/m 3 to 1250 Bq/m 3 .The spring ood of April-May 1987 was simulated to compare the THREETOX results with the available data of measurements.The simulated distribution of 137 Cs total concentration C + SC s Bq/m 3

at water
The mathematical modelling of radionuclide transport 47 surface on 30 May is depicted in gure 2a.Deviation of the simulated data from the measured ones did not exceed 10 % at the end of the rst month and 31 % at the end of the second month.This is in the range of uncertainties of the measured data.The vertical variability is signi cant mainly for the radionuclide concentration adherent to suspended sediments, as it is shown in gure 2b where the simulated results (SC s ) are presented for the reservoir vertical cross-section.Seasonal and short-term predictions are in reasonable agreement with the measured data for spring oods, rainstorm oods, consequences of radionuclide releases from the Prypiat oodplain, as results of the ice jams in the winters of 1991 and 1993 2,3].To combine analyses of the local e ects with the whole Dnieper river modelling the above presented results of the 2-D simulation of dam constructions on the Prypiat oodplain were followed by the box modelling of the whole Dnieper reservoir cascade ( gure 4).
To estimate the collective dose for the population of Ukraine from the consumption of the Dnieper water during 70 years after the accident (till 2057) the WATOX code was used in the version based on three-month averaged input data.

Conclusions
A set of models of di erent scales of resolution has been developed to simulate near range and long range radionuclide dispersion from the vicinity of the Chornobyl Nuclear Power Plant via the surface water system Prypiat River-Dnieper reservoirs-Dnieper-Bug Estuary.The models have been validated on the basis of the monitoring data sets.The model set was used for the prediction of radionuclide concentration in water for di erent scenarios of hydrological events and radionuclide releases from the contaminated watersheds of the Chornobyl zone.The modelling results were used as a technical background for the construction of a oodplain protection dike in the town of Prypiat.Simulation of the processes of radionuclide-sediment interaction continues to be the most complicated problem of these studies.It is necessary to provide further investigations of the balance between a complexity of the model (e.g.including two-step kinetic processes submodels) and its applicability to the adequate description of radionuclide transport in water bodies.The developed model set can be used in the decision support systems for o -site emergency management.

Figure 1 .Figure 2 .Figure 3 .
Figure 1.Depth-averaged 90 Sr concentration (Bq/m 3 ) during the ood (Q = 2000 m 3 /s) on the Prypiat oodplain upstream the Chornobyl NPP before water protection measures (upper chart) and after the design of two dams.The ow direction is from the left to the right.4.2.River branch scale Bottom traps were constructed to settle down the contaminated sediments in the Prypiat River channel.The traps are regions (length 250-500 m) of increasing the water depth (from 5-6 m to 12-16 m) and the canal width.The e ciency of the traps was simulated by the 2-D vertical VERTOX model.The computed results demonstrate the selectivity of the traps for pollutant deposition depending on radionuclide distribution through the size fractions of the sediment grains. 137Cs is transported mainly by the nest suspended sediments that makes such measures ine cient.Less then 1 % of the horizontal ux of ne sediments was settled down in the traps.4.3.Reservoir scaleThe Prypiat River in uxes into the north-western part of the Kyiv reservoir 30 km downstream the Chornobyl nuclear power plant.The lake-like part of the reservoir is 60 km long and 16 km wide, with maximum depth 15 m.The hydropower plant dam at the southern end of the reservoir is very close to the northern part of Kyiv suburbs.After the Chornobyl accident the bottom deposition of the reservoir was heavily contaminated and contained

Figure 4 .
Figure 4.The simulated 90 Sr distribution in the Dnieper reservoirs during the spring ood with the maximum discharge 2000 m 3 /s for di erent scenarios of measures on the Prypiat oodplain.For most practical problems of radionuclide dispersion evaluation in the Kyiv reservoir the 2-D model also gives reasonable results.In the strati ed Dnieper-Bug estuary the radionuclide transport can be simulated only on the basis of the 3-D model 10].