tgcd (version 2.7)

simqOTOR: Thermoluminescence glow peak simulation

Description

Simulating glow peaks according to the one trap-one recombination center (OTOR) model using the quasi-equilibrium approximation.

Usage

simqOTOR(temps, n0, Nn, Ah, An, ff, ae, 
         hr, outfile = NULL, plot = TRUE)

Value

Return an invisible list containing the following elements:

temps

a vector of temperature values

tl

values of the thermoluminescence intensity

n

variation of concentration of trapped electrons with temperature

sp

parameters used for describing the shape of a glow peak, see function simPeak for details

Arguments

temps

vector(required): temperature values (K) where the values of the thermoluminescence intensity will be computed, it needs to be sorted increasingly

n0

numeric(required): initial concentration of trapped electrons (1/cm^3)

Nn

numeric(required): total concentration of the traps in the crystal (1/cm^3)

Ah

numeric(optional): probability coefficient of electron recombining with holes in the recombination center (cm^3/s)

An

numeric(optional): probability coefficient of electron retrapping in the traps (cm^3/s)

ff

numeric(required): the frequency factor (1/s)

ae

numeric(required): the activation energy (eV)

hr

numeric(with default): the linear heating rate (K/s)

outfile

character(optional): if specified, simulated intensities of glow peaks will be written to a file named "outfile" in CSV format and saved to the current work directory

plot

logical(with default): draw a plot according to the simulated glow peak or not

Details

Function simqOTOR simulates a synthetic glow peak according to the OTOR model using the quasi-equilibrium approximation. This function may be used to simulating glow peaks of first-, second-, and general-order, depending on the given kinetic parameters. The approximate equation of the OTOR model derived using the quasi-equilibrium approximation can be described by (Pagonis et al., 2006):

\(\frac{d_n}{d_T}=\frac{-A_hn^2Sexp(-\frac{E}{kT})}{[nA_h+(N_n-n)A_n]\beta}\)

where \(n\) is the concentration of trapped electrons, \(\frac{d_n}{d_T}\) the rate of change of the concentration of trapped electrons, \(S\) the frequency factor, \(E\) the activation energy, \(T\) the absolute temperature, \(k\) the Boltzmann constant, \(N_n\) the total concentration of the traps in the crystal, \(A_h\) the probability coefficient of electron recombining with holes in the recombination center, \(A_n\) the probability coefficient of electron retrapping in the traps, and \(\beta\) the linear heating rate.

The ordinary equation is solved by the Fortran 77 subroutine lsoda (original version written by Linda R. Petzold and Alan C. Hindmarsh available at Netlib: https://www.netlib.org/odepack/, modified version by R. Woodrow Setzer from the R package deSolve (Soetaert et al., 2010) available at CRAN: https://CRAN.R-project.org/package=deSolve).

References

Pagonis V, Kitis G, Furetta C, 2006. Numerical and practical exercises in thermoluminescence. Springer Science & Business Media.

Soetaert K, Petzoldt T, Setzer RW, 2010. Solving Differential Equations in R: Package deSolve. Journal of Statistical Software, 33(9): 1-25.

See Also

tgcd; simPeak

Examples

    # Synthesizing a glow curve consisting of five glow peaks.
     temps <- seq(330, 730, by=0.5)
     peak1 <- simqOTOR(temps, n0=0.7e10, Nn=1e10, Ah=1e-3, An=1e-7,
       ff=1e14, ae=1.5, hr=1, outfile = NULL, plot = TRUE)
     peak2 <- simqOTOR(temps, n0=0.5e10, Nn=1e10, Ah=1e-7, An=1e-7,
       ff=1e17, ae=1.9, hr=1, outfile = NULL, plot = TRUE)
     peak3 <- simqOTOR(temps, n0=0.2e10, Nn=1e10, Ah=1e-5, An=1e-7,
       ff=1e15, ae=1.45, hr=1, outfile = NULL, plot = TRUE)
     peak4 <- simqOTOR(temps, n0=0.2e10, Nn=1e10, Ah=1e-5, An=1e-7,
       ff=1e9, ae=0.85, hr=1, outfile = NULL, plot = TRUE)
     peak5 <- simqOTOR(temps, n0=0.3e10, Nn=1e10, Ah=1e-7, An=1e-7,
       ff=1e11, ae=1.4, hr=1, outfile = NULL, plot = TRUE)
     peaks <- cbind(peak1$tl, peak2$tl, peak3$tl, peak4$tl, peak5$tl, 
       peak1$tl+peak2$tl+peak3$tl+peak4$tl+peak5$tl)
     matplot(temps, y=peaks, type="l", lwd=2, lty="solid", 
       xlab="Temperature (K)", ylab="TL intensity (counts)")