Determining radium-226 concentration from radon-222 emanation in building materials : a theoretical model

It was developed an improved theoretical model capable to estimate the radium concentration in building materials solely measuring the radon-222 concentration in a confined atmosphere. This non-destructive technique is not limited by the size of the samples, and it intrinsically includes back diffusion. The resulting equation provides the exact solution for the concentration of radon-222 as a function of time and distance in one dimension. The effective concentration of radium-226 is a fit parameter of this equation. In order to reduce its complexity, this equation was simplified considering two cases: low diffusion in the building material compared to the air, and a building material initially saturated with radon-222. These simplified versions of the exact one dimension solution were used to fit experimental data. Radon-222 concentration was continuously measured for twelve days with an AlphaGUARDTM detector, located at the Laboratory of Applied Nuclear Physics at Universidade Tecnológica Federal do Paraná (UTFPR). This model was applied to two different materials: cement mortar and concrete, which results were respectively (15.7±8.3)Bq/kg and (10.5±2.4)Bq/kg for the radium-226 effective concentration. This estimation was confronted with the direct measurements of radium in the same materials (same sources) using gamma-ray spectrometry, fulfilled at Centro de Desenvolvimento da Tecnologia Nuclear (CDTN), which results were respectively (13.81±0.23)Bq/kg and (12.61±0.22)Bq/kg.


INTRODUCTION
Building materials such as concrete, cement, brick, and gypsum usually present small amounts of radioactive elements.These naturally occurring radioactive materials (NORM) contribute to individual exposure to ionizing radiation.Exposure can occur in two ways: external and internal.External Barreto, et al. • Braz.J. Rad.Sci.• 2019 2 exposure is due to sources in the surroundings, outside the individual's body.This is the case of radium-226 present in a building material, which decays by alpha emission into radon-222, subsequently emitting a gamma ray.Radon is a noble gas, thus it is chemically inert.It means radon-222 is colorless, odorless, and tasteless, so it can be unknowingly inhaled.The inhalation of radon-222 is an example of internal exposure.Radon-222 is an alpha emitter with half-life of 3.8 days, which decay products are post-transition alpha and gamma radioactive metals.Therefore, these products can be absorbed by lung tissue, readily producing DNA damage, which can lead to pulmonary cancer [1].The World Health Organization estimates that between 3% to 14% of lung cancer deaths are associated with radon [2].The transport mechanisms of radon-222 in solid materials are emanation, convection (diffusion and advection), and exhalation [3,4].Emanation occurs when the mineral grain releases the gas molecule.It means emanation is the rate of flow of radon in the material, and it is proportional to the activity of radium-226.The emanation coefficient is the fraction of radon that reaches the pore space, and it is dependent on the crystalline structure, grain size and moisture [4,5].Its reported values vary from 1% to 50% over a wide range of materials, conditions, definitions and measurement methods [6,7].Convection is the movement of groups of molecules driven by concentration gradient (diffusion) and pressure gradient (advection).The most common approaches to describe convection in fluid mixtures through porous media are the advective-diffusive model (ADM) and the dusty-gas model (DGM) [8].The ADM is the straightforward union of ordinary diffusion using Fick's Law and advection using Darcy's Law.In the ADM, slip effects and Knudsen diffusion are included through a Klinkenberg correction [9].The DGM is a sophisticate model based on kinetic theory and the Stefan-Maxwell equations [10].Although both ADM and DGM have the contribution of the pressure and concentration gradients, in DGM the multiplying factors of the gradients are depend on the mole fractions.However, in the trace gas limit these factors become constant, and ADM and DGM take the same form, differing only by the meaning of the coefficients of the gradients [8].Since radon is usually a trace gas, its convection in porous media in the absence of weather is dominated by diffusion [3,4].The escape of radon from the material into the open air is called exhalation.When radon escapes into a confined atmosphere, such as a room, the term exhalation is sometimes substituted by entry [3].The exhalation rate is the flux density in the surface of the material.Once the gas can accumulate inside dwellings, it is necessary to determinate the indoor radon-222 activity [11][12][13][14][15], its exhalation from soil [16][17][18][19][20], building materials [21][22][23][24][25][26][27], and objects rich in radium-226 from private and museum collections [28][29][30].There are several different techniques to determine the emanation factor and the effective radium-226 activity from radon-222 exhalation [18,25,[31][32][33].Among them, there are two main approaches that use a time varying mean radon-222 concentration approximation in a confined atmosphere [31,32].In 2002, Ferry et al. proposed a non-destructive technique to determine the emanation factor in small rocks using accumulation chambers [31].In this method, the samples must be smaller than the diffusion length.Using this condition, it is possible to derive an equation to describe the concentration of radon with time.In this method, the decay rate is effective and a fitting parameter.In 2004, Jang et al. proposed a mathematical model to estimate the radon-222 mean concentration in a confined atmosphere [32].This model uses the equilibrium concentration along the external atmosphere to derive the exhalation rate as a function of concentration.By using the exhalation rate, it is possible to derive an equation for the mean concentration as a function of time.The expression presented by Jang is equivalent to the one proposed by Ferry.But it places an improved interpretation of the parameters, since it considers the back diffusion, which is proportional to the concentration in the chamber.In 2006, Shafi-ur-Rehman et al. explored Jang's model to study the exhalation rate of soil and sand [33].In the same year, Rabi et al. explored Jang's model in 2D numerical simulations studies [21].In 2008, Faheem et al. used Jang's model to study the influence of moisture on the exhalation of soil and build materials samples [34].In 2013, Perrier and Girault used Ferry's model to design a method to determine the effective radium concentration in samples using small scintilation flasks [35].In 2016, Perrier et al. revised this methodology for high-sensitivity measurements of radium-226 in liquid samples through radon-222 emanation [36].In 2017, Girault et al. used Ferry's technique to measure the radon-222 loss from meteorites [30].In the same year, Campos et al. used Jang's equation to study the radon exhalation in recycled phosphogypsum [37].
In this work, it is proposed an improved equation to describe the mean radon-222 concentration in a confined atmosphere as a function of time.By using this mathematical model, it is possible to estimate the radium-226 effective concentration in materials using experimental data concerning the radon-222 activity.

