R-Lum / RLumModel
1
#' sequence step irradiation
2
#'
3
#' This function simulates the irradiaton of quartz in the energy-band-model.
4
#'
5
#' @param temp \code{\link{cnumeric}} (\bold{required}): temperature [deg. C] at which the dose should be applied
6
#'
7
#' @param dose \code{\link{numeric}} (\bold{required}): dose to apply in Gray
8
#'
9
#' @param dose_rate \code{\link{numeric}} (\bold{required}): Dose rate in Gy/s
10
#'
11
#' @param n \code{\link{numeric}} or \code{\linkS4class{RLum.Results}} (\bold{required}):
12
#' concentration of electron-/holetraps, valence- and conduction band
13
#' from step before. This is necessary to get the boundary condition for the ODEs.
14
#'
15
#' @param parms \code{\linkS4class{RLum.Results}} (\bold{required}): The specific model parameters are used to simulate
16
#' numerical quartz luminescence results.
17
#'
18
#' @return This function returns an RLum.Results object.
19
#'
20
#' @section Function version: 0.1.2 [2017-11-20]
21
#'
22
#' @author Johannes Friedrich, University of Bayreuth (Germany),
23
#'
24
#' @references
25
#'
26
#' Bailey, R.M., 2001. Towards a general kinetic model for optically and thermally stimulated
27
#' luminescence of quartz. Radiation Measurements 33, 17-45.
28
#'
29
#' @seealso \code{\link{model_LuminescenceSignals}}, \code{\link{simulate_RL}}
30
#'
31
#' @examples
32
#'
33
#' #so far no example available
34
#'
35
#' @noRd
36
.simulate_irradiation <- function(
37
  temp,
38
  dose,
39
  dose_rate,
40
  n,
41
  parms
42
){
43

44
# check input arguments ---------------------------------------------------
45

46
  ##check if temperature is > 0 K (-273 degree celsius)
47 5
  if(temp < -273){
48 5
    stop("\n [.simulate_irradiation()] Argument 'temp' has to be > 0 K!")
49
  }
50
  ##check if dose_rate is a positive number
51 5
  if(dose_rate < 0){
52 5
    stop("\n [.simulate_irradiation()] Argument 'dose_rate' has to be a positive number!")
53
  }
54

55
  ##check if dose is a positive number
56 5
  if(dose < 0){
57 5
    stop("\n [.simulate_irradiation()] Argument 'dose' has to be a positive number!")
58
  }
59

60
  ##check if n is a RLum object
61 5
  if(class(n) != "RLum.Results"){
62 5
    n <- n
63
  } else {
64 5
    n <- n$n
65
  }
66

67
# solving ODE ---------------------------------------------------
68
  
69 5
  if(dose != 0){
70
  
71
  ##============================================================================##
72
  # SETTING PARAMETERS FOR IRRADIATION
73
  #
74
  # R: electron-hole-production-rate (in Bailey 2004: 2.5e10, Bailey 2002: 3e10, else: 5e7)
75
  # P: Photonflux (in Bailey 2004: wavelength [nm])
76
  # b: heating rate [deg. C/s]
77
  ##============================================================================##
78
  ## check if R is given in customized parameter sets
79 5
  if("R" %in% names(parms) && parms$R != 0){
80
    
81 0
    R <- dose_rate*parms$R
82
    
83
  } else {
84

85 5
    if(parms$model == "Bailey2004"){
86 5
      R <- dose_rate*2.5e10
87
    } else {
88
      
89 5
      if(parms$model == "Bailey2002"){
90 5
        R <- dose_rate*3e10
91
      } else {
92
        
93 5
        if(parms$model == "Friedrich2018"){
94 0
          R <- dose_rate*6.3e7
95
        } else {
96
        
97 5
          R <- dose_rate*5e7  # all other simulations
98
        }
99
      }
100
    }
101
  }
102

103 5
  P <- 0
104 5
  b <- 0
105

106
  ##============================================================================##
107
  # SETTING PARAMETERS FOR ODE
108
  ##============================================================================##
109

110 5
  times   <- seq(0, dose/(dose_rate), by = (dose/dose_rate)/100)
111 5
  parameters.step <- .extract_pars(parameters.step = list(
112 5
    R = R,
113 5
    P = P,
114 5
    temp = temp,
115 5
    b = b,
116 5
    times = times,
117 5
    parms = parms))
118
  
119
  ##============================================================================##
120
  # SOLVING ODE (deSolve requiered)
121
  ##============================================================================##
122
  
123 5
  out <- deSolve::lsoda(y = n, times = times, parms = parameters.step, func = .set_ODE_Rcpp, rtol = 1e-6, atol = 1e-6);
124

125
  ##============================================================================##
126
  # TAKING THE LAST LINE OF "OUT" TO COMMIT IT TO THE NEXT STEP
127
  ##============================================================================##
128

129 5
  return(Luminescence::set_RLum(class = "RLum.Results",
130 5
                  data = list(
131 5
                    n = out[length(times),-1],
132 5
                    temp = temp
133
                  )))
134

135
  } else {
136
    
137 5
  return(Luminescence::set_RLum(class = "RLum.Results",
138 5
                                data = list(
139 5
                                  n = n,
140 5
                                  temp = temp
141
                                )))
142
  
143
  
144
}
145
}

Read our documentation on viewing source code .

Loading