/
run_MC_CW_IRSL_TUN.R
153 lines (143 loc) · 5.68 KB
/
run_MC_CW_IRSL_TUN.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
#' @title Run Monte-Carlo Simulation for CW-IRSL (tunnelling transitions)
#'
#' @description Runs a Monte-Carlo (MC) simulation of continuous wave infrared stimulated luminescence
#' (CW-IRSL) using the model for tunnelling transitions. Tunnelling refers to quantum mechanical
#' tunnelling processes from the excited state of the trap,
#' into a recombination centre.
#'
#' @details
#'
#' **The model**
#'
#' \deqn{
#' I_{TUN}(r',t) = -dn/dt = A * exp(-(\rho')^{-1/3} * r')* n (r',t)
#' }
#'
#'Where in the function: \cr
#' A := effective optical excitation rate for the tunnelling process (s^-1) \cr
#' r' := the dimensionless tunnelling radius \cr
#' \eqn{\rho}' := `rho'` the dimensionless density of recombination centres (see Huntley (2006)) \cr
#' t := time (s) \cr
#' n := the instantaneous number of electrons corresponding to the radius r' at time t
#'
#' @param A [numeric] (**required**): The effective optical excitation rate for the tunnelling process
#' (`s^-1`).
#'
#' @param rho [numeric] (**required**): The density of recombination centres
#' (defined as \eqn{\rho}' in Huntley 2006) (dimensionless).
#'
#' @param times [numeric] (**required**): The sequence of time steps within the simulation (s).
#'
#' @param clusters [numeric] (*with default*): The number of created clusters for the MC runs. The input can be the output of [create_ClusterSystem]. In that case `n_filled` indicate absolute numbers of a system.
#'
#' @param N_e [numeric] (*width default*): The total number of electron traps available (dimensionless).
#' Can be a vector of `length(clusters)`, shorter values are recycled.
#'
#' @param r_c [numeric] (*with default*): Critical distance (>0) that must be provided if the
#' sample has been thermally and/or optically pretreated. This parameter expresses the fact
#' that electron-hole pairs within a critical radius `r_c` have already recombined.
#'
#' @param delta.r [numeric] (*with default*): Increments of the dimensionless distance parameter r'
#'
#' @param method [character] (*with default*): Sequential `'seq'` or parallel `'par'`processing. In
#' the parallel mode the function tries to run the simulation on multiple CPU cores (if available) with
#' a positive effect on the computation time.
#'
#' @param output [character] (*with default*): output is either the `'signal'` (the default) or
#' `'remaining_e'` (the remaining charges/electrons in the trap)
#'
#' @param \dots further arguments, such as `cores` to control the number of used CPU cores or `verbose` to silence the terminal
#'
#' @return This function returns an object of class `RLumCarlo_Model_Output` which
#' is a [list] consisting of an [array] with dimension length(times) x length(r) x clusters
#' and a [numeric] time vector.
#'
#' @section Function version: 0.2.0
#'
#' @author Johannes Friedrich, University of Bayreuth (Germany),
#' Sebastian Kreutzer, Institute of Geography, Heidelberg University (Germany)
#'
#' @references
#' Huntley, D.J., 2006. An explanation of the power-law decay of luminescence.
#' Journal of Physics: Condensed Matter, 18(4), 1359.
#'
#' Pagonis, V., Friedrich, J., Discher, M., Müller-Kirschbaum, A., Schlosser, V., Kreutzer, S.,
#' Chen, R. and Schmidt, C., 2019. Excited state luminescence signals from a random distribution of defects:
#' A new Monte Carlo simulation approach for feldspar.
#' Journal of Luminescence 207, 266–272. \doi{10.1016/j.jlumin.2018.11.024}
#'
#' **Further reading**
#'
#' Aitken, M.J., 1985. Thermoluminescence dating. Academic Press.
#'
#' Jain, M., Guralnik, B., Andersen, M.T., 2012. Stimulated luminescence emission from
#' localized recombination in randomly distributed defects.
#' Journal of Physics: Condensed Matter 24, 385402.
#'
#' Chen, R., McKeever, S.W.S., 1997. Theory of Thermoluminescence and Related Phenomena.
#' WORLD SCIENTIFIC. \doi{10.1142/2781}
#'
#' @examples
#' run_MC_CW_IRSL_TUN(
#' A = 0.8,
#' rho = 1e-4,
#' times = 0:50,
#' r_c = 0.05,
#' delta.r = 0.1,
#' method = "seq",
#' clusters = 10,
#' output = "signal") %>%
#' plot_RLumCarlo(norm = TRUE, legend = TRUE)
#'
#' @keywords models data
#' @encoding UTF-8
#' @md
#' @export
run_MC_CW_IRSL_TUN <- function(
A,
rho,
times,
clusters = 10,
r_c = 0,
delta.r = 0.1,
N_e = 200,
method = "seq",
output = "signal",
...){
# Integrity checks ----------------------------------------------------------------------------
if(!output %in% c("signal", "remaining_e"))
stop("[run_MC_CW_IRSL_TUN()] Allowed keywords for 'output' are either 'signal' or 'remaining_e'!",
call. = FALSE)
# Register multi-core back end ----------------------------------------------------------------
cl <- .registerClusters(method, ...)
on.exit(parallel::stopCluster(cl))
# Setting parameters --------------------------------------------------------------------------
r <- seq(abs(r_c[1]), 2, abs(delta.r[1]))
# Enable dosimetric cluster system -----------------------------------------
if(class(clusters)[1] == "RLumCarlo_ClusterSystem"){
N_e <- .distribute_electrons(
clusters = clusters,
N_system = N_e[1])[["e_in_cluster"]]
clusters <- clusters$cl_groups
}
# Expand parameters -------------------------------------------------------
N_e <- rep(N_e, length.out = max(clusters))
# Run model -----------------------------------------------------------------------------------
temp <- foreach(
c = 1:max(clusters),
.packages = 'RLumCarlo',
.combine = 'comb_array',
.multicombine = TRUE
) %dopar% {
results <- MC_C_CW_IRSL_TUN(
times = times,
N_e = N_e[c],
r = r,
rho = rho[1],
A = A[1]
)
return(results[[output]])
} # end c-loop
# Return --------------------------------------------------------------------------------------
.return_ModelOutput(signal = temp, time = times)
}