Mathematical Development
Let n( r,t) be the radon-222 density (concentration, in atoms per unit volume) at some point of the space at a given time.The time rate of change of the concentration has the contribution of three terms: the radon-222 emanation per unit volume (w), the negative of current density divergence ( ∇ • j), and the negative of radon-222 activity density (λ n, where λ is the decay constant).
The radon-222 emanation per unit volume w is the production rate density, which is equal to the product of the radium-226 effective activity density α (in decays per unit time per unit mass) and the material density ρ (in unit mass per unit volume).In a time-scale of days, this activity is approximately constant.The emanation is assumed to be uniform inside a homogeneous material, and zero at any point of the atmosphere.For convection in a porous media, in first approximation, it is possible to write the current density adding the Fick's first law of diffusion and the Darcy's law, j = −D ∇n − n k p µ ∇p.The first term containing the concentration gradient ∇n drives the diffusion of radon-222 in relation to the gas mixture.In this term, D is an effective diffusion coefficient, which is affected by porosity, tortuosity and moisture [3,4].The second term contains the pressure gradient ∇p, the intrinsic permeability of the medium k p and the gas viscosity µ.This term is responsible for driving the gas mixture as a whole (not only radon-222).Since the same equation must be applied to the background atmosphere of concentration n air , it is possible to write the pressure gradient in terms of j air and ∇n air .However, on near stationary equilibrium, in the absence of temperature and pressure fluctuations, the background atmosphere is assumed to be static, j air = 0.So it is possible to equate the radon-222 current density and the concentration gradients, j = −D ∇n + n n air D air ∇n air .Therefore, in the limit of trace gas (n/n air → 0), the second term in the right side of this equation can be neglected, and it is assumed the convection is dominated by diffusion, which is a common practice [3,4].By replacing the Fick's first law into equation 1, it is obtained a second order partial differential equation.
This equation can be simplified if it is used an auxiliary density function q, which is also a function of time and position, n = w/λ + qe −λt .Finally, by replacing this equation into equation 2, we obtain the Fick's second law for q, ∂ q ∂t = D∇ 2 q.

General One-Dimensional Solution
In one dimension, the independent variables are only x and t.By using the method of separation of variables, it is derived the general solution to this equation, which is a term that varies linearly with distance, a constant, and a combination of exponentials.
In equation 3, a κ and κ depends on the initial and boundary conditions.The use of prime in a κ is for convenience (it will be clear along the development).It is important to realize that κ 2 must be real, otherwise q(x,t) would oscillate in time.An oscillating density has no physical meaning in this problem, since Fick's second law is only first order with respect to time.Therefore, κ is either a real or a pure imaginary number.This means the solution can include trigonometric and hyperbolic functions of x.To apply the general solution (eq.3), it is considered a cylindrical solid of length L I .
This solid is filling the left side of a tube of length L and area S (fig.1).In such a way that the right side is filled with an atmosphere of length L E = L − L I .The ends of the tube are closed.The values of D are different in the solid material and the external atmosphere.It means that a κ and κ possibly vary as well.From this point forward, it is used the indices I (internal) and E (external) to denote the material and the external atmosphere, respectively.One must note that by using this symmetry and boundaries, the gas necessarily interacts with the lateral wall of the cylinder (perpendicular to x direction), and the concentration is dependent on both x and r (two dimensions).
The solution in r are Bessel functions of zero order, but this will not be explored in this development to keep it simple.

