Waste Biomass Valor (2013) 4:197–211
DOI 10.1007/s12649-012-9156-0

ORIGINAL PAPER

Environmental Impact Assessment of Liquid Waste Ponds
in Uranium Milling Installations
C. P. Naveira-Cotta • E. M. Pontedeiro •R. M. Cotta • J. Su • M. T. van Genuchten

Received: 1 March 2012 / Accepted: 30 July 2012 / Published online: 19 August 2012
Ó Springer Science+Business Media B.V. 2012

Abstract A detailed environmental impact assessment isrequired in many countries when considering the disposalof wastes containing radioactive materials. In this paper we
present hybrid numerical-analytical solutions for the subsurface transport of radioactive contaminant decay chainsthat may be used for such an assessment. The model
involves the advective–dispersive transport of multiple
radionuclide species within separate but coupled saturated
and unsaturated soil domains. The resulting partial differential equations were solved using the Generalized Integral
Transform Technique to yield analytical expressions forthe concentration distributions versus distance, and analytical or numerical solutions as a function of time. The
potential of the hybrid modeling approach is illustrated bymeans of an environmental impact assessment of an uranium milling liquid waste pond near Caetité, Brazil. Calculated radionuclide concentration distributions were forthis purpose used in subsequent radiation dose calculations.
Several scenarios were analyzed, including operational,
human intrusion, farmer/road builder, and leaching scenarios. Numerical results show that for all scenarios

C. P. Naveira-Cotta  R. M. Cotta (&)  M. T. van Genuchten
Mechanical Engineering Department, Center for Analysis andSimulations in Environmental Engineering, CASEE,POLI&COPPE, Universidade Federal do Rio de Janeiro,
Cidade Universitária, Cx. Postal 68503, Rio de Janeiro,
RJ 21945-970, Brazil
e-mail: cotta@mecanica.coppe.ufrj.br
E. M. Pontedeiro  J. Su
Nuclear Engineering Department, POLI&COPPE,
Universidade Federal do Rio de Janeiro,
Cidade Universitária, Cx. Postal 68503,
Rio de Janeiro, RJ 21945-970, Brazil

analyzed the total dose is less than the Brazilian regulatory
dose limit of 0.3 mSv/year above natural background for aperiod of 10,000 years.
Keywords Radionuclide transport Environmental impact assessment  Hybrid methods Integral transforms
List of symbols
Ci Concentration of i-th radionuclide in decay chain[Bq/m3]
C* Filtered concentrations [Bq/m3]
D
Dispersion coefficient [m2/year]fi
Concentration of the i-th radionuclide at the inlet[Bq/m3]hMass transfer coefficient [1/year]i
Subscript, refers to the i-th radionuclide in the decay
chain
Kd Distribution coefficient [m3/kg]
L
Length of the transport domain [m]Nr Total of species in the radionuclide decay chainq
Recharge rate given by the Darcy-Buckingham
equation [m/year]
Ri Retardation factor of i-th radionuclide [-]tTime [year]
V
Pore-water velocity [m/year]x
Distance [m]Greek letters
a Dispersivity [m]k Eigenvalue
h Variable water content in the unsaturated zone [m3/m3]
l Decay rate [1/year]
li Decay rate of the i-th radionuclide [1/year]u Eigenfunction

123

198

Introduction
Regulatory bodies in many countries require a detailed
environmental impact assessment when considering thedisposal of solid and liquid wastes from industrial activities involving radioactive materials [1–3]. The assessment
typically includes an analysis of the transport of radionuclide decay chains through engineered and natural barriers
into underlying and downgradient soil and aquifer systems,
for both design and operation licensing purposes. Industrial
activities producing radioactive materials include nuclear
power production, electricity generation, petroleum andnatural gas production, and mining and milling operations[4–6]. Part of the environmental impact assessment generally pertains to the long-term effects to human health[7, 8] of radionuclides in the geosphere through all possible
pathways (soil, water, air). In order to simulate the transport of radionuclide decay chains in soils, one must solve
the coupled system of transport equations for each species
in the chain for the hydrologic and geochemical conditionsof the specific site being considered.
Several powerful and well-tested analytical and numericalsimulation tools are publicly or commercially available for therequired calculations. While many field problems may requirecomprehensive one- or multi-dimensional numerical models
simulating transient water, heat and/or contaminant transport[9, 10], some transport problems are better addressed using
simplified analytical or quasi-analytical models. As pointed
out in many studies (e.g., [11–15]), analytical models areuseful for initial or approximate analyses of alternative pollution scenarios, conducting sensitivity analyses to investigate
the effects of various parameters or processes on transport,
extrapolating results over large times and spatial scales where
numerical solutions become impractical, serving as screening
models, providing benchmark solutions for more complex
transport processes that cannot be solved analytically, and fortesting more comprehensive numerical solutions of the governing transport equations.
Especially attractive for long-term environmental
impact assessments are classical integral transform methods, which during the past several decades have slowly
evolved into very flexible hybrid numerical-analytical
approaches that offer user-controlled accuracy and more
efficient computational performance for a wide variety ofa priori non-transformable problems [16, 17], including
nonlinear mass transfer and heat and fluid flow applications
related to environmental modeling [18]. The hybrid method
in various forms has been applied over the years to a large
number of environmental problems [18–24]. The approachis particularly well suited for benchmarking purposes
because of its automatic error control features, and byretaining the same characteristics as a purely analyticalsolution. In addition to error control and estimation, an

123

Waste Biomass Valor (2013) 4:197–211

attractive aspect of the method is its straightforward
extension to multi-dimensional situations with only amoderate increase in computational effort as compared to
one-dimensional applications [25]. This is due to the hybrid
nature of the solution process since the analytical part isemployed for all but one independent variable, while thenumerical part is used only for integration of a system ofordinary differential equations (ODE) in one single independent variable, usually the time variable.
In this study we used the Generalized Integral Transform
Technique (GITT) in a hybrid numerical-analytical solutionof the set of advection–dispersion equations describing the
migration of radioactive decay chains in finite porous media,
as implemented within the Mathematica 5.2 platform [26].
The one-dimensional transient formulation is applicable tochains of any prescribed length and for varying retardation(sorption) and decay properties among the radionuclides
involved. The method can be implemented using either prescribed initial concentrations within the repository as part of
the transport domain, or by invoking inlet concentrations
obtained from solutions of the Bateman equations [27]assuming a well-mixed repository releasing its radionuclides
over time [28, 29]. First- or third-type boundary conditions
can be used for the leaching (inlet) boundary, and second- or
third-type conditions for the downgradient (outlet) side of thetransport domain. A simple filtering strategy can further be
employed to homogenize the leaching boundary condition,
whereas the NDSolve routine [26] is used to numerically solve
the transformed ODE system. Each transformed governing
equation for the individual radionuclides can be employedalso separately in a recursive implementation of the numericalsolution of the ODE system. Finally, an integral balance
strategy is employed to enhance convergence behavior. In thisstudy we tested the accuracy of the solution against previous
results by van Genuchten [29] and Lung et al. [30], who used
analytical solutions derived for far more restrictive transportconditions.
The hybrid solution will be applied to four different
scenarios for the migration of radioactive contaminantsfrom liquid waste ponds at the Uranium mining facility of
INB (Industrias Nucleares do Brasil) located near themunicipality of Caetité, Bahia, Brazil. The solution procedure was selected because of its unique advantages foranalysis of long-term contaminant transport problems(which for certain applications may run for tens of thousands of years), and because the approach is one of therecommended methods for radioactive waste safetyassessment by the International Atomic Energy Agency[31, 32]. Calculated radionuclide concentrations are further
used in subsequent radiation dose calculations. Safety ordose calculations for long time periods are essential for the
management of long-lived radioactive wastes [3, 33]. Since
institutional controls may not last much longer than a few

