SCM

SCM Repository

[rmetrics] Annotation of /pkg/randtoolbox/R/pseudoRNG.R
ViewVC logotype

Annotation of /pkg/randtoolbox/R/pseudoRNG.R

Parent Directory Parent Directory | Revision Log Revision Log


Revision 4229 - (view) (download)

1 : dutangc 4229 ##
2 :     # @file pseudoRNG.R
3 :     # @brief R file for all pseudo RNGs
4 :     #
5 :     # @author Christophe Dutang
6 :     # @author Petr Savicky
7 :     #
8 :     #
9 :     # Copyright (C) 2009, Christophe Dutang,
10 :     # Petr Savicky, Academy of Sciences of the Czech Republic.
11 :     # All rights reserved.
12 :     #
13 :     # The new BSD License is applied to this software.
14 :     # Copyright (c) 2009 Christophe Dutang, Petr Savicky.
15 :     # All rights reserved.
16 :     #
17 :     # Redistribution and use in source and binary forms, with or without
18 :     # modification, are permitted provided that the following conditions are
19 :     # met:
20 :     #
21 :     # - Redistributions of source code must retain the above copyright
22 :     # notice, this list of conditions and the following disclaimer.
23 :     # - Redistributions in binary form must reproduce the above
24 :     # copyright notice, this list of conditions and the following
25 :     # disclaimer in the documentation and/or other materials provided
26 :     # with the distribution.
27 :     # - Neither the name of the Academy of Sciences of the Czech Republic
28 :     # nor the names of its contributors may be used to endorse or promote
29 :     # products derived from this software without specific prior written
30 :     # permission.
31 :     #
32 :     # THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
33 :     # "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
34 :     # LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
35 :     # A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
36 :     # OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
37 :     # SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
38 :     # LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
39 :     # DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
40 :     # THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
41 :     # (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
42 :     # OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
43 :     #
44 :     #
45 : dutangc 3644 ### pseudo random generation
46 :     ###
47 :     ### R functions
48 :     ###
49 :    
50 :    
51 :     ### set the seed ###
52 :    
53 :     setSeed <- function(seed)
54 :     invisible( .Call("doSetSeed", seed) )
55 :    
56 :    
57 :     ### pseudo random generation ###
58 :    
59 :     congruRand <- function(n, dim = 1, mod = 2^31-1, mult = 16807, incr = 0, echo)
60 :     {
61 :     if(!is.numeric(n) || any(n <=0))
62 :     stop("invalid argument 'n'")
63 :     if(!is.numeric(dim) || length(dim) !=1 || any(dim <= 0))
64 :     stop("invalid argument 'dim'")
65 :     if(!is.numeric(mod) || length(mod) !=1)
66 :     stop("invalid argument 'mod'")
67 :     if(!is.numeric(mult) || length(mult) != 1 || mult > mod || mult < 0)
68 :     stop("invalid argument 'mult'")
69 :     if(!is.numeric(incr) || length(incr) != 1 || incr > mod || incr < 0)
70 :     stop("invalid argument 'incr'")
71 :    
72 :     if(missing(echo))
73 :     echo <- FALSE
74 :    
75 :     if(length(n) > 1)
76 :     res <- .Call("doCongruRand", length(n), dim, mod, mult, incr, echo)
77 :     else
78 :     res <- .Call("doCongruRand", n, dim, mod, mult, incr, echo)
79 :     if(dim == 1)
80 :     as.vector(res)
81 :     else
82 :     as.matrix(res)
83 :     }
84 :    
85 :     SFMT <- function(n, dim = 1, mexp = 19937, usepset = TRUE, withtorus = FALSE, usetime = FALSE)
86 :     {
87 :     if(n <0 || is.array(n))
88 :     stop("invalid argument 'n'")
89 :     if(dim < 0 || length(dim) >1)
90 :     stop("invalid argument 'dim'")
91 :     if(!is.logical(withtorus) && !is.numeric(withtorus))
92 :     stop("invalid argument 'withtorus'")
93 :     if(!is.numeric(mexp))
94 :     stop("invalid argument 'mexp'")
95 :     if(!is.logical(usepset))
96 :     stop("invalid argument 'usepset'")
97 :    
98 :     authorizedParam <- c(607, 1279, 2281, 4253, 11213, 19937, 44497, 86243, 132049, 216091)
99 :    
100 :     if( !(mexp %in% authorizedParam) )
101 :     stop("'mexp' must be in {607, 1279, 2281, 4253, 11213, 19937, 44497, 86243, 132049, 216091}. ")
102 :    
103 :    
104 :     if(!is.logical(withtorus))
105 :     {
106 :     if(0 < withtorus && withtorus <= 1)
107 :     nbTorus <- floor( withtorus * n )
108 :     if(withtorus <=0 || withtorus > 1)
109 :     stop("invalid argument 'withtorus'")
110 :     }
111 :     if(is.logical(withtorus))
112 :     {
113 :     if(!withtorus)
114 :     nbTorus <- 0
115 :     else
116 :     stop("invalid argument 'withtorus'")
117 :     }
118 :    
119 :     if(nbTorus == 0)
120 :     {
121 :     if(length(n) > 1)
122 :     res <- .Call("doSFMersenneTwister", length(n), dim, mexp, usepset)
123 :     else
124 :     res <- .Call("doSFMersenneTwister", n, dim, mexp, usepset)
125 :     }
126 :     else
127 :     {
128 :     restorus <- torus(nbTorus, dim, mixed = FALSE, usetime = usetime)
129 :    
130 :     if(length(n) > 1)
131 :     res <- .Call("doSFMersenneTwister", length(n) - nbTorus, dim, mexp, usepset)
132 :     else
133 :     res <- .Call("doSFMersenneTwister", n- nbTorus, dim, mexp, usepset)
134 :    
135 :     res <- rbind(res, as.matrix( restorus, nbTorus, dim) )
136 :     }
137 :    
138 :     if(dim == 1)
139 :     as.vector(res)
140 :     else
141 :     as.matrix(res)
142 :     }
143 :    
144 :     WELL <- function(n, dim = 1, order = 512, temper = FALSE, version = "a")
145 :     {
146 :     if(n <0 || is.array(n))
147 :     stop("invalid argument 'n'")
148 :     if(dim < 0 || length(dim) >1)
149 :     stop("invalid argument 'dim'")
150 :     if(!is.numeric(order))
151 :     stop("invalid argument 'order'")
152 :     if( !(order %in% c(512, 521, 607, 800, 1024, 19937, 21701, 23209, 44497) ) )
153 :     stop("'order' must be in {512, 521, 607, 800, 1024, 19937, 21071, 23209, 44497}.")
154 :     if( !(version %in% c("a", "b") ) )
155 :     stop("'version' must be either 'a' or 'b'.")
156 :    
157 :     if(!is.logical(temper))
158 :     stop("invalid argument 'temper'")
159 :     if(temper && order %in% c(512, 521, 607, 1024))
160 :     stop("tempering impossible")
161 :    
162 :     zeversion <- 0
163 :     if(version == "a")
164 :     zeversion <- 1
165 :     if(version == "b")
166 :     zeversion <- 2
167 :     if(zeversion == 0)
168 :     stop("wrong version for WELL RNG")
169 :     if(version == "b" && order %in% c(512, 21701) )
170 :     stop("this WELL RNG does not have a 'b' version")
171 :    
172 :     if(length(n) > 1)
173 :     res <- .Call("doWELL", length(n), dim, order, temper, zeversion)
174 :     else
175 :     res <- .Call("doWELL", n, dim, order, temper, zeversion)
176 :    
177 :    
178 :    
179 :     if(dim == 1)
180 :     as.vector(res)
181 :     else
182 :     as.matrix(res)
183 :     }
184 :    
185 :     knuthTAOCP <- function(n, dim = 1)
186 :     {
187 :     if(n <0 || is.array(n))
188 :     stop("invalid argument 'n'")
189 :     if(dim < 0 || length(dim) >1)
190 :     stop("invalid argument 'dim'")
191 :    
192 :     if(length(n) > 1)
193 :     res <- .Call("doKnuthTAOCP", length(n), dim)
194 :     else
195 :     res <- .Call("doKnuthTAOCP", n, dim)
196 :    
197 :     if(dim == 1)
198 :     as.vector(res)
199 :     else
200 :     as.matrix(res)
201 :     }
202 :    

root@r-forge.r-project.org
ViewVC Help
Powered by ViewVC 1.0.0  
Thanks to:
Vienna University of Economics and Business Powered By FusionForge