Applying the Boundary Conditions at the Extremities
At the ends of the tube, the current density must be zero, j I (0,t) = j E (L,t) = 0.
These relations are satisfied for any value of t only if b = 0 and the terms with ±κ cancel in pairs.This happens if a κ I = a −κ I and a κ E e κ E L = a −κ E e −κ E L .If κ is a real number, then these conditions result in hyperbolic functions of x for q(x,t).Otherwise, if κ is purely imaginary, it will result in trigonometric functions (κ = ik, where k is real).For convenience, we will rewrite equation 3 in terms of trigonometric functions, assuming κ is purely imaginary. (5) In these equations, a k = 2a k .It is used k ≥ 0 in the sum to denote that the negative values have already been taken into account.Despite this, one must note that k becomes purely imaginary when κ is real.In this case, hyperbolic functions will be used.The appropriate procedure was to expand the summation, discriminating both cases.But this would duplicate the size of the equations, turning them cumbersome to read.Therefore, we are neglecting the mathematical rigor for the sake of clarity.

Applying the Boundary Conditions at the Interface
In three dimensions, the boundary conditions at the interface material-atmosphere is determined by applying both Gauss and Stokes Theorems with the current density.These result in two conditions that must be satisfied at the interface.

Current density at the interface
The current density must be the same at both sides on the interface, j I (L I ,t) = j E (L I ,t).This leads to the following relation.
By using k E = β k I , one can rewrite a k E in terms of a k I .
By replacing equation 9 into equation 6, one can rewrite q E in terms of k I , a k I and β .Therefore, for convenience, we will suppress the index I in k I from this point forward.

Radon-222 concentration at the interface
The product of the diffusion coefficient and the concentration must be the same at both sides on the interface, D I n I (L I ,t) = D E n E (L I ,t).This leads to the following relation.
Over again, this equation is satisfied for any value of t only if the terms with the same exponentials are equal.The time independent terms equate each other when k is purely imaginary, k 2 = −λ /D I .
For convenience, we will represent it as κ λ = λ /D I , which is recognized as the inverse of the diffusion length [4].This is the unique instance in which hyperbolic functions appear.
By using this equation, it is possible to express a κ λ in terms of the other constants.
All other terms in the summation equate independently of a k , which depend solely on the initial conditions.These terms provide a transcendental equation, which define a discrete set of real values for k.
The terms in which k = 0 provide one equation for c I and c E , while a second equation will be provided by the initial conditions.

Applying the Initial Conditions to the Number of Radon-222 Atoms
By integrating the time rate of change of the concentration (eq. 1) in the system volume, it is obtained the time rate of change of the number of radon-222 atoms, Ṅ = W − λ N (in Newton's notation).In this, W represents the net radon-222 production rate, which is equal to the radium-226 activity, while λ N is the radon-222 activity.The second term in equation 1 (current density divergence) vanishes because the current density is zero at the external boundaries.The solution to the number of radon-222 atoms in the system is straightforward, where N 0 represents the initial number of atoms in the system.
The direct integration of the concentration in the system volume must provide the same result.In the one-dimensional case, the integration is limited to the length of the tube.
In this equation, the first term in the right side is N I (t), while the second term is N E (t).Both functions can be obtained with the aid of equations 5 and 6. ) These equations have the same summation term, differing only by the sign, in such a way they vanish when N I (t) and N E (t) are added.By doing so, and subtracting from equation 15, it is derived the second equation relating c I and c E .
By using equations 14 and 19, one finds c I and c E .
In these equations, the volume V is equal the product of the cross section area S and the length L.
So, V I = SL I is the volume of the solid, and V E = SL E is the volume of the external atmosphere.