Waste Biomass Valor (2013) 4:197–211

decades after site closure, the wastes are unlikely to remain
safe. Natural leaching of radionuclides from the repository
into the adjacent biosphere is almost certainly to occur
within a few centuries or even earlier. Since mining wastes
containing low-level natural radioactive materials are typically deposited in industrial landfills, regulations in Brazil
and elsewhere require a safety assessment of the disposal
facility using a leaching and off-site contaminant transportscenario.

199

oCi
oCi
o2 C iþV¼ D 2  li Ri Ci ðx; tÞ þ li1 Ri1 Ci1 ðx; tÞ;otoxox
0  x  L; t [ 0; i ¼ 1; . . .; NrRi

ð1ÞThe associated initial and boundary conditions are
Ci ðx; 0Þ ¼ 0;

0  x  L;

i ¼ 1; ::; NroCi
þ VCi ¼ Vfi ðtÞ orox
i ¼ 1; . . .; NrD
Problem Formulation and Solution Methodology
Our conceptual model for the site involves a liquid waste
pond, an underlying unsaturated soil zone through whichwater flows vertically downward, a permanently saturated
granular aquifer through which water flows laterally, and acorresponding mixing zone where water from the unsaturated zone mixes with the saturated zone (Fig. 1). The
recharge rate through the unsaturated zone into groundwater is defined by the long-term hydrologic characteristicsof the site, thus disregarding operational conditions and therelatively short time period needed to fill the pond. Lateralflow through the aquifer is directly obtained from thehydraulic gradient and the saturated hydraulic conductivity.
The mixing process of unsaturated recharge water and its
dissolved constituents with water laterally moving throughthe saturated zone was modeled using the EPA mixing
zone model [34] based on a penetration depth concept.
To solve the radionuclide transport problem with the hybrid
method [16–18], one-dimensional transient advection–dispersion equations are assumed for each of the Nr species in theradionuclide decay chain. The solution procedure is firstdescribed for the saturated medium. The governing transport
equations for the individual species are given by:DoCi
þ hCi ¼ 0;
ox

ð2ÞCi ¼ fi ðtÞ;x ¼ 0;

t [ 0;
ð3Þ

x ¼ L;

t [ 0;i ¼ 1; . . .; Nrð4Þ

where fi(t) is the concentration of the i-th radionuclide at the
inlet. The boundary conditions at x = 0 can be of a third or
first type, while the exit boundary condition can be of a third or
second type depending upon the selected value for the masstransfer coefficient h. Before applying the Generalized
Integral Transform Technique (GITT), an analytical filtering
approach is used first to homogenize the boundary conditions
at x = 0, which will improve convergence of the proposed
eigenfunction expansions. While several filtering methods areavailable for this purpose [35, 36], here we use a very simple
filter that employs the boundary source term concentration, as
follows:
Ci ðx; tÞ ¼ Ci ðx; tÞ þ fi ðtÞ;i ¼ 1; . . .; Nr0  x  L;t [ 0;

ð5Þ

The transient advection–dispersion problem withhomogeneous boundaries can now be rewritten in the form
oCi
oC 
o2 C 
þ V i ¼ D 2i  li Ri Ci þ li1 Ri1 Ci1þ gi ðtÞ;otoxox
0  x  L; t [ 0; i ¼ 1; . . .; NrRið6Þdfi ðtÞ
 li Ri fi ðtÞ þ li1 Ri1 fi1 ðtÞ;dt
i ¼ 1; . . .; Nrgi ðtÞ ¼ Ri

Ci ðx; 0Þ ¼ fi ð0Þ;

0  x  L;oCi
þ VCi ¼ 0; x ¼ 0;
ox
t [ 0; i ¼ 1; . . .; NrDt [ 0;
ð7Þ

i ¼ 1; . . .; Nr

ð8Þ

oCi
¼ 0;
ox

ð9Þ

x ¼ L;

To illustrate the remainder of the solution procedure, we use
the above third type inlet and second type exit boundary
conditions. The following auxiliary eigenvalue problem
results:
Fig. 1 Schematic of the adopted conceptual model for verticalinfiltration into residual unsaturated soil coupled to lateral flow in ahorizontal granular aquiferD

d2 un ðxÞþ k2n un ðxÞ ¼ 0;
dx2

0\x\L

ð10Þ

123

200

Waste Biomass Valor (2013) 4:197–211

dun ðxÞ
þ Vun ðxÞ ¼ 0;
dx
x¼L

D

dun ðxÞ
¼ 0;
dx

x ¼ 0;

and with the transformed initial conditions given by
ð11Þ

h

ð0Þ ¼ fi ð0Þfn ¼ Ci;nThis Sturm-Louiville problem is readily solved for the
eigenfunctions and the related eigenvalues and norms: V sin pknﬃﬃﬃ x
kn
D
pﬃﬃﬃﬃ
un ðxÞ ¼ cos pﬃﬃﬃﬃ x þ
ð12Þ
kn D
D



 V cos pknﬃﬃﬃ L
kn
D
pﬃﬃﬃﬃ
kn sin pﬃﬃﬃﬃ L ¼
ð13Þ
D
D

V  V cosðpLlﬃﬃDﬃÞ þ
n

i
pﬃﬃﬃﬃ
D sinðpLlﬃﬃDﬃÞln fi ð0Þ
n

l2n

ð20Þ
The above transformed system, while amenable ofanalytical treatment using the corresponding coefficientmatrix eigensystem analysis, is most flexibly handled using
the NDSolve routine of Mathematica 5.2 [26] to generalize
the solution procedure to unsaturated systems as will be
discussed later. Once the transformed concentrations are
obtained numerically using Mathematica, an interpolating




2kn


V 2 þ Dk2n sin pﬃﬃﬃﬃ L
2Lkn V 2 þ Dk2n
2kn
D
ZL
pﬃﬃﬃﬃ
þ
2Vkn  2Vkn cos pﬃﬃﬃﬃ L þDDD
normn ¼ u2n ðxÞdx ¼
4k3n






ð14Þ

0

The integral transform pair is defined next in terms of an
eigenfunction expansion for the filtered concentrations,Ci ðx; tÞ, by taking the eigenfunctions, un ðxÞ, as the basis
for the expansions:

ðtÞ ¼
Ci;n

ZL

un ðxÞCi ðx; tÞdx!

Transform

ð15Þ

0Ci ðx; tÞ ¼

1
X

un ðxÞ C ðtÞ
normn i;nn¼1!

Inversefunction is automatically generated to provide the time
behavior in a continuous fashion. The inverse can then berecalled to explicitly provide the spatial dependency of theconcentration fields.
For the more general situation of an unsaturated soil, thegoverning transport equations are given by

ð16Þ

By implementing the above filtered problem using
0 un ðxÞdx, the space variable x is eliminated, leading to
the following transformed problem involving an infinite
system of ordinary differential equations in time:o½ðh þ qKdi ÞCi  oðqCi Þ
o
oCi
þ
¼ ðhD
Þ
ot
ox
ox
ox
 li ðh þ qKdi ÞCi ðx; tÞ þ li1 ðh þ qKdi1 ÞCi1 ðx; tÞ;
0  x  L; t [ 0; i ¼ 1; . . .; Nr
ð21Þ

RL

Ri


1
X
dCi;nðtÞþVan;m Ci;m
ðtÞ
dt
m¼1



ðtÞ  li Ri Ci;nðtÞ þ li1 Ri1 Ci1;nðtÞ þ gi ðtÞfn ;¼ k2n Ci;nt [ 0; i ¼ 1; . . .; Nr ; n ¼ 1; 2; . . .

ð17Þ
where

For long-term simulations such as those considered here
(i.e., in the order of thousands of years), a time invariable
recharge rate can be assumed, in which case the watercontent will be a function of position only. The formulation
then simplifies tooCioCi
o2 C iþq¼ aq 2  li ½hðxÞ
otoxoxþ qKdi Ci ðx; tÞ þ li1 ½hðxÞ þ qKdi1 Ci1 ðx; tÞ;
0  x  L; t [ 0; i ¼ 1; . . .; Nr½hðxÞ þ qKdi ð22Þ
an;m ¼

1Nm

ZL

un ðxÞu0m ðxÞdx

ð18Þ

0

fn ¼

ZL

un ðxÞdx

0

123

ð19Þ

where hðxÞ is used to indicate that the water content is afunction of the spatial coordinate, x. The corresponding
initial and boundary conditions are
Ci ðx; 0Þ ¼ 0;

0  x  L;

i ¼ 1; . . .; Nrð23ÞWaste Biomass Valor (2013) 4:197–211

Ci ¼ fi ðtÞ;x ¼ 0;

t [ 0;
ð24Þ

a

oCi
þ h Ci ¼ 0;
ox

x ¼ L;

t [ 0;i ¼ 1; . . .; Nrð25Þ

in which we used the definition of the dispersion coefficientas the product of the dispersivity, a, and the recharge porewater velocity, V (i.e., D = aV).
The integral transforms of the system of equations given by
Eqs. (22) to (25) are very similar to those of the saturated
medium, and as such will not be repeated here. The only
difference is that the integral transformation in this casegenerates other coupling terms, including the coupling matrix
for the transient terms of the transformed concentrations. For
any arbitrary water content distribution, h(x), we may adopt asemi-analytical integration procedure to generate the coefficient matrix with substantial savings in computational effort[37]. Also, the input profile for the water content depends on
the solution of the nonlinear flow problem for the unsaturated
zone as given by the Richards equation for steady-state flow.
That initial profile was obtained with the HYDRUS-1D software package [38] based on a numerical solution of the
Richards equation [39] for the soil profile being consideredassuming a steady infiltration rate as dictated by the long-time
recharge rate in the study area.

Numerical Verification
The above hybrid numerical-analytical procedure was
implemented in the Mathematica 5.2 system [26], including all of the symbolic calculations of the solution process.
The numerical accuracy of the resulting algorithm wastested against two codes previously developed for radionuclide decay chains: The LBL (Lawrence Berkeley Laboratory) computer code of Lung et al. [30] and the CHAIN
code of van Genuchten [29]. Solutions employed in thesecodes were derived using Laplace transform techniques asapplied to semi-infinite media.
Comparison with the Lawrence Berkeley Laboratory
Code of Lung et al. [30]
The LBL code deals with purely advective transport (no dispersion) of radionuclides in saturated semi-infinite media toallow for explicit analytical inversion of the Laplace transforms. We consider here the three element LBL chain
example formed by 234U, 230Th, and 226Ra because of its
similarity with the chain to be analyzed later for the Caetité

site. The following parameter values were used: D = 50 m2/
year, h = 3 9 10-3, V = 1 m/year, R1 = 120, R2 = 1,500,
R3 = 300, l1 = 2.806 9 10-6 year-1, l2 = 8.664 9 10-6
year-1, and l3 = 4.332 9 10-4 year-1. The normalized
concentrations at boundary x = 0 are given by f1 = 1, f2 = 1,
and f3 = 10. Numerical results were obtained for the sametime values (10 and 1,000 years) as used by Lung et al.
(1987), while we used a domain length L = 200 m to avoid
any effects of the exit boundary. Figure 2 compares calculated radionuclide concentration distributions for the twosolutions at t = 1,000 years. Accounting for some loss inprecision from directly reading the graphs in the LBL report[30], the agreement is very good. While the GITT results
were obtained for system truncation orders N \ 150,
graphical convergence was attained much earlier. Figure 2
shows relatively rapid advancement of the contamination
plumes for 234U and 226Ra, which have much lower retardation coefficients in comparison with 230Th. Table 1
illustrates the convergence behavior of the eigenfunction
expansion for 234U at t = 1,000 years. The selected truncation orders were N = 50, 75, 100, 125 and 150 for thewhole solution domain. Convergence is overall very satisfactory, yielding five or six significant digits in the results.
The GITT solution hence was found to be precise in predicting concentrations that are orders-of-magnitude smaller than those at the inlet. This is important since utilities areusually required to satisfy environmental health criteria at
the edge of their property, which in general involve much
lower concentrations than at the source of contamination.
Comparison with the CHAIN Code of van Genuchten[29]
The CHAIN code considers the advective–dispersive
transport of a radionuclide decay chain of up to four species. For the present test we employed the four element
1
234

U

230

Th

226

Ra

0.8

C/Co

oCi
þ Ci ¼ fi ðtÞ orox
i ¼ 1; . . .; Nr

a

201

0.6
0.4
0.2

20

40

60

80

100

x (m)
Fig. 2 Comparison of concentration distributions at t = 1,000 years
for three radionuclides (234U, 230Th, and 226Ra) as obtained in thisstudy using GITT (solid lines) and by Lung et al. [30] using Laplacetransforms (symbols)

123

202

Waste Biomass Valor (2013) 4:197–211

Table 1 Convergence of the
GITT solution for 234U
concentration distribution at
t = 1,000 years for the LBL test
case of Lung et al. [30]involving 234U, 230Th and 226Ra(N is the truncation limit)

Table 2 Convergence of the
GITT solution for the relatively
concentration distribution of264
Ra at t = 10,000 years forthe test case of van Genuchten[29] involving 238Pu, 234U,
230
Th, and 226Ra (N is thetruncation order)x [m]N = 50N = 75N = 100

N = 125

N = 15010.980879

0.980935

0.980955

0.980962

0.980963

100.797314

0.797306

0.797301

0.797300

0.797300200.585814

0.585807

0.585811

0.585810

0.58581030

0.393687

0.393696

0.393693

0.393694

0.393694400.240582

0.2405780.240579

0.240579

0.240579

60

0.0663583

0.0663609

0.0663615

0.0663618

0.0663619

80

0.0119769

0.0119747

0.0119757

0.0119756

0.0119755

100

0.00139205

0.00139003

0.00139092

0.00139085

0.00139068

x [m]N = 50N = 100

N = 150

N = 200

N = 250CHAIN

101.2737 E-6

1.5689 E-6

1.5630 E-6

1.5629 E-6

1.5663 E-6

1.5111 E-6

203.5447 E-5

3.5923 E-5

3.5927 E-5

3.5927 E-5

3.5934 E-5

3.5898 E-530

8.9891 E-5

9.0126 E-5

9.0136 E-5

9.0135 E-5

9.0138 E-5

9.0103 E-5401.4340 E-4

1.4347 E-4

1.4346 E-4

1.4346 E-4

1.4346 E-4

1.4343 E-450

1.9501 E-4

1.9494 E-4

1.9494 E-4

1.9494 E-4

1.9494 E-4

1.9491 E-460

2.4041 E-4

2.4031 E-4

2.4031 E-4

2.4031 E-4

2.4031 E-4

2.4029 E-4802.6106 E-4

2.6112 E-4

2.6112 E-4

2.6112 E-4

2.6112 E-4

2.6109 E-4100

2.4736 E-4

2.4739 E-4

2.4739 E-4

2.4739 E-4

2.4739 E-4

2.4736 E-4

chain 238Pu, 234U, 230Th, and 226Ra as in [29] using inlet
boundary conditions obtained from the exact solution of
the Bateman equations for radioactive decay in andleaching from a waste site. We used the following
parameter values: D = 10 m2/year, V = 100 m/year, R1 =10,000, R2 = 14,000, R3 = 50,000, R4 = 500, l1 = 7.9 9
10-3 1/year, l2 = 2.8 9 10-6 1/year, l3 = 8.7 9 10-6 1/
year and l4 = 4.3 9 10-4 1/year. Numerical results wereobtained at t = 10,000 years as provided also in [29],Relative Concentration1

t = 10,000 y

0.1
234

U
Th
226
Ra
230

0.01

while a domain length, L, of 300 m was sufficient to honor
the semi-infinite medium assumption in our calculations.
Table 2 illustrates the convergence behavior for 226Ra forthis strongly advective case. The various columns show
GITT results for different truncation orders using increasing values of N. Also shown are the CHAIN results (lastcolumn) as obtained directly from the CHAIN software
output.
The results in Table 2 indicate that the convergence rate of
the eigenfunction expansions is again very good, yielding four
digit agreement with the CHAIN benchmark solution, except
near the inlet boundary. Figure 3 shows the concentrations of
the four elements at t = 10,000 years. Notice that the 238Pu
concentrations are not visible within the invoked log concentration range, and that 226Ra moves by far the fastest
through the system because of its relatively low retardation
factor.

Application to Uranium Waste Disposal Site
1E-3

1E-4
1

10

100

Distance, x (m)Fig. 3 Calculated concentration distributions versus distance for fourradionuclides, 238Pu, 234U (top curve), 230Th (middle curve), 226Ra(bottom curve) at t = 10,000 years as obtained with GITT for the test
case of van Genuchten [29]123

Our code using the GITT approach is now demonstrated
for the environmental impact assessment of an industrial
installation that deals with liquid and solid waste containing natural radioactive contaminants. The case study concerns possible pollution from an uranium mining andmilling facility operated by INB (Indústrias Nucleares do
Brasil) in Caetité (Bahia, Brazil). Liquid wastes resulting
from the uranium extraction process were deposited in
ponds within the mining installation. Each cell of the ponds

Waste Biomass Valor (2013) 4:197–211

203

Table 3 Scenarios considered for the migration of radionuclides from an uranium mining waste disposal pond
Case

Zones considered

L1/La2 [m]Infiltrationa

Water content[m3/m3]

Kd0 s

1

Unsaturated/saturated

6.9/1,150

q = P - ET/Darcy

h = 0.308/0.475

Table 42

Unsaturated/saturated

4.2/1,150

q = P - ET/Darcy

h = 0.308/0.475

Table 4

3

Unsaturated/saturated

6.9/1,150

q = P/Darcy

h = 0.338/0.475

Table 4

4

Unsaturated/saturated

6.9/1,150

q = P/Darcy

h = 0.338/0.475

Claya/Table 4a

L1 = depth of unsaturated zone, L2 = length of saturated zone, P = annual precipitation rate, ET = annual evapotranspiration rate,
Darcy = flow rate estimated using Darcy’s Law, Clay = Kd0 s of a clay soil

has an estimated operation period of 4 years. The Brazilian Nuclear Energy Commission (CNEN) required ananalysis of the long-term transport of radioactive contaminants from the ponds into the subsurface, including
an assessment of the total radioactive doses resulting frompossible exposure to the radionuclides. Details of theanalysis are given elsewhere [40], including model
development, code construction, identification of flow andtransport properties, and dose calculations. Here we focus
on the contaminant transport simulations of four cases(Table 3) representing the most important possibilities foroff-site migration.
The base case (Case 1) assumes a 6.9-m deep unsaturated zone below the pond, with recharge estimated from ahydrologic water balance of the Caetité region, and withwater contents calculated with the HYDRUS-1D software
package [38] using soil hydraulic properties [41] estimatedfrom collected soil samples [40, 42]. The unsaturated zone
was coupled to a horizontal granular aquifer of 1,150 m
length. The flow rate through the horizontal aquifer was
estimated from the hydraulic gradient and the hydraulicconductivity using Darcy’s law.
The other three cases (2, 3, 4) considered here are relatively simple perturbations of the base case. Of these,
Case 2 assumes a shallower unsaturated zone (4.2 m),
typical of wet seasons in the area, while keeping theremaining parameters unchanged. Case 3 considers a muchhigher recharge rate through the unsaturated zone by considering only precipitation and neglecting evapotranspiration. Finally, Case 4 illustrates the importance of accurately
estimating the partitioning coefficients, Kd0 s, of the radionuclides. In this case reported literature values for clay
soils (stronger retardation) are used rather than the measured values listed in Table 4. Table 4 also gives estimatesof the radioactive inventory of the simulated waste pond
[40].
Calculated concentration distributions for the various
cases are presented in Figs. 4 through 8. In each plot
showing distributions versus distance, x, we employed
black lines for 238U and 234U (graphs for these two radionuclides always coincided for the time periods considered

Table 4 Partitioning coefficients, Kd [ml/g], and radioactive inventories, A [Bq]Element

Kd [ml/g]

A [Bq]

238

2,985.6

5.34 9 1011

230

44,279.1

1.68 9 1010

226

198.315.51 9 1010

210

6,772.8

6.16 9 1011

234

U/

Th
Ra
Pb

U

in this study), blue lines for 230Th, red lines for 226Ra, andgreen lines for 210Pb.
Results for the base case (Case 1) were compared firstwith a very conservative scenario, not further referred to in
Table 3, which completely neglects the unsaturated zone
below the pond. Figure 4 show for this scenario the calculated concentration distributions along the horizontal
aquifer after 500 years. Clearly, including the unsaturated
zone considerably slows down the transport process, as
expected, but also significantly lowers the maximum concentrations. These effects are most pronounced for lead,radium and thorium, which for our soils had much larger
values of the retardation factors as compared to uranium.
Figures 5 shows for Case 1 the radionuclide concentration distributions in the aquifer at times t = 1,000,
2,000, 3,000, 4,000, and 10,000 years. Results again
illustrate the relatively rapid migration of uranium (lower
Kd values), but also very low concentrations characteristicof thorium transport. At the same time, radium and leadbecome more abundant at later times despite their relatively high retardation factors.
Figure 6 provides an interesting comparison between
the base case and Case 2 (shallower unsaturated zone). The
top plots show for the two cases concentration profiles of
the radionuclides within the unsaturated zone at
t = 100 years, while the plots at the bottom illustrate theirmigration along the horizontal aquifer at t = 500 years.
Results again indicate that a thicker unsaturated zone (6.9
vs. 4.2 m) only minimally decreases the transport rate in
the aquifer, but significantly reduces the maximum concentrations in the aquifer at similar total migration times.

123

204

Waste Biomass Valor (2013) 4:197–211
t = 500 y1

t = 500 y
0.1
0.01

238

U, 234U
230
Th
226
Ra
210
Pb

0.01

238

U,
Th
226
Ra
210
Pb

234

U

230

1E-3

Ci (Bq/L)

Ci (Bq/L)

0.1

1E-4
1E-5

1E-3
1E-6
1E-7

1E-4
1

10

100

1

1000

10

x (m)

100

1000

x (m)

Fig. 4 Calculated radionuclide concentrations versus distance x in the horizontal aquifer for base Case without (left) and with (right) accounting
for the unsaturated zone

0.00060.5

1,000 y2,000 y3,000 y4,000 y10,000 y

0.3

0.2

0.00041,000 y2,000 y3,000 y4,000 y10,000 y

0.0003230

238

U (Bq/L)

0.4

Th (Bq/L)

0.00050.00020.1

0.00010.0

0.0000
0

200

400

600

800

1000

1200

0

200

400

x (m)

600

800

1000

1200

x (m)

0.012

0.025

226

0.020

0.008

0.006

Pb

1,000 y2,000 y3,000 y4,000 y10,000 y

1,000 y2,000 y3,000 y4,000 y10,000 y

0.015

210

Ra (Bq/L)

0.010

0.010

0.004
0.005

0.002

0.000

0.000

0

200

400

600

x (m)

800

1000

1200

0

200

400

600

800

1000

1200

x (m)

Fig. 5 Calculated radionuclide concentrations versus distance x in the horizontal aquifer for base Case 1 (with unsaturated zone). Shown are the
distributions for 238U (top left), 230Th (top right), 226Ra (bottom left) and 210Pb (bottom right)

123

Waste Biomass Valor (2013) 4:197–211

205

t = 100 y

t = 100 y

1

1

238

0.01

238

0.1

230

234

U, U
Th
226
Ra
210
Pb
230

Ci (Bq/L)

Ci (Bq/L)

234

U, U
Th
226
Ra
210
Pb

0.1

1E-3

0.01

1E-3

1E-4
0.01

1E-4

0.1

1

0.1

1

1

x (m)

x (m)

t = 500 y

t = 500 y

1

0.1
0.1

0.01

238

U,
Th
226
Ra
210
Pb

234

U

238

230

Ci (Bq/L)

Ci (Bq/L)

1E-3
1E-4

234

U, U
Th
226
Ra
210
Pb

0.01

230

1E-3

1E-4

1E-5

1E-5

1E-6

1E-6

1E-7
1

10

100

1000

1

10

x (m)

100

1000

x (m)

Fig. 6 Calculated radionuclide concentrations versus depth in the unsaturated zone (top) and the horizontal aquifer for base Case 1 (6.9 m deep
vadose zone, plots on the left) and Case 2 (4.2 m vadose zone; plots on the right)

The effect of a higher recharge rate (Case 3; 1.104 m/year)
relative to Case 1 (0.471 mm/year) is shown in Fig. 7. The top
plots show concentrations in the unsaturated zone at
t = 500 years, while the plots at the bottom show radionuclide transport in the aquifer at t = 10,000 years. The higher
recharge rate for Case 3 (plots on the right) causes more rapidtransport of the contaminants, which now move much further
into to aquifer at similar elapsed times, and with higher
concentrations.
Figure 8 next compares Case 1 (plots on the left side)with Case 4 for no evapotranspiration (higher rechargerate) and higher Kd values (plots on the right). The top
plots present 238U and 234U concentrations in the unsaturated zone for times up to 1,000 years, while the bottom
plots show the radionuclide distributions in the aquifer at
t = 500 years. The higher Kd values typical of clay soils

lead to much higher retardation factors and hence tomarkedly different distributions for the two cases. For
example, the U concentrations after 1,000 years in theunsaturated zone are much lower for Case 4 and, as
expected, barely reach the horizontal aquifer after
500 years.

Performance Assessment
The main results of the performance and dose calculations
are presented in this section. The primary objective of aperformance or safety assessment is to provide evidence
that the public and the environment are protected as muchas possible in the future. This is possible by carrying out
simulations and using the results to judge if the proposed

123

206

Waste Biomass Valor (2013) 4:197–211
t = 500 y
t = 500 y
1
1

238

238

234

U, U
230
Th
226
Ra
210
Pb

0.01

Ci (Bq/L)

Ci (Bq/L)

0.1

1E-3

234

U, U
Th
226
Ra
210
Pb

0.1

230

0.01

1E-3

1E-4

1E-4
0.1

0.1

1

1

x (m)

x (m)
1

t = 10,000 y
238

238

234

U, U
Th
226
Ra
210
Pb

1E-3

0.1

Ci (Bq/L)

0.01

234

U, U
Th
226
Ra
210
Pb
230

230

0.1

Ci (Bq/L)

t = 10,000 y

1

0.01

1E-3

1E-4
1

10

100

1000

x (m)

1E-4
1

10

100

1000

x (m)

Fig. 7 Calculated radionuclide concentrations versus distance, x: in the unsaturated zone (plots on top) and horizontal aquifer (bottom) for Case1
having a low recharge rate (left) and Case 3 having a high recharge rate (right)

disposal site meets regulatory standards for long-term
protection of human health and the environment. A pathway analysis with different scenarios permits a systematicway of evaluating potential routes by which the public
could be exposed to radiation [6, 8, 43]. Because of theuncertain nature of the long-term integrity of disposal
systems, as well as the unpredictability of human endeavors, all post-closure scenarios are by their very nature
hypothetical.
The safety assessment assumes that the dose projection
for members of the public cannot exceed 0.3 mSv/year
above natural background. Assuming a risk factor of 5 %per Sievert (Sv) as recommended by ICRP [44], this leads
to an assessed radiological risk of fatality to a potentially
exposed group of 5 9 10-5 (i.e., 5 in one hundred thousand per year). The acceptable fraction, termed the dose or
risk constraint, is generally determined by the regulatory123

agency involved. For example, Brazilian regulations use adose constraint of 0.3 mSv/y above natural background,which used in this study.
Consistent with calculations in the previous section, ourperformance assessment considered vertical flow throughthe waste site (pond) and the 6.9-m deep unsaturated zone,
and subsequently horizontal flow through the underlying
aquifer until water was intercepted by a well or discharged
to a stream outside the URA. We assumed that the well or
stream at that location was the only source of water for
local residents for human consumption, irrigation or for use
by farm animals. The radiological dose was calculated foran exposed critical group of individuals receiving thehighest dose or risk. Data regarding human behavior werebased on the dietary habits of adults from that group (e.g.,water and food consumption, breathing rate, and number ofhours being indoor or outdoor). Residents received doses

Waste Biomass Valor (2013) 4:197–211

207

t = 500 y

t = 500 y
1

1

238

234

U, U
Th
226
Ra
210
Pb
230

238

0.1

234

U, U
Th
226
Ra
210
Pb
230

0.01

Ci (Bq/L)

Ci (Bq/L)

0.1

1E-3

0.01

1E-3

1E-4

1E-4

0.1

1

x (m)

0.1

1

x (m)

Fig. 8 Calculated 238U radionuclide concentrations (Bq/L) versus depth, x, and time, t (top plots) and versus depth at t = 500 y (bottom) in theunsaturated zone for Case 1 (left plots) and Case 4 having a higher recharge rate and higher Kd values (right)

from external radiation, ingestion of contaminated waterand food, and inhalation of airborne radionuclides transported from the waste site or suspended in runoff water
following irrigation. Following previous studies [6, 32], thefollowing 10 specific processes in the biosphere wereconsidered: (a) ingestion of water from the well, (b) irrigation, (c) re-suspension or particulates and inhalation,(d) external radiation exposure, (e) consumption of homegrown produce, (f) consumption of contaminated meat,(g) ingestion of contaminated milk, (h) accidental ingestionof contaminated soil, (i) inhalation of radon and decay
products from soil, and (j) contact with surface water,
including transfer to humans. However, we did consider theeating habits and current agricultural activities of the critical group for the dose calculations. For this reason we did
not include exposure through the ingestion of fish since thelocal population has an abundant supply of fish from outside the region. The dose calculations also assumed current
percentages of the local versus outside production of milk

and meat. The diet of the critical group was assumed to be
based on manioc, rice, beans, potatoes and meat (bothchicken and beef), while their fruits consisted of banana,orange and lemon, and their vegetable of carrot, bitter
melon (chuchu), beet and some lettuce. Ingestion, inhalation and external dose conversion factors (Sv Bq-1), as
well as transfer coefficients and concentration factors weretaken from IAEA [43].
We considered four scenarios: a Leaching Scenario(CEN1), a Human Intrusion Scenario (CEN2), a FarmerDweller Intrusion Scenario (CEN3) and a Road Builder
Scenario (CEN4). Radionuclide concentrations calculated
in the previous section were used for both the dose
calculations of Leaching Scenario CEN1. The dose calculations presented here were carried out using methodologies developed by IAEA [43] as implemented byPontedeiro [3].
For Leaching Scenario CEN1 we assumed that radionuclides migrated through an 1,850 m long horizontal

123

208

Waste Biomass Valor (2013) 4:197–211
0.020

Water Ingestion Dose (mSv/y)

0.06

Dose (mSv/y)

0.015

Particulate
Vegetables
Meat
Milk
Water

0.010

0.005

2000

4000

238

234

U, U
Th
226
Ra
210
Pb

0.05

230

0.04

0.03

0.02

0.01

6000

1000

2000

3000

t (years)

4000

5000

6000

7000

t (years)

0.02000.10

Total Dose (mSv/y)0.0175

Dose (mSv/y)

0.0150External
Inhalation
Radon
Ingestion

0.01250.01000.0075

0.08

0.06

0.04

0.00500.02
0.0025

0

1000

2000

3000

4000

5000

6000

7000

0

2000

4000

t (years)

8000

10000

t (years)Fig. 10 Dose calculations for Human Intrusion Scenario CEN2:
Doses due to ingestion of contaminated water (top) and total dose(bottom)

0.020

0.015

0.25
0.010

0.20
0.005

0

1000

2000

3000

4000

5000

6000

7000

t (years)Fig. 9 Dose calculations for Leaching Scenario CEN1: Ingestion bydifferent paths (top), exposure to different pathways (middle) andtotal dose (bottom)

aquifer until reaching a small stream or well outside the
URA boundary. Contaminated water was then used for
irrigation and consumption by cattle and the farmer and his

123

Total Dose (mSv/y)Total Dose (mSv/y)

6000

0.15

0.10

0.05

0

50

100

150

200

250

t (years)Fig. 11 Total dose for Farmer-Dweller Scenario CEN3300

Waste Biomass Valor (2013) 4:197–211

209

Total Dose (mSv/y)

0.005

0.004

0.003

0.002

0.001

0

50

100

150

200

250

300

t (years)Fig. 12 Total dose for Road Builder Scenario CEN4

family. Dose calculations for CEN1 are presented in Fig. 9.
Results indicate that the dose due to the ingestion of water
was dominant as compared to the doses resulting from theother exposure pathways. The maximum total dose
remained less than 0.0185 mSv/year during the first
10,000 years after closure of the repository. This maximum
is far less than the Brazilian regulatory limit of 0.3 mSv/
year above natural background. It is important to note thatthe total dose of the Leaching Scenario was almost zerofor the first several thousands of years, being only 1.14 9
10-7 mSv/year at t = 2,000 years.
For the dose calculations of the Human Intrusion Scenario (CEN2), shown in Fig. 10, we used calculated
radionuclide concentrations leaving the unsaturated zone
below the original waste layer at a depth of 6.9 m. We
assumed that intrusion occurred in the URA restricted area
after institutional control ended, and that water was
removed from the aquifer exactly under the embankmentfor ingestion. This water for ingestion was the only exposure pathway considered for CEN2. Doses of the five radionuclides resulting from this ingestion are presented in
the top plot of Fig. 10. Notice that the 238U and 234U peaks
occurred at about 1,000 years, while the dose due toPb-210 increased during the entire 10,000 years consideredhere. The total dose of the Human Intrusion Scenario isshown in the bottom plot. The plot indicates a first peak
due to 238U and 234U leaching, followed by increases later
due to 210Pb. The total dose value was 0.1007 mSv/year at
10,000 years, again below the Brazilian regulatory limit of
0.3 mSv/year above natural background.
For the Farmer-Dweller Intrusion Scenario (CEN3) we
considered the presence in the area of a small rural settlement. The scenario assumes some dilution of theradioactive waste with nonradioactive materials (localsoil). For the calculations we assumed that 10 % of thewaste was brought to the surface, and that the average
specific activity of all material used by the residents wasequal to the total activity of the waste divided by the total
mass of the material used. For this scenario we included the
10 exposure pathways [7] listed earlier, except any ingestion of water by a well and exposure to surface water,
including transfer to fish and humans. Figure 11 shows thatthe total dose for the CEN3 Framer-Dweller Scenario
declined rapidly from a relatively high initial value of0.28 mSv/year to very small values later.
Finally, for the Road Builders Scenario (CEN4) we
considered inadvertent intrusion, and access to or removal
of waste at or near the surface. The contaminated soilpresumably was used for the construction of a road throughthe property. The road was assumed to be constructed at an
average speed of 20 km per year, with the workers being
present 8 h per day, 20 days per month. The contaminatedmaterial itself was assumed to be spread along a pathequivalent to 1 km, with the road builders working with thecontaminated material for a total of 90 h per year working.
Expose pathways involved (a) external exposure to contaminated soil, (b) dust inhalation, and (c) accidental
ingestion of contaminated soil. As shown in Fig. 12, theinitial value of the total dose was calculated to be
0.004 mSv/year. This is far below the 0.3 mSv/year Brazilian regulatory limit. The total dose again declined very
quickly with time, similarly as for scenario CEN3.

Concluding Remarks
Our results illustrate the considerable power of hybrid
numerical-analytical solutions for evaluating relatively
complex subsurface contaminant transport problems, inthis case the long-term transport of radioactive contaminant
decay chains in a coupled saturated and unsaturated soilaquifer system, followed in this example by a risk assessment study. The generalized integral transform technique(GITT) produced analytical expressions for the concentration distributions versus distance, and analytical or
numerical estimates of the concentration as a function oftime. Especially attractive features of the hybrid method
for long-time environmental assessment studies are its
user-controlled accuracy and very efficient computational
performance for non-transformable problems such as
nonlinear water flow and contaminant transport formulations. Also, the ability to accurately solve the equations forvery small concentrations (compared to the initial concentration) is an important attribute of the GITT method. In
addition to straightforward error control and estimation, an
attractive aspect of the method is its easy extension tomulti-dimensional situations. This is due to the hybrid
nature of the solution process since the analytical part is

123

210

employed for all but one independent variable, while thenumerical part is used only for integration of a system ofordinary differential equations in one single independent
variable, usually the time variable. As such, the approach
can be extended in a relatively straightforward manner also
to non-ideal transport processes such as contaminanttransport in structured (fractured or macroporous) media[18], with concomitant problems of preferential flow. We
finally note that the GITT approach is able to solve serial
reactions in which the sequential reaction rate is equal tothe previous reaction rate, which is a limiting feature that
plagues most or all fully analytical solutions for sequential
decay chain reactions [45]. In fact, since the reaction terms
are not here included in the auxiliary eigenvalue problems,such singularities are not carried into the analytical portionof the solution, while the numerical solution of the transformed equations proceeds as normal. We conclude that themethodology followed in this paper is well suited for
simulations of a broad range of especially long-term contaminant transport problems in soils and groundwater.

References
1. CNEN: Licensing of Radioactive Installations, Report No.
CNEN-NE-6.02, Brazilian Nuclear Energy Commission, CNEN,
Rio de Janeiro, Brazil (1998)
2. IAEA: Extent of environmental contamination by naturally occurring radioactive material (NORM) and technological options for
mitigation, Technical Report Series No. 419, Int. Atomic Energy
Agency (IAEA), Vienna (2003) (available at http://www-pub.
iaea.org/MTCD/publications/Pubdetails.asp?pubId=6789)
3. Pontedeiro, E.M.: Avaliação de Modelos de Impacto Ambiental
para Deposição de Resı́duos Sólidos Contendo Radionuclı́deos
Naturais em Instalações Minero-Industriais, DSc Thesis, PEM/
COPPE/UFRJ, Brazil (2006)
4. Heilbron, P.F.L., Perez-Guerrero, J.S., Pontedeiro, E.M., RupertiJr., N., Cotta, R.M.: Analytical model for environmental impact
due to solid waste disposal from electricity generation. Hybrid
Meth. Eng. 4(1–2), 1–26 (2002)
5. Landa, E.R.: Naturally occurring radionuclides from industrial
sources: characteristics and fate in the environment. Radioact.
Environ. 10, 211–237 (2007)
6. Pontedeiro, E.M., Heilbron, P.F.L., Cotta, R.M.: Assessment ofmineral industry NORM/TENORM disposal in hazardous landfills. J. Hazard. Mater. 139, 563–568 (2007)
7. Kennedy, J.D., Strenge, D.L.: Residual radioactive contaminationfrom decommissioning: technical basis for translating contamination levels to annual total effective dose equivalent. NUREG
CR-5512, PNL-7994, U.S. Nuclear Regulatory Commission,
Washington, DC (1992)
8. Yu, C., Gnanapragasam, E.,Zielen, Biwer, B,M., Kamboi, S., Cheng,
J.-J., Yuan, Klett, T., LePoire, D., Ziele, A.J., Chen, S.Y., Williams,
W.A., Wallo, A., Domotor, S., Mo, T., Schwartzman, A.: User’s
manual for RESRAD-OFFSITE, Version 2. ANL/EVS/TM/07-1,
Argonne National Laboratory, Argonne, IL (2007)
9. Simunek, J., van Genuchten, M.T., Sejna, M.: Development andapplications of the HYDRUS and STANMOD software packages
and related codes. Vadose Zone J. 7(2), 587–600 (2008)

123

Waste Biomass Valor (2013) 4:197–211
10. Finsterle, S., Doughty, C., Kowalsky, M.B., Moridis, M.B., Pan,
L., Xu, T., Zhang, Y., Pruess, K.: Advanced vadose zone simulations using TOUGH. Vadose Zone J. 7(2), 601–609
11. Javandel, I., Doughty, C., Tsang, C.-F.: Groundwater Transport:
Handbook of Mathematical Models. Water Resour. Monograph
10, Am. Geophys. Union, Washington, D.C. (1984)
12. Leij, F.J., Skaggs, T.H., van Genuchten, M.T.: Analytical solutions for solute transport in three-dimensional semi-infinite porous media. Water Resour. Res. 27(10), 2719–2733 (1991)
13. Bauer, P., Attinger, S., Kinzelbach, W.: Transport of a decay
chain in homogenous porous media: analytical solutions. J. Contam. Hydrol. 49, 217–239 (2001)
14. Vanderborght, J., Kasteel, R., Vereecken, H.: Stochastic continuum transport equations for field-scale solute transport: Overview
of theoretical and experimental results. Vadose Zone J. 5,
184–203 (2006)
15. Pérez Guerrero, J.S., Pimentel, L.C.G., Skaggs, T.H., Van Genuchten, M.T.: Analytical solution of the advection-diffusiontransport equation using a change-of-variable and integral transform technique. Int. J. Heat Mass Transfer 52, 3297–3304 (2009)
16. Cotta, R.M.: Integral Transforms in Computational Heat andFluid Flow. CRC Press, Boca Raton (1993)
17. Cotta, R.M., Mikhailov, M.D.: Hybrid methods and symbolic
computations. In: Minkowycz, W.J., Sparrow, E.M., Murthy, J.Y.
(eds.) Handbook of Numerical Heat Transfer, 2nd edition,
Chapter 16. Wiley, New York (2006)
18. Cotta, R.M., Ungs, M.J., Mikhailov, M.D.: Contaminant transport
in finite fractured porous medium: integral transforms andlumped-differential formulations. Ann. Nucl. Energy 30,
261–285 (2003)
19. Cotta, R.M., Mikhailov, M.D., Ruperti Jr., N.J.: Analysis ofradioactive waste contamination in soils: solution via symbolic
manipulation, Proceedings of the International Conference on theRadiological Accident of Goiânia; 10 Years Later, Brazil,
October 1997 (CD-ROM); also in IAEA Conf. Proceedings
Series, Vienna, IAEA-GOCP, pp. 298–308 (1998)
20. Liu, C., Szecsody, J.E., Zachara, J.M., Ball, W.P.: Use of the generalized integral transform method for solving equations of solutetransport in porous media. Adv. Water Resour. 23, 483–492 (2000)
21. Leal, M.A., Ruperti Jr, N.J.: A hybrid solution for simulation of
2-D contaminant transport in heterogeneous subsurface systems.
Hybrid Methods Eng. 3(2), 111–129 (2001)
22. Barros, F.P.J., Mills, W.B., Cotta, R.M.: Integral transform
solution of two-dimensional model for contaminants dispersion inrivers and channels with spatially variable coefficients. Environ.
Model. Softw. 21, 699–709 (2006)
23. Skaggs, T.H., Jarvis, N.J., Pontedeiro, E.M., van Genuchten,
M.T., Cotta, R.M.: Analytical advection-dispersion model fortransport and plant uptake of contaminants in the root zone.
Vadose Zone J. 6, 890–898 (2007)
24. Pérez Guerrero, J.S., Skaggs, T.H., van Genuchten, M.T.: Analytical solution for multi-species contaminant transport in finite
media with time-varying boundary conditions. Transp. PorousMedia 85, 171–188 (2010)
25. Cotta, R.M., Quaresma, J.N.N., Sphaier, L.A., Naveira-Cotta, C.P.:
Unified integral transform approach in the hybrid solution of multidimensional nonlinear convection-diffusion problems. 14th Int.
Heat Transfer Conf., Washington, DC, USA, August 2010
26. Wolfram, S.: Mathematica 5.2, Cambridge, UK, Wolfram Media(2005)
27. Bateman, H.: The solution of a system of differential equations
occurring in the theory of radioactive transformations. Proc.
Cambridge Philos. Soc. 15, 423–427 (1910)
28. Higashi, K., Pigford, T.H.: Analytical models for migration ofradionuclides in geologic sorbing media. J. Nucl. Sci. Technol.
17(9), 700–709 (1980)

Waste Biomass Valor (2013) 4:197–211
29. van Genuchten, M.T.: Convective-dispersive transport of solutes
involved in sequential first-order decay reactions. Comput.
Geosci. 11(2), 129–147 (1985)
30. Lung, H.C., Chambré, P.L., Pigford, T.H., Lee, W.W.L.: Transport of radioactive decay chains in finite and semi-infinite porous
media, Earth Sciences Division, Lawrence Berkeley Laboratory,
Report LBL23987 (1987)
31. IAEA: Safety Assessment Methodologies for Near Surface Disposal Facilities, Vol. 1, Review and Enhancement of Safety
Assessment Approaches and Tools, International Atomic Energy
Commission (IAEA), Vienna, Austria (2004)
32. Heilbron, P.F.L., Xavier, A.M., Pontedeiro, E.M., Ferreira, R.S.:
Segurança Nuclear e Proteção do Meio Ambiente. Brazil, Editora
E-Papers, Rio de Janeiro (2004)
33. Pontedeiro, E.M., van Genuchten, M.T., Cotta, R.M., Simunek,
J.: The effects of preferential flow and soil texture on risk
assessments of a NORM waste disposal site. J. Hazard. Mater.
174(1–3), 648–655 (2010)
34. EPA: Soil Screening Guidance: Technical Background Document, EPA/540/R-95/128, PB96–963502. Office of Solid Waste
and Emergency Response, U.S. Environmental Protection
Agency, Washington, DC (1996)
35. Almeida, A.R., Cotta, R.M.: A comparison of convergenceacceleration schemes for eigenfunction expansions of partial
differential equations. Int. J. Numer. Method. Heat Fluid Flow
6(no. 6), 85–97 (1996)
36. Cotta, R.M., Mikhailov, M.D.: Heat Conduction: Lumped Analysis, Integral Transforms, Symbolic Computation. Wiley, NewYork (1997)
37. Sphaier, L.A., Cotta, R.M., Naveira-Cotta, C.P., Quaresma,
J.N.N.: The UNIT algorithm for solving one-dimensional convection-diffusion problems via integral transforms. Int. Commun.
Heat Mass Transf. 38(5), 565–571 (2011)

211
38. Simunek, J., van Genuchten, M.T., Sejna, M.: The HYDRUS-1D
Software Package for Simulating the One-Dimensional Movement of Water, Heat And Multiple Solutes in Variably-Saturated
Media, Version 3.0. HYDRUS Software Series 1, Department ofEnvironmental Sciences, University of California, Riverside,
240 p (2005)
39. Richards, L.A.: Capillary conduction of fluid through porous
mediums. Physics 1, 318–333 (1931)
40. Orlande, H.R.B., Cotta, R.M., Su, J., Couto, P., Naveira, C.P.,
Moreira, P.H.S., Sabino, V.: Estudo da Dispersão de Rejeitos
Radioativos das Células dos Depósitos de Rejeitos Lı́quidos
Tratados (Ponds) na Unidade de Concentrado de Urânio na INBCaetité, BA, Relatório final COPPETEC Projeto PEM 7552
(2006)
41. van Genuchten, M.T.: A closed-form equation for predicting thehydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J.
44, 892–898 (1980)
42. Orlande, H.R.B., van Genuchten, M. T., Cotta, R.M., Moreira,P.H.: Bayesian estimation of hydraulic and solute transportparameters from laboratory soil columns experiments. In: Rebay,
M. (ed.), Proc. Int. Symp. on Convective Heat and Mass Transfer
in Sustainable Energy, CONV-09, April 26–May 1, Hammamet,
Tunisia (2009)
43. IAEA: Safety assessment methodologies for near surface disposal
facilities, Vol. 2, Test cases, Int. Atomic Energy Commission(IAEA), Vienna, Austria (2004)
44. ICRP: Radionuclide Transformations Energy and Intensity ofEmissions. International Commission on Radiological Protection,
ICRP Publ. 38. Pergamon Press, Oxford (1983)
45. Sun, Y., Petersen, J.N., Clement, T.P., Skeen, R.S.: Developmentof analytical solutions for multispecies transport with serial andparallel reactions. Water Resour. Res. 35(1), 185–190 (1999)

123