Radon-222 Density and Number of Radon-222 Atoms
To obtain the radon-222 density formula, one start by separating the summation in q(x,t) (eq. 5 and 6) in three terms: the first with k purely imaginary (k 2 = −κ 2 λ ), the second with k = 0, and the third converted into a new summation with k real and positive (k > 0).By substituting c (eq. 20 and 21) into q(x,t), one finds that the a 0 terms vanishes.The final step is to replace the κ λ term with a κ λ from equation 12.
In these, the set of k is determined by the transcendental equation 13, while the a k values are determined by the initial density distribution along the volume.The n I (x, ∞) and n E (x, ∞) terms are the equilibrium values of the concentration in each volume. ) By using equations 17 and 18, one can also derive the expressions for the number of radon-222 atoms in the solid and in the external atmosphere.
The set of A k coefficients is dependent on the boundary conditions and the set of a k , which is dependent on the initial conditions.
After an infinite time period, the number of atoms in each volume will converge to N I (∞) and N E (∞). )

Mean Radon-222 Concentration in the Atmosphere
Usually it is measured the radon-222 mean activity density in the external atmosphere, λ nE .By dividing equation 27 by V E , one obtains an expression for the mean radon-222 concentration, nE In this equation, γ = V E /V I , āk = A k /V E , and the equilibrium mean radon-222 concentration is given by nE The first term in equation 31 is determined by the asymptotic behavior of the experimental data.This first term depends on β , D I , and w.The second term in equation 31 depends on N 0 .This second term would not contribute if the solid material were initially isolated, because in this case N 0 /V I = w/λ .

Discussion
The terms in the k summation are very sensitive to the initial conditions, and are crucial to determine the shape of the function.The exponential arguments determine D I and β , while the set of āk determine the initial mean concentration, nE (0).The first term depends on w and β .By finding w and β , using the second term, it is possible to determine N 0 .Therefore, it is necessary at least the first four terms from equation 31 to estimate all relevant unknown system parameters and initial conditions.

Low diffusion in the solid material
If the diffusion coefficient in the solid is much lower than the coefficient in the external atmosphere, D I << D E , then β is significantly small.In this case, it is worth to expand the first and the second terms of equation 31 in a Maclaurin series.By neglecting the higher order terms, it is derived a simplified version of equation 31.

Solid material saturated with radon-222
The second term in equation 31 is a source of uncertainty because it depends on three unknowns (w, β and N 0 ).This term is specially problematic if the solid is near the saturation, N 0 /V I ≈ w/λ .But if the experiment is designed to start with a material saturated (N I (0) = N I (∞)) and an atmosphere with no radon-222 (N E (0) = 0), then the initial number of atoms in the system will be equal to the equilibrium value in the solid, N 0 = N I (∞).By adding equations 30 and 29, and dividing by V I , it is obtained an expression for the equilibrium mean radon-222 concentration in the solid, nI (∞) = N I (∞)/V I .By substituting this expression into equation 31, it is obtained a less sensitive second term.This term now depends only on β , since nE (∞) is determined by the asymptotic behavior.

Low diffusion coefficient and saturated solid material
By joining the previous two approximations, it is possible to write a very simplified equation for the mean concentration in the external atmosphere.
By using this equation, at principle, it would be necessary only the first two terms to determine w.Applying the one-dimensional model to this setup demands caution, since the area changes along the circuit and the air is pumped through the confined atmosphere.Therefore, both the solid and external atmosphere lengths (L I and L E ) are "effective".They could be used as adjustable parameters to tune the method, turning it a semi-empirical approach.But this analysis demands extra effort and it will be done in further developments.In this analysis, it was used as "effective lengths" simply dividing the volumes by the exposed area.So the effective length of the solid material was L I = (1.11± 0.05) cm, while the effective length of the external atmosphere was L E = (4.7 ± 0.1) cm.

Fitting Procedure
One should be capable of fitting the experimental results using the expression for the mean concentration nE (eq.31).After a short stage of testing, it was realized that the second term of this expression is cumbersome.This term is important in defining the shape of the curve in the "middle range", in the range of 50 hours to 100 hours, where the initial distribution connects to the asymptotic behavior.But it depends on an unknown initial condition (N 0 ) and two other fit parameters (w and β ).It implies that without knowing N 0 , the second term is useless to determine the unknown parameters of the system.
The low diffusion approximation for nE (eq.33) presents the same issue of equation 31.However, the solid material saturated with radon-222 approximation for nE (eq.34) does not depend on N 0 .
Therefore, after the initial tests, it was decided to fit the experimental data using both "saturated solid" approximations, named as Model A (eq. 36) and Model B (eq. 37), derived from equations 34 and 35 respectively.
In these equations, λ 1 = D I k 2 1 + λ and λ 2 = D I k 2 2 + λ , where k 1 and k 2 are the first and second roots of the transcendental equation 13.The program used to the fit the data was GNUPLOT, which can be programmed to interactively solve the transcendental equation for k during the fitting procedure.fits are in good agreement with the direct measurement of radium-226 using gamma-ray spectrometry.Model A resulted in (17.5 ± 7.4) Bq/kg and (13.0 ± 3.2) Bq/kg, for cement mortar and concrete respectively.Model B provided (15.7 ± 8.3) Bq/kg for cement mortar and (10.5 ± 2.4) Bq/kg for concrete.These results overlap the direct measurement, which was (13.81 ± 0.23) Bq/kg for cement mortar, and (12.61 ± 0.22) Bq/kg for concrete.The value of β relates the ratio of the effective diffusion coefficients (D I and D E ), and it determines the secular equilibrium of radon-222 in the atmosphere.Despite being an important fit parameter, it is difficult to compare it to experimental data.

ACKNOWLEDGMENT
The authors are thankful to CNEN, CAPES, CNPq and Fundação Araucaria of Paraná State for their support.

( 7 )
This equation is satisfied for any value of t only if the exponentials are equal.It means D I k 2 I = D E k 2 E , and consequently k E = β k I where β = D I D E .It implies the factors that multiply the exponentials equate on both sides.
Lets assume L I , L E , S and λ are known parameters, and one wishes to use equation 31 to fit experimental data to determine w.The w is the most important parameter because it provides the radon-222 emanation, which is used to determine the radium-226 activity density α in the material.The unknown parameters are β , D I , and w, which determine the set of k.The unknown initial condition is the radon-222 concentration along the volume at t = 0, which is determined by N 0 and the set of a k .

Radon- 222
activity concentration was continuously measured for twelve days with an AlphaGUARD TM detector.The experimental setup is showed in reference[38, p. 50, fig.16].The detector measures activity concentration [Bq/m 3 ] which is equal to the product of the radon-222 decay constant and the mean concentration in a confined atmosphere.The half-life of radon-222 used was 3.82 days.When converted to decay constant, it results in a value of 2.10 × 10 −6 s −1 .An air pump circulated the air inside the atmosphere along a circuit of tubes connected to the detector.The measurements were made on two different materials: cement mortar and concrete.The testing samples were four cylinders with (10.0 ± 0.1) cm height and (5.0 ± 0.1) cm diameter, resulting in a total material volume of V I = (785 ± 32) cm3 .The total mass of the testing samples was (1.704 ± 0.001) kg for cement mortar and (1.756 ± 0.001) kg for concrete.The cylinders were placed in a jar connected to the circuit.The external atmosphere volume was V E = (2515 ± 50) cm 3 (the volume of the chamber subtracted by the volume of the samples).The exposed area of each set of four cylinders was S = (707 ± 13) cm 2 .

Figure 2
Figure2presents the experimental data and the fit using equations 36 and 37 for cement mortar and concrete, respectively.

Figure 2 :
Figure 2: Experimental results using cement mortar and concrete samples.The fit results for models A and B are presented.

Table 1
presents the results of the fit.Despite the large experimental uncertainty, the results of the Barreto, et al. • Braz.J. Rad.Sci.• 2019 13

Table 1 :
4esults of the fit using the saturated solid models Bq/kg) 17.5 ± 7.413.0± 3.1 15.7 ± 8.3 10.5 ± 2.4 β 0.17 ± 0.04 0.20 ± 0.02 0.15 ± 0.04 0.19 ± 0.02 5. CONCLUSIONS In this work, it was developed a detailed mathematical model to describe the mean radon-222 concentration in a chamber as a function of time.Using this model, it is possible to estimate the effective radium-226 concentration in samples, independent of its size.This non-destructive technique intrinsically includes back diffusion.It uses an AlphaGUARD TM detector, which is portable and low cost, or any other similar device.The mathematical model presents improvements over the previous approaches because it is exact in one dimension in the absence of convection.By controlling the initial conditions, the simplified versions of the model can provide good results.The best estimation for the radium-226 effective concentration using this model for cement mortar and concrete samples were respectively (15.7 ± 8.3) Bq/kg and (10.5 ± 2.4) Bq/kg.These results are in good agreement to direct measurements using gamma-ray spectrometry, which results were respectively (13.81 ± 0.23) Bq/kg and (12.61 ± 0.22) Bq/